| 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 | |