1- // ===-- Implementation of atanh(x) function -------------------------------===//
1+ // ===-- Half-precision atanh(x) function --- -------------------------------===//
22//
33// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
44// See https://llvm.org/LICENSE.txt for license information.
1010#include " explogxf.h"
1111#include " hdr/errno_macros.h"
1212#include " hdr/fenv_macros.h"
13+ #include " src/__support/FPUtil/FEnvImpl.h"
1314#include " src/__support/FPUtil/FPBits.h"
1415#include " src/__support/FPUtil/PolyEval.h"
1516#include " src/__support/FPUtil/cast.h"
2122namespace LIBC_NAMESPACE_DECL {
2223
2324LLVM_LIBC_FUNCTION (float16, atanhf16, (float16 x)) {
24- using FPBits = typename fputil::FPBits<float16>;
25+ using FPBits = fputil::FPBits<float16>;
2526
2627 FPBits xbits (x);
2728 Sign sign = xbits.sign ();
2829 uint16_t x_abs = xbits.abs ().uintval ();
2930
31+ // |x| >= 1
3032 if (LIBC_UNLIKELY (x_abs >= 0x3c00U )) {
31- if (xbits.is_nan ()) {
33+ if (xbits.is_nan ())
3234 return x;
33- }
35+
3436 // |x| == 1.0
3537 if (x_abs == 0x3c00U ) {
3638 fputil::set_errno_if_required (ERANGE);
3739 fputil::raise_except_if_required (FE_DIVBYZERO);
3840 return FPBits::inf (sign).get_val ();
39- } else {
40- fputil::set_errno_if_required (EDOM);
41- fputil::raise_except_if_required (FE_INVALID);
42- return FPBits::quiet_nan ().get_val ();
4341 }
42+ // |x| > 1.0
43+ fputil::set_errno_if_required (EDOM);
44+ fputil::raise_except_if_required (FE_INVALID);
45+ return FPBits::quiet_nan ().get_val ();
4446 }
4547
4648 // For |x| less than approximately 0.10
@@ -52,35 +54,25 @@ LLVM_LIBC_FUNCTION(float16, atanhf16, (float16 x)) {
5254 // atanh(x) ≈ x + (1/3)*x^3
5355 if (LIBC_UNLIKELY (x_abs < 0x0100U )) {
5456 return static_cast <float16>(
55- LIBC_UNLIKELY (x_abs == 0 ) ? x : (x + 0x1 .555556p-2 * x * x * x));
57+ LIBC_UNLIKELY (x_abs == 0 ) ? x : (x + 0x1 .555556p-2f * x * x * x));
5658 }
5759
5860 // For 0x0100U <= |x| <= 0x2e66U:
5961 // Let t = x^2.
6062 // Define P(t) ≈ (1/3)*t + (1/5)*t^2 + (1/7)*t^3 + (1/9)*t^4 + (1/11)*t^5.
61- // The coefficients below were derived using Sollya:
62- // > display = hexadecimal;
63- // > round(1/3, SG, RN);
64- // > round(1/5, SG, RN);
65- // > round(1/7, SG, RN);
66- // > round(1/9, SG, RN);
67- // > round(1/11, SG, RN);
68- // This yields:
69- // 0x1.555556p-2
70- // 0x1.99999ap-3
71- // 0x1.24924ap-3
72- // 0x1.c71c72p-4
73- // 0x1.745d18p-4f
63+ // Coefficients (from Sollya, RN, hexadecimal):
64+ // 1/3 = 0x1.555556p-2, 1/5 = 0x1.99999ap-3, 1/7 = 0x1.24924ap-3,
65+ // 1/9 = 0x1.c71c72p-4, 1/11 = 0x1.745d18p-4f
7466 // Thus, atanh(x) ≈ x * (1 + P(x^2)).
7567 float xf = x;
7668 float x2 = xf * xf;
7769 float pe = fputil::polyeval (x2, 0 .0f , 0x1 .555556p-2f , 0x1 .99999ap-3f ,
7870 0x1 .24924ap-3f , 0x1 .c71c72p -4f , 0x1 .745d18p-4f );
79- return static_cast <float16>(fputil::multiply_add (xf, pe, xf));
71+ return fputil::cast <float16>(fputil::multiply_add (xf, pe, xf));
8072 }
8173
8274 float xf = x;
83- return static_cast <float16>(0.5 * log_eval ((xf + 1.0 ) / (xf - 1.0 )));
75+ return fputil::cast <float16>(0.5 * log_eval ((xf + 1.0 ) / (xf - 1.0 )));
8476}
8577
8678} // namespace LIBC_NAMESPACE_DECL
0 commit comments