libc/src/math/generic/atanhf.cpp

   1 //===-- Single-precision atanh function -----------------------------------===//
   2 //
   3 // Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
   4 // See https://llvm.org/LICENSE.txt for license information.
   5 // SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
   6 //
   7 //===----------------------------------------------------------------------===//
   8
   9 #include "src/math/atanhf.h"
  10 #include "src/__support/FPUtil/FPBits.h"
  11 #include "src/__support/macros/optimization.h" // LIBC_UNLIKELY
  12 #include "src/math/generic/explogxf.h"
  13
  14 namespace __llvm_libc {
  15
  16 LLVM_LIBC_FUNCTION(float, atanhf, (float x)) {
  17   using FPBits = typename fputil::FPBits<float>;
  18   FPBits xbits(x);
  19   bool sign = xbits.get_sign();
  20   uint32_t x_abs = xbits.uintval() & FPBits::FloatProp::EXP_MANT_MASK;
  21
  22   // |x| >= 1.0
  23   if (LIBC_UNLIKELY(x_abs >= 0x3F80'0000U)) {
  24     if (xbits.is_nan()) {
  25       return x;
  26     }
  27     // |x| == 1.0
  28     if (x_abs == 0x3F80'0000U) {
  29       fputil::set_errno_if_required(ERANGE);
  30       fputil::raise_except_if_required(FE_DIVBYZERO);
  31       return FPBits::inf(sign).get_val();
  32     } else {
  33       fputil::set_errno_if_required(EDOM);
  34       fputil::raise_except_if_required(FE_INVALID);
  35       return FPBits::build_quiet_nan(0);
  36     }
  37   }
  38
  39   // |x| < ~0.10
  40   if (LIBC_UNLIKELY(x_abs <= 0x3dcc'0000U)) {
  41     // |x| <= 2^-26
  42     if (LIBC_UNLIKELY(x_abs <= 0x3280'0000U)) {
  43       return static_cast<float>(LIBC_UNLIKELY(x_abs == 0)
  44                                     ? x
  45                                     : (x + 0x1.5555555555555p-2 * x * x * x));
  46     }
  47
  48     double xdbl = x;
  49     double x2 = xdbl * xdbl;
  50     // Pure Taylor series.
  51     double pe = fputil::polyeval(x2, 0.0, 0x1.5555555555555p-2,
  52                                  0x1.999999999999ap-3, 0x1.2492492492492p-3,
  53                                  0x1.c71c71c71c71cp-4, 0x1.745d1745d1746p-4);
  54     return static_cast<float>(fputil::multiply_add(xdbl, pe, xdbl));
  55   }
  56   double xdbl = x;
  57   return static_cast<float>(0.5 * log_eval((xdbl + 1.0) / (xdbl - 1.0)));
  58 }
  59
  60 } // namespace __llvm_libc