1 | //===-- Half-precision atanf16(x) 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/atanf16.h" |
10 | #include "hdr/errno_macros.h" |
11 | #include "hdr/fenv_macros.h" |
12 | #include "src/__support/FPUtil/FEnvImpl.h" |
13 | #include "src/__support/FPUtil/FPBits.h" |
14 | #include "src/__support/FPUtil/PolyEval.h" |
15 | #include "src/__support/FPUtil/cast.h" |
16 | #include "src/__support/FPUtil/except_value_utils.h" |
17 | #include "src/__support/FPUtil/multiply_add.h" |
18 | #include "src/__support/FPUtil/sqrt.h" |
19 | #include "src/__support/macros/optimization.h" |
20 | |
21 | namespace LIBC_NAMESPACE_DECL { |
22 | |
23 | // Generated by Solly using the following command: |
24 | // > round(pi/2, SG, RN); |
25 | static constexpr float PI_2 = 0x1.921fb6p0; |
26 | |
27 | #ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS |
28 | static constexpr size_t N_EXCEPTS = 6; |
29 | |
30 | static constexpr fputil::ExceptValues<float16, N_EXCEPTS> ATANF16_EXCEPTS{{ |
31 | // (input, RZ output, RU offset, RD offset, RN offset) |
32 | {0x2745, 0x2744, 1, 0, 1}, |
33 | {0x3099, 0x3090, 1, 0, 1}, |
34 | {0x3c6c, 0x3aae, 1, 0, 1}, |
35 | {0x466e, 0x3daa, 1, 0, 1}, |
36 | {0x48ae, 0x3ddb, 1, 0, 0}, |
37 | {0x5619, 0x3e3d, 1, 0, 1}, |
38 | }}; |
39 | #endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS |
40 | |
41 | LLVM_LIBC_FUNCTION(float16, atanf16, (float16 x)) { |
42 | using FPBits = fputil::FPBits<float16>; |
43 | FPBits xbits(x); |
44 | |
45 | uint16_t x_u = xbits.uintval(); |
46 | uint16_t x_abs = x_u & 0x7fff; |
47 | bool x_sign = x_u >> 15; |
48 | float sign = (x_sign ? -1.0 : 1.0); |
49 | |
50 | // |x| >= +/-inf |
51 | if (LIBC_UNLIKELY(x_abs >= 0x7c00)) { |
52 | if (xbits.is_nan()) { |
53 | if (xbits.is_signaling_nan()) { |
54 | fputil::raise_except_if_required(FE_INVALID); |
55 | return FPBits::quiet_nan().get_val(); |
56 | } |
57 | return x; |
58 | } |
59 | |
60 | // atanf16(+/-inf) = +/-pi/2 |
61 | return fputil::cast<float16>(sign * PI_2); |
62 | } |
63 | |
64 | float xf = x; |
65 | float xsq = xf * xf; |
66 | #ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS |
67 | // Handle exceptional values |
68 | if (auto r = ATANF16_EXCEPTS.lookup_odd(x_abs, x_sign); |
69 | LIBC_UNLIKELY(r.has_value())) |
70 | return r.value(); |
71 | #endif |
72 | |
73 | // |x| <= 0x1p0, |x| <= 1 |
74 | if (x_abs <= 0x3c00) { |
75 | // atanf16(+/-0) = +/-0 |
76 | if (LIBC_UNLIKELY(x_abs == 0)) |
77 | return x; |
78 | |
79 | // Degree-14 minimax odd polynomial of atan(x) generated by Sollya with: |
80 | // > P = fpminimax(atan(x)/x, [|0, 2, 4, 6, 8, 10, 12, 14|], [|SG...|], |
81 | // [0, 1]); |
82 | float result = fputil::polyeval( |
83 | xsq, 0x1.fffffcp-1f, -0x1.55519ep-2f, 0x1.98f6a8p-3f, -0x1.1f0a92p-3f, |
84 | 0x1.95b654p-4f, -0x1.e65492p-5f, 0x1.8c0c36p-6f, -0x1.32316ep-8f); |
85 | return fputil::cast<float16>(xf * result); |
86 | } |
87 | |
88 | // If |x| > 1 |
89 | // y = atan(x) = sign(x) * atan(|x|) |
90 | // atan(|x|) = pi/2 - atan(1/|x|) |
91 | // Recall, 1/|x| < 1 |
92 | float x_inv_sq = 1.0f / xsq; |
93 | float x_inv = fputil::sqrt<float>(x_inv_sq); |
94 | |
95 | // Degree-14 minimax odd polynomial of atan(x) generated by Sollya with: |
96 | // > P = fpminimax(atan(x)/x, [|0, 2, 4, 6, 8, 10, 12, 14|], [|SG...|], |
97 | // [0, 1]); |
98 | float interm = |
99 | fputil::polyeval(x_inv_sq, 0x1.fffffcp-1f, -0x1.55519ep-2f, |
100 | 0x1.98f6a8p-3f, -0x1.1f0a92p-3f, 0x1.95b654p-4f, |
101 | -0x1.e65492p-5f, 0x1.8c0c36p-6f, -0x1.32316ep-8f); |
102 | |
103 | return fputil::cast<float16>(sign * |
104 | fputil::multiply_add(x_inv, -interm, PI_2)); |
105 | } |
106 | |
107 | } // namespace LIBC_NAMESPACE_DECL |
108 | |