1//===----------------------------------------------------------------------===//
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/// \file
10/// Common utilities and tables for single-precision powf(x, y).
11///
12//===----------------------------------------------------------------------===//
13
14#ifndef LLVM_LIBC_SRC___SUPPORT_MATH_POWF_UTILS_H
15#define LLVM_LIBC_SRC___SUPPORT_MATH_POWF_UTILS_H
16
17#include "hdr/errno_macros.h"
18#include "hdr/fenv_macros.h"
19#include "src/__support/CPP/bit.h"
20#include "src/__support/CPP/optional.h"
21#include "src/__support/FPUtil/FEnvImpl.h"
22#include "src/__support/FPUtil/FPBits.h"
23#include "src/__support/FPUtil/double_double.h"
24#include "src/__support/FPUtil/nearest_integer.h"
25#include "src/__support/FPUtil/rounding_mode.h"
26#include "src/__support/FPUtil/triple_double.h"
27#include "src/__support/common.h"
28#include "src/__support/macros/attributes.h"
29#include "src/__support/macros/config.h"
30#include "src/__support/macros/optimization.h"
31#include "src/__support/macros/properties/cpu_features.h"
32#include "src/__support/math/common_constants.h"
33
34#if defined(LIBC_TARGET_CPU_HAS_FPU_FLOAT) || \
35 !defined(LIBC_MATH_HAS_SMALL_TABLES)
36#include "src/__support/FPUtil/sqrt.h"
37#endif // LIBC_TARGET_CPU_HAS_FPU_FLOAT || !LIBC_MATH_HAS_SMALL_TABLES
38
39#ifndef LIBC_MATH_HAS_SMALL_TABLES
40#include "src/__support/math/exp10f.h"
41#include "src/__support/math/exp2f.h"
42#endif // !LIBC_MATH_HAS_SMALL_TABLES
43
44namespace LIBC_NAMESPACE_DECL {
45namespace math {
46namespace powf_internal {
47
48using fputil::DoubleDouble;
49using fputil::FloatFloat;
50using fputil::TripleDouble;
51
52// Table of -log2(R[i]) represented in triple-double precision {lo, mid, hi},
53// where R[i] = 2^-8 * ceil(2^8 * (1 - 2^-8) / (1 + i * 2^-7)) for i = 0..127.
54//
55// The least significant bit of the high part `hi` is chosen to be >= 2^-20,
56// so that the total precision of (e_x + hi) is at most 8 + 20 = 28 bits
57// (since |e_x| <= 150 < 2^8 fits in 8 bits).
58// Because the single-precision exponent y has 24 bits of mantissa, the product
59// y * (e_x + hi) spans at most 24 + 28 = 52 bits <= 53 bits, making it exact
60// in double precision.
61// The mid part provides 53 bits for fast evaluation, and the lo part provides
62// another 53 bits for the accurate fallback, yielding >= 127 bits of precision.
63//
64// Generated by Sollya with:
65// > display = hexadecimal;
66// > prec = 200;
67// > for i from 0 to 127 do {
68// if (i == 0) then {
69// print("{0.0, 0.0, 0.0},");
70// } else if (i == 127) then {
71// print("{0.0, 0.0, 1.0},");
72// } else {
73// r = 2^(-8) * ceil( 2^8 * (1 - 2^(-8)) / (1 + i*2^(-7)) );
74// l = -log2(r);
75// h = nearestint(l * 2^20) * 2^(-20);
76// m = round(l - h, D, RN);
77// lo = round(l - h - m, D, RN);
78// print("{", lo, ", ", m, ", ", h, "},");
79// };
80// };
81LIBC_INLINE_VAR constexpr TripleDouble LOG2_R_TD[128] = {
82 {.lo: 0.0, .mid: 0.0, .hi: 0.0},
83 {.lo: 0x1.84a2c615b70adp-79, .mid: -0x1.177c23362928cp-25, .hi: 0x1.72c8p-7},
84 {.lo: -0x1.f27b820fd03eap-76, .mid: -0x1.179e0caa9c9abp-22, .hi: 0x1.744p-6},
85 {.lo: -0x1.f27ef487c8f34p-77, .mid: -0x1.c6cea541f5b7p-23, .hi: 0x1.184cp-5},
86 {.lo: -0x1.e3f80fbc71454p-76, .mid: -0x1.66c4d4e554434p-22, .hi: 0x1.773ap-5},
87 {.lo: -0x1.9f8ef14d5f6eep-79, .mid: -0x1.70700a00fdd55p-24, .hi: 0x1.d6ecp-5},
88 {.lo: 0x1.452bbce7398c1p-77, .mid: 0x1.53002a4e86631p-23, .hi: 0x1.1bb3p-4},
89 {.lo: -0x1.990555535afdp-81, .mid: 0x1.fcd15f101c142p-25, .hi: 0x1.4c56p-4},
90 {.lo: 0x1.447e30ad393eep-78, .mid: 0x1.25b3eed319cedp-22, .hi: 0x1.7d6p-4},
91 {.lo: 0x1.b7759da88a2dap-76, .mid: -0x1.4195120d8486fp-22, .hi: 0x1.960dp-4},
92 {.lo: 0x1.cee7766ece702p-78, .mid: 0x1.45b878e27d0d9p-23, .hi: 0x1.c7b5p-4},
93 {.lo: -0x1.a55c745ecdc2fp-77, .mid: 0x1.770744593a4cbp-22, .hi: 0x1.f9c9p-4},
94 {.lo: 0x1.f7ec992caa67fp-77, .mid: 0x1.c673032495d24p-22, .hi: 0x1.097ep-3},
95 {.lo: -0x1.433638c6ece3ep-77, .mid: -0x1.1eaa65b49696ep-22, .hi: 0x1.22dbp-3},
96 {.lo: 0x1.58f27b6518824p-76, .mid: 0x1.b2866f2850b22p-22, .hi: 0x1.3c6f8p-3},
97 {.lo: -0x1.86bdcfdfd4a4cp-79, .mid: 0x1.8ee37cd2ea9d3p-25, .hi: 0x1.494f8p-3},
98 {.lo: -0x1.ff7044a68a7fap-80, .mid: 0x1.7e86f9c2154fbp-24, .hi: 0x1.633a8p-3},
99 {.lo: -0x1.aa21694561327p-81, .mid: 0x1.8e3cfc25f0ce6p-26, .hi: 0x1.7046p-3},
100 {.lo: -0x1.d209f2d4239c6p-87, .mid: 0x1.57f7a64ccd537p-28, .hi: 0x1.8a898p-3},
101 {.lo: -0x1.a55e97e60e632p-76, .mid: -0x1.a761c09fbd2aep-22, .hi: 0x1.97c2p-3},
102 {.lo: 0x1.261179225541ep-76, .mid: 0x1.24bea9a2c66f3p-22, .hi: 0x1.b26p-3},
103 {.lo: -0x1.08fa30510fca9p-82, .mid: -0x1.60002ccfe43f5p-25, .hi: 0x1.bfc68p-3},
104 {.lo: -0x1.63ec8d56242f9p-76, .mid: 0x1.69f220e97f22cp-22, .hi: 0x1.dac2p-3},
105 {.lo: 0x1.8bcdaf0534365p-76, .mid: -0x1.6164f64c210ep-22, .hi: 0x1.e858p-3},
106 {.lo: 0x1.1003282896056p-78, .mid: -0x1.0c1678ae89767p-24, .hi: 0x1.01d9cp-2},
107 {.lo: 0x1.01bcc7025fa92p-78, .mid: -0x1.f26a05c813d57p-22, .hi: 0x1.08bdp-2},
108 {.lo: -0x1.fe8a8648e9ebcp-80, .mid: 0x1.4d8fc561c8d44p-24, .hi: 0x1.169cp-2},
109 {.lo: 0x1.08dfb23650c75p-79, .mid: -0x1.362ad8f7ca2dp-22, .hi: 0x1.1d984p-2},
110 {.lo: -0x1.f8d5a89861a5ep-79, .mid: 0x1.2b13cd6c4d042p-22, .hi: 0x1.249ccp-2},
111 {.lo: -0x1.a1c872983511ep-76, .mid: -0x1.1c8f11979a5dbp-22, .hi: 0x1.32cp-2},
112 {.lo: 0x1.e8e21bff3336bp-77, .mid: 0x1.c2ab3edefe569p-23, .hi: 0x1.39de8p-2},
113 {.lo: 0x1.fd1994fb2c4a1p-80, .mid: 0x1.7c3eca28e69cap-26, .hi: 0x1.4106p-2},
114 {.lo: 0x1.6b94b51cf76b1p-80, .mid: -0x1.34c4e99e1c6c6p-24, .hi: 0x1.4f6fcp-2},
115 {.lo: -0x1.31d55da1d0f66p-76, .mid: -0x1.194a871b63619p-22, .hi: 0x1.56b24p-2},
116 {.lo: -0x1.378b22691e28bp-77, .mid: 0x1.e3dd5c1c885aep-23, .hi: 0x1.5dfdcp-2},
117 {.lo: 0x1.99e302970e411p-83, .mid: -0x1.6ccf3b1129b7cp-23, .hi: 0x1.6552cp-2},
118 {.lo: 0x1.20164a049664dp-82, .mid: -0x1.2f346e2bf924bp-23, .hi: 0x1.6cb1p-2},
119 {.lo: -0x1.d14aac4d864c3p-77, .mid: -0x1.fa61aaa59c1d8p-23, .hi: 0x1.7b8ap-2},
120 {.lo: 0x1.496ab4e4b293fp-79, .mid: 0x1.90c11fd32a3abp-22, .hi: 0x1.8304cp-2},
121 {.lo: -0x1.d209f2d4239c6p-86, .mid: 0x1.57f7a64ccd537p-27, .hi: 0x1.8a898p-2},
122 {.lo: 0x1.eae3326327babp-81, .mid: 0x1.249ba76fee235p-27, .hi: 0x1.9218p-2},
123 {.lo: 0x1.fa05bddfded8cp-77, .mid: -0x1.aad2729b21ae5p-23, .hi: 0x1.99b08p-2},
124 {.lo: -0x1.624140d175ba2p-77, .mid: 0x1.71810a5e1818p-22, .hi: 0x1.a8ff8p-2},
125 {.lo: 0x1.f1c5160c515c1p-81, .mid: -0x1.6172fe015e13cp-27, .hi: 0x1.b0b68p-2},
126 {.lo: -0x1.86a6204eec8cp-79, .mid: 0x1.5ec6c1bfbf89ap-24, .hi: 0x1.b877cp-2},
127 {.lo: 0x1.718f761dd3915p-78, .mid: 0x1.678bf6cdedf51p-24, .hi: 0x1.c0438p-2},
128 {.lo: -0x1.d4ee66c3700e4p-76, .mid: 0x1.c2d45fe43895ep-22, .hi: 0x1.c819cp-2},
129 {.lo: -0x1.7d14533586306p-77, .mid: -0x1.9ee52ed49d71dp-22, .hi: 0x1.cffbp-2},
130 {.lo: 0x1.5ce9fb5a7bb5bp-81, .mid: 0x1.5786af187a96bp-27, .hi: 0x1.d7e6cp-2},
131 {.lo: -0x1.ae6face57ad3bp-77, .mid: 0x1.3ab0dc56138c9p-23, .hi: 0x1.dfdd8p-2},
132 {.lo: 0x1.5ac93b443d55fp-78, .mid: 0x1.fe538ab34efb5p-22, .hi: 0x1.e7df4p-2},
133 {.lo: 0x1.f1753e0ae1e8fp-76, .mid: -0x1.e4fee07aa4b68p-22, .hi: 0x1.efec8p-2},
134 {.lo: 0x1.cdfd4c297069bp-76, .mid: -0x1.172f32fe67287p-22, .hi: 0x1.f804cp-2},
135 {.lo: 0x1.97a0e8f3ba742p-79, .mid: -0x1.9a83ff9ab9cc8p-22, .hi: 0x1.00144p-1},
136 {.lo: -0x1.800450f5b2357p-78, .mid: -0x1.68cb06cece193p-22, .hi: 0x1.042bep-1},
137 {.lo: -0x1.a839041241fe7p-78, .mid: 0x1.8cd71ddf82e2p-22, .hi: 0x1.08494p-1},
138 {.lo: 0x1.ed0b8eeccca86p-78, .mid: 0x1.5e18ab2df3ae6p-22, .hi: 0x1.0c6cap-1},
139 {.lo: 0x1.3dd41df9689b3p-79, .mid: 0x1.5dee4d9d8a273p-25, .hi: 0x1.1096p-1},
140 {.lo: -0x1.990555535afdp-82, .mid: 0x1.fcd15f101c142p-26, .hi: 0x1.14c56p-1},
141 {.lo: -0x1.1773d02c9055cp-77, .mid: -0x1.2474b0f992ba1p-23, .hi: 0x1.18faep-1},
142 {.lo: -0x1.4aeef330c53c1p-78, .mid: 0x1.4b5a92a606047p-24, .hi: 0x1.1d368p-1},
143 {.lo: 0x1.8e6ff749ebacbp-77, .mid: 0x1.16186fcf54bbdp-22, .hi: 0x1.21786p-1},
144 {.lo: 0x1.c09d761c548ebp-84, .mid: 0x1.18efabeb7d722p-27, .hi: 0x1.25c0ap-1},
145 {.lo: 0x1.aaa73a428e1e4p-78, .mid: -0x1.e5fc7d238691dp-24, .hi: 0x1.2a0f4p-1},
146 {.lo: -0x1.af2f3d8b63fbap-79, .mid: 0x1.f5809faf6283cp-22, .hi: 0x1.2e644p-1},
147 {.lo: -0x1.af2f3d8b63fbap-79, .mid: 0x1.f5809faf6283cp-22, .hi: 0x1.2e644p-1},
148 {.lo: 0x1.78de359f2bb88p-77, .mid: 0x1.c6e1dcd0cb449p-22, .hi: 0x1.32bfep-1},
149 {.lo: -0x1.415ae1a715618p-76, .mid: 0x1.76e0e8f74b4d5p-22, .hi: 0x1.37222p-1},
150 {.lo: -0x1.4991b5375621fp-79, .mid: -0x1.cb82c89692d99p-24, .hi: 0x1.3b8b2p-1},
151 {.lo: -0x1.827d37deb2236p-76, .mid: -0x1.63161c5432aebp-22, .hi: 0x1.3ffaep-1},
152 {.lo: 0x1.9576edac01c78p-77, .mid: 0x1.458104c41b901p-22, .hi: 0x1.44716p-1},
153 {.lo: 0x1.9576edac01c78p-77, .mid: 0x1.458104c41b901p-22, .hi: 0x1.44716p-1},
154 {.lo: -0x1.05a27b81e2219p-77, .mid: -0x1.cd9d0cde578d5p-22, .hi: 0x1.48efp-1},
155 {.lo: 0x1.237616778b4bap-82, .mid: 0x1.b9884591add87p-26, .hi: 0x1.4d738p-1},
156 {.lo: 0x1.3b7d7e5d148bbp-76, .mid: 0x1.c6042978605ffp-22, .hi: 0x1.51ff2p-1},
157 {.lo: -0x1.cc3f936a5977cp-79, .mid: -0x1.fc4c96b37dcf6p-22, .hi: 0x1.56922p-1},
158 {.lo: 0x1.20164a049664dp-83, .mid: -0x1.2f346e2bf924bp-24, .hi: 0x1.5b2c4p-1},
159 {.lo: 0x1.20164a049664dp-83, .mid: -0x1.2f346e2bf924bp-24, .hi: 0x1.5b2c4p-1},
160 {.lo: -0x1.a212919a92f7ap-77, .mid: 0x1.c4e4fbb68a4d1p-22, .hi: 0x1.5fcdcp-1},
161 {.lo: -0x1.b64b03f7230ddp-77, .mid: -0x1.9d499bd9b3226p-23, .hi: 0x1.6476ep-1},
162 {.lo: -0x1.1ec6379e6e3b9p-77, .mid: -0x1.f89b355ede26fp-23, .hi: 0x1.69278p-1},
163 {.lo: -0x1.1ec6379e6e3b9p-77, .mid: -0x1.f89b355ede26fp-23, .hi: 0x1.69278p-1},
164 {.lo: -0x1.4ba44c03bfbbdp-78, .mid: 0x1.53c7e319f6e92p-24, .hi: 0x1.6ddfcp-1},
165 {.lo: -0x1.c36fc650d030fp-77, .mid: -0x1.b291f070528c7p-22, .hi: 0x1.729fep-1},
166 {.lo: -0x1.69e5693a7f067p-80, .mid: 0x1.2967a451a7b48p-25, .hi: 0x1.7767cp-1},
167 {.lo: -0x1.69e5693a7f067p-80, .mid: 0x1.2967a451a7b48p-25, .hi: 0x1.7767cp-1},
168 {.lo: 0x1.6598aae91499ap-76, .mid: 0x1.244fcff690fcep-22, .hi: 0x1.7c37ap-1},
169 {.lo: 0x1.99d61ec432837p-77, .mid: 0x1.46fd97f5dc572p-23, .hi: 0x1.810fap-1},
170 {.lo: 0x1.99d61ec432837p-77, .mid: 0x1.46fd97f5dc572p-23, .hi: 0x1.810fap-1},
171 {.lo: 0x1.855c42078f81bp-76, .mid: -0x1.f3a7352663e5p-22, .hi: 0x1.85efep-1},
172 {.lo: -0x1.59408e815107p-77, .mid: 0x1.b3cda690370b5p-23, .hi: 0x1.8ad84p-1},
173 {.lo: -0x1.59408e815107p-77, .mid: 0x1.b3cda690370b5p-23, .hi: 0x1.8ad84p-1},
174 {.lo: 0x1.33b318085e50ap-78, .mid: 0x1.3226b211bf1d9p-23, .hi: 0x1.8fc92p-1},
175 {.lo: 0x1.343fe7c9cb4aep-79, .mid: 0x1.d24b136c101eep-23, .hi: 0x1.94c28p-1},
176 {.lo: 0x1.343fe7c9cb4aep-79, .mid: 0x1.d24b136c101eep-23, .hi: 0x1.94c28p-1},
177 {.lo: -0x1.d19522e56fe6p-76, .mid: 0x1.7c40c7907e82ap-22, .hi: 0x1.99c48p-1},
178 {.lo: -0x1.23b9d8ea55c3ep-77, .mid: -0x1.e81781d97ee91p-22, .hi: 0x1.9ecf6p-1},
179 {.lo: -0x1.23b9d8ea55c3ep-77, .mid: -0x1.e81781d97ee91p-22, .hi: 0x1.9ecf6p-1},
180 {.lo: 0x1.829440c24aeb6p-78, .mid: -0x1.6a77813f94e01p-22, .hi: 0x1.a3e3p-1},
181 {.lo: -0x1.624140d175ba2p-76, .mid: -0x1.1cfdeb43cfdp-22, .hi: 0x1.a8ffap-1},
182 {.lo: -0x1.624140d175ba2p-76, .mid: -0x1.1cfdeb43cfdp-22, .hi: 0x1.a8ffap-1},
183 {.lo: 0x1.afa6f024fb045p-77, .mid: -0x1.f983f74d3138fp-23, .hi: 0x1.ae256p-1},
184 {.lo: -0x1.603ad3a5d326dp-78, .mid: -0x1.e278ae1a1f51fp-23, .hi: 0x1.b3546p-1},
185 {.lo: -0x1.603ad3a5d326dp-78, .mid: -0x1.e278ae1a1f51fp-23, .hi: 0x1.b3546p-1},
186 {.lo: -0x1.0c1e0e5855d6ap-77, .mid: -0x1.97552b7b5ea45p-23, .hi: 0x1.b88ccp-1},
187 {.lo: -0x1.0c1e0e5855d6ap-77, .mid: -0x1.97552b7b5ea45p-23, .hi: 0x1.b88ccp-1},
188 {.lo: 0x1.c817ad56baa16p-78, .mid: -0x1.19b4f3c72c4f8p-24, .hi: 0x1.bdceap-1},
189 {.lo: 0x1.44c47ac1bf62bp-77, .mid: 0x1.f7402d26f1a12p-23, .hi: 0x1.c31a2p-1},
190 {.lo: 0x1.44c47ac1bf62bp-77, .mid: 0x1.f7402d26f1a12p-23, .hi: 0x1.c31a2p-1},
191 {.lo: -0x1.69b9465eae1e6p-78, .mid: -0x1.2056d5dd31d96p-23, .hi: 0x1.c86f8p-1},
192 {.lo: -0x1.69b9465eae1e6p-78, .mid: -0x1.2056d5dd31d96p-23, .hi: 0x1.c86f8p-1},
193 {.lo: -0x1.24a6d9d1d1904p-79, .mid: -0x1.6e46335aae723p-24, .hi: 0x1.cdcecp-1},
194 {.lo: -0x1.3826144575ac4p-76, .mid: -0x1.beb244c59f331p-22, .hi: 0x1.d3382p-1},
195 {.lo: -0x1.3826144575ac4p-76, .mid: -0x1.beb244c59f331p-22, .hi: 0x1.d3382p-1},
196 {.lo: 0x1.dbc96b3b12b25p-81, .mid: 0x1.16c071e93fd97p-27, .hi: 0x1.d8abap-1},
197 {.lo: 0x1.dbc96b3b12b25p-81, .mid: 0x1.16c071e93fd97p-27, .hi: 0x1.d8abap-1},
198 {.lo: 0x1.68a8ccdbd1f33p-77, .mid: 0x1.d8175819530c2p-22, .hi: 0x1.de298p-1},
199 {.lo: 0x1.68a8ccdbd1f33p-77, .mid: 0x1.d8175819530c2p-22, .hi: 0x1.de298p-1},
200 {.lo: 0x1.e586711df5ea1p-79, .mid: 0x1.51bd552842c1cp-23, .hi: 0x1.e3b2p-1},
201 {.lo: 0x1.e586711df5ea1p-79, .mid: 0x1.51bd552842c1cp-23, .hi: 0x1.e3b2p-1},
202 {.lo: -0x1.bc25adf042483p-79, .mid: 0x1.914e204f19d94p-22, .hi: 0x1.e9452p-1},
203 {.lo: -0x1.bc25adf042483p-79, .mid: 0x1.914e204f19d94p-22, .hi: 0x1.e9452p-1},
204 {.lo: 0x1.d7d82b65c5686p-76, .mid: 0x1.c55d997da24fdp-22, .hi: 0x1.eee32p-1},
205 {.lo: 0x1.d7d82b65c5686p-76, .mid: 0x1.c55d997da24fdp-22, .hi: 0x1.eee32p-1},
206 {.lo: -0x1.3f108c0857ca3p-77, .mid: -0x1.685c2d2298a6ep-22, .hi: 0x1.f48c4p-1},
207 {.lo: -0x1.3f108c0857ca3p-77, .mid: -0x1.685c2d2298a6ep-22, .hi: 0x1.f48c4p-1},
208 {.lo: -0x1.bd800bca7a221p-78, .mid: 0x1.7a4887bd74039p-22, .hi: 0x1.fa406p-1},
209 {.lo: 0.0, .mid: 0.0, .hi: 1.0},
210};
211
212// Lookup table for 5-bit (32 entries) range reduction for float_eval and
213// double_eval fast path:
214// r(0) = 1.0f
215// r(i) = 2^-6 * ceil(2^6 * (1 - 2^-6) / (1 + i * 2^-5)), i = 1..31.
216// The constants are chosen so that dx = r * m_x - 1 is exact in
217// single precision (with hardware FMA or 14-bit Sterbenz split), and
218// -2^-6 <= dx < 2^-5.
219//
220// Generated by Sollya with:
221// > display = hexadecimal;
222// > for i from 0 to 31 do {
223// if (i == 0) then {
224// print("0x1.0p+0f,");
225// } else {
226// r = 2^(-6) * ceil( 2^6 * (1 - 2^(-6)) / (1 + i*2^(-5)) );
227// print(r @ "f,");
228// };
229// };
230LIBC_INLINE_VAR constexpr float R_32[32] = {
231 0x1.0p+0f, 0x1.fp-1f, 0x1.ep-1f, 0x1.dp-1f, 0x1.cp-1f, 0x1.b8p-1f,
232 0x1.bp-1f, 0x1.ap-1f, 0x1.98p-1f, 0x1.9p-1f, 0x1.8p-1f, 0x1.78p-1f,
233 0x1.7p-1f, 0x1.68p-1f, 0x1.6p-1f, 0x1.58p-1f, 0x1.5p-1f, 0x1.5p-1f,
234 0x1.48p-1f, 0x1.4p-1f, 0x1.38p-1f, 0x1.38p-1f, 0x1.3p-1f, 0x1.28p-1f,
235 0x1.2p-1f, 0x1.2p-1f, 0x1.18p-1f, 0x1.18p-1f, 0x1.1p-1f, 0x1.1p-1f,
236 0x1.08p-1f, 0x1.0p-1f,
237};
238
239// R_32 represented in double precision for double_eval.
240LIBC_INLINE_VAR constexpr double R_32_D[32] = {
241 0x1.0p+0, 0x1.fp-1, 0x1.ep-1, 0x1.dp-1, 0x1.cp-1, 0x1.b8p-1, 0x1.bp-1,
242 0x1.ap-1, 0x1.98p-1, 0x1.9p-1, 0x1.8p-1, 0x1.78p-1, 0x1.7p-1, 0x1.68p-1,
243 0x1.6p-1, 0x1.58p-1, 0x1.5p-1, 0x1.5p-1, 0x1.48p-1, 0x1.4p-1, 0x1.38p-1,
244 0x1.38p-1, 0x1.3p-1, 0x1.28p-1, 0x1.2p-1, 0x1.2p-1, 0x1.18p-1, 0x1.18p-1,
245 0x1.1p-1, 0x1.1p-1, 0x1.08p-1, 0x1.0p-1,
246};
247
248// Table of -log2(R_32[i]) represented in double precision for double_eval.
249//
250// Generated by Sollya with:
251// > display = hexadecimal;
252// > prec = 200;
253// > for i from 0 to 31 do {
254// if (i == 0) then {
255// print("0x0.0000000000000p+0,");
256// } else if (i == 31) then {
257// print("0x1.0000000000000p+0,");
258// } else {
259// r = 2^(-6) * ceil( 2^6 * (1 - 2^(-6)) / (1 + i*2^(-5)) );
260// print(round(-log2(r), D, RN), ",");
261// };
262// };
263LIBC_INLINE_VAR constexpr double LOG2_R_32[32] = {
264 0x0.0000000000000p+0, 0x1.77394c9d958d5p-5, 0x1.7d60496cfbb4cp-4,
265 0x1.22dadc2ab3497p-3, 0x1.8a8980abfbd32p-3, 0x1.bfc67a7fff4ccp-3,
266 0x1.f5fd8a9063e35p-3, 0x1.32bfee370ee68p-2, 0x1.4f6fbb2cec598p-2,
267 0x1.6cb0f6865c8eap-2, 0x1.a8ff971810a5ep-2, 0x1.c819dc2d45fe4p-2,
268 0x1.e7df5fe538ab3p-2, 0x1.042bd4b9a7c99p-1, 0x1.14c560fe68af9p-1,
269 0x1.25c0a0463bebp-1, 0x1.37222bb70747cp-1, 0x1.37222bb70747cp-1,
270 0x1.48eef19317991p-1, 0x1.5b2c3da19723bp-1, 0x1.6ddfc2a78fc63p-1,
271 0x1.6ddfc2a78fc63p-1, 0x1.810fa51bf65fdp-1, 0x1.94c287492c4dbp-1,
272 0x1.a8ff971810a5ep-1, 0x1.a8ff971810a5ep-1, 0x1.bdce9dcc96187p-1,
273 0x1.bdce9dcc96187p-1, 0x1.d338120a6dd9dp-1, 0x1.d338120a6dd9dp-1,
274 0x1.e9452c8a71028p-1, 0x1.0000000000000p+0,
275};
276
277// Table of -log2(R_32[i]) represented in FloatFloat precision {lo, hi} for
278// float_eval.
279// We choose the precision of the high part to be 24 - 8 = 16 bits, so that
280// e_xf + LOG2_R_FF_32[i].hi
281// is exact in float for |e_x| <= 150.
282// For i = 31, R_32[31] = 0.5, so -log2(R_32[31]) = 1.0f is exact.
283//
284// Generated by Sollya with:
285// > display = hexadecimal;
286// > prec = 200;
287// > for i from 0 to 31 do {
288// if (i == 0) then {
289// print("{0x0.000000p+0f, 0x0.000000p+0f},");
290// } else if (i == 31) then {
291// print("{0x0.000000p+0f, 0x1.000000p+0f},");
292// } else {
293// r = 2^(-6) * ceil( 2^6 * (1 - 2^(-6)) / (1 + i*2^(-5)) );
294// l = -log2(r);
295// h = round(1 + l, 17, RN) - 1;
296// lo = round(l - h, SG, RN);
297// print("{" @ lo @ "f, " @ h @ "f},");
298// };
299// };
300LIBC_INLINE_VAR constexpr FloatFloat LOG2_R_FF_32[32] = {
301 {.lo: 0x0.000000p+0f, .hi: 0x0.000000p+0f}, {.lo: -0x1.acd89ap-19f, .hi: 0x1.774p-5f},
302 {.lo: 0x1.25b3eep-22f, .hi: 0x1.7d6p-4f}, {.lo: 0x1.6e155ap-18f, .hi: 0x1.22d8p-3f},
303 {.lo: 0x1.80abfcp-19f, .hi: 0x1.8a88p-3f}, {.lo: -0x1.858p-19f, .hi: 0x1.bfc8p-3f},
304 {.lo: -0x1.3ab7cep-18f, .hi: 0x1.f6p-3f}, {.lo: -0x1.1c8f12p-22f, .hi: 0x1.32cp-2f},
305 {.lo: -0x1.134c4ep-20f, .hi: 0x1.4f7p-2f}, {.lo: 0x1.ed0cbap-19f, .hi: 0x1.6cbp-2f},
306 {.lo: -0x1.a39fbep-20f, .hi: 0x1.a9p-2f}, {.lo: 0x1.dc2d46p-18f, .hi: 0x1.c818p-2f},
307 {.lo: -0x1.40358ep-19f, .hi: 0x1.e7ep-2f}, {.lo: -0x1.5a32c2p-20f, .hi: 0x1.042cp-1f},
308 {.lo: -0x1.3e032ep-18f, .hi: 0x1.14c6p-1f}, {.lo: 0x1.408c78p-18f, .hi: 0x1.25c0p-1f},
309 {.lo: 0x1.5db83ap-20f, .hi: 0x1.3722p-1f}, {.lo: 0x1.5db83ap-20f, .hi: 0x1.3722p-1f},
310 {.lo: 0x1.e3263p-18f, .hi: 0x1.48eep-1f}, {.lo: 0x1.ed0cbap-20f, .hi: 0x1.5b2cp-1f},
311 {.lo: -0x1.eac382p-20f, .hi: 0x1.6de0p-1f}, {.lo: -0x1.eac382p-20f, .hi: 0x1.6de0p-1f},
312 {.lo: -0x1.6b9026p-19f, .hi: 0x1.811p-1f}, {.lo: 0x1.0e9258p-18f, .hi: 0x1.94c2p-1f},
313 {.lo: -0x1.a39fbep-19f, .hi: 0x1.a9p-1f}, {.lo: -0x1.a39fbep-19f, .hi: 0x1.a9p-1f},
314 {.lo: 0x1.3b992cp-18f, .hi: 0x1.bdcep-1f}, {.lo: 0x1.3b992cp-18f, .hi: 0x1.bdcep-1f},
315 {.lo: 0x1.20a6dep-21f, .hi: 0x1.d338p-1f}, {.lo: 0x1.20a6dep-21f, .hi: 0x1.d338p-1f},
316 {.lo: -0x1.a6eb1ep-18f, .hi: 0x1.e946p-1f}, {.lo: 0x0.000000p+0f, .hi: 0x1.000000p+0f},
317};
318
319// Upper bound for y = 150 / |log2(1 - 2^-24)|, generated by Sollya:
320// > y = round(-150 / log2(1 - 2^(-24)), SG, RU);
321// > y;
322// 0x1.9fe368p30
323// > printsingle(y);
324// 0x4ecff1b4
325constexpr uint32_t Y_UPPER_BOUND = 0x4ecf'f1b4;
326// Lower bound for y = 2^-25 / 150, generated by Sollya:
327// > y = round(2^(-25) / 150, SG, RD);
328// > y;
329// 0x1.b4e81ap-33
330// > printsingle(y);
331// 0x2f5a740d
332constexpr uint32_t Y_LOWER_BOUND = 0x2f5a'740d;
333
334// Check if x is an odd integer: the lowest set bit must be at the unit
335// position:
336// x_e + lsb == UNIT_EXPONENT.
337LIBC_INLINE bool is_odd_integer(float x) {
338 using FPBits = fputil::FPBits<float>;
339 FPBits xbits(x);
340 uint32_t x_u = xbits.uintval();
341 unsigned x_e = static_cast<unsigned>(xbits.get_biased_exponent());
342 unsigned lsb =
343 static_cast<unsigned>(cpp::countr_zero(value: x_u | FPBits::EXP_MASK));
344 constexpr unsigned UNIT_EXPONENT =
345 static_cast<unsigned>(FPBits::EXP_BIAS + FPBits::FRACTION_LEN);
346 return (x_e + lsb == UNIT_EXPONENT);
347}
348
349// Check if x is an integer: the lowest set bit must be at or above the unit
350// position:
351// x_e + lsb >= UNIT_EXPONENT.
352LIBC_INLINE bool is_integer(float x) {
353 if (x == 0.0f)
354 return true;
355 using FPBits = fputil::FPBits<float>;
356 FPBits xbits(x);
357 uint32_t x_u = xbits.uintval();
358 unsigned x_e = static_cast<unsigned>(xbits.get_biased_exponent());
359 unsigned lsb =
360 static_cast<unsigned>(cpp::countr_zero(value: x_u | FPBits::EXP_MASK));
361 constexpr unsigned UNIT_EXPONENT =
362 static_cast<unsigned>(FPBits::EXP_BIAS + FPBits::FRACTION_LEN);
363 return (x_e + lsb >= UNIT_EXPONENT);
364}
365
366LIBC_INLINE float set_overflow(Sign sign = Sign::POS) {
367 fputil::set_errno_if_required(ERANGE);
368 fputil::raise_overflow_except_if_required<float>();
369 using FPBits = fputil::FPBits<float>;
370#ifndef LIBC_MATH_HAS_ASSUME_ROUND_NEAREST_ONLY
371 int rounding = fputil::quick_get_round();
372 if (rounding == FE_TOWARDZERO)
373 return FPBits::max_normal(sign).get_val();
374 if (rounding == FE_DOWNWARD)
375 return sign.is_neg() ? FPBits::inf(sign: Sign::NEG).get_val()
376 : FPBits::max_normal(sign: Sign::POS).get_val();
377 if (rounding == FE_UPWARD)
378 return sign.is_neg() ? FPBits::max_normal(sign: Sign::NEG).get_val()
379 : FPBits::inf(sign: Sign::POS).get_val();
380#endif // LIBC_MATH_HAS_ASSUME_ROUND_NEAREST_ONLY
381 return FPBits::inf(sign).get_val();
382}
383
384LIBC_INLINE float set_underflow(Sign sign = Sign::POS) {
385 fputil::set_errno_if_required(ERANGE);
386 fputil::raise_underflow_except_if_required<float>();
387 using FPBits = fputil::FPBits<float>;
388#ifndef LIBC_MATH_HAS_ASSUME_ROUND_NEAREST_ONLY
389 int rounding = fputil::quick_get_round();
390 if (rounding == FE_UPWARD && sign.is_pos())
391 return FPBits::min_subnormal(sign: Sign::POS).get_val();
392 if (rounding == FE_DOWNWARD && sign.is_neg())
393 return FPBits::min_subnormal(sign: Sign::NEG).get_val();
394#endif // LIBC_MATH_HAS_ASSUME_ROUND_NEAREST_ONLY
395 return FPBits::zero(sign).get_val();
396}
397
398// Fast checks for special inputs:
399// x = 0, +-1, 2^k, 10, +- Inf
400// y = 0, +-1, 2, 0.5
401LIBC_ALWAYS_INLINE cpp::optional<float> check_special_inputs(float x, float y) {
402 using FPBits = fputil::FPBits<float>;
403 FPBits xbits(x), ybits(y);
404
405 bool x_sign = xbits.sign() == Sign::NEG;
406 bool y_sign = ybits.sign() == Sign::NEG;
407
408 FPBits x_abs = xbits.abs();
409 FPBits y_abs = ybits.abs();
410
411 // If x or y is signaling NaN
412 if (x_abs.is_signaling_nan() || y_abs.is_signaling_nan()) {
413 fputil::raise_except_if_required(FE_INVALID);
414 return FPBits::quiet_nan().get_val();
415 }
416
417 if (x == 1.0f || y == 0.0f)
418 return 1.0f;
419
420 if (x == 0.0f) {
421 if (y_abs.is_nan())
422 return y;
423 if (y_abs.is_inf())
424 return y_sign ? FPBits::inf().get_val() : 0.0f;
425 bool out_is_neg = x_sign && is_odd_integer(x: y);
426 if (y_sign) {
427 // pow(0, negative number) = Inf
428 fputil::set_errno_if_required(EDOM);
429 fputil::raise_except_if_required(FE_DIVBYZERO);
430 return FPBits::inf(sign: out_is_neg ? Sign::NEG : Sign::POS).get_val();
431 }
432 // pow(0, positive number) = 0
433 return out_is_neg ? -0.0f : 0.0f;
434 }
435
436 if (y == 1.0f)
437 return x;
438
439 if (y == 2.0f)
440 return x * x;
441
442#ifdef LIBC_TARGET_CPU_HAS_FPU_FLOAT
443 if (y == 0.5f && !x_sign)
444 return fputil::sqrt<float>(x);
445#endif // LIBC_TARGET_CPU_HAS_FPU_FLOAT
446
447 // Exact power of 2:
448 // pow(2^k, y) = 2^(k * y).
449 // This is exact if k * y is an integer.
450 uint32_t y_a = y_abs.uintval();
451 if (x_abs.is_normal() && xbits.get_mantissa() == 0 && y_a > Y_LOWER_BOUND &&
452 y_a < Y_UPPER_BOUND) {
453 if (!x_sign || is_integer(x: y)) {
454 Sign out_sign = (x_sign && is_odd_integer(x: y)) ? Sign::NEG : Sign::POS;
455 int e_x = xbits.get_exponent();
456 double ey = static_cast<double>(e_x) * static_cast<double>(y);
457 double ey_int = fputil::nearest_integer(x: ey);
458 if (ey_int == ey) {
459 if (ey > 127.0)
460 return set_overflow(out_sign);
461 if (ey < -149.0)
462 return set_underflow(out_sign);
463 int k = static_cast<int>(ey);
464 if (k >= -126)
465 return FPBits::create_value(sign: out_sign, biased_exp: static_cast<uint32_t>(k + 127),
466 mantissa: 0)
467 .get_val();
468 return FPBits::create_value(sign: out_sign, biased_exp: 0, mantissa: 1U << (k + 149)).get_val();
469 }
470 }
471 }
472
473#ifndef LIBC_MATH_HAS_SMALL_TABLES
474 if (x == 2.0f)
475 return math::exp2f(x: y);
476
477 if (x == 10.0f)
478 return math::exp10f(x: y);
479
480 // For x = 2^(+- 2^n):
481 // x^y = 2^(+- 2^n * y).
482 if (!x_sign && x_abs.is_normal() && xbits.get_mantissa() == 0 &&
483 y_a > Y_LOWER_BOUND && y_a < Y_UPPER_BOUND) {
484 int e_x = xbits.get_exponent();
485 uint32_t abs_e = static_cast<uint32_t>(e_x > 0 ? e_x : -e_x);
486 if (cpp::has_single_bit(value: abs_e)) {
487 float hi = static_cast<float>(e_x) * y;
488 if (hi >= 128.0f)
489 return set_overflow();
490 if (hi < -150.0f)
491 return set_underflow();
492 if (hi >= -126.0f)
493 return math::exp2f(x: hi);
494 }
495 }
496#endif // !LIBC_MATH_HAS_SMALL_TABLES
497
498 return cpp::nullopt;
499}
500
501// Filters out extreme input ranges, infinities, NaNs, normalizes denormal
502// inputs, and handles negative bases. Returns a value if the result is
503// determined, or nullopt if regular evaluation should proceed (in which case x,
504// y, ex, sign may be updated).
505template <typename SignType = uint64_t>
506LIBC_ALWAYS_INLINE cpp::optional<float>
507check_exceptional_cases(float &x, float &y, int &ex, SignType &sign) {
508 using FloatBits = fputil::FPBits<float>;
509 FloatBits xbits(x), ybits(y);
510 bool x_sign = xbits.sign() == Sign::NEG;
511 bool y_sign = ybits.sign() == Sign::NEG;
512
513 FloatBits x_abs = xbits.abs();
514 FloatBits y_abs = ybits.abs();
515
516 uint32_t x_a = x_abs.uintval();
517 uint32_t y_a = y_abs.uintval();
518
519 // If x or y is signaling NaN
520 if (x_abs.is_signaling_nan() || y_abs.is_signaling_nan()) {
521 fputil::raise_except_if_required(FE_INVALID);
522 return FloatBits::quiet_nan().get_val();
523 }
524
525 // Extreme |y| >= Y_UPPER_BOUND: |y * log2(x)| >= 150 for all x != 1.
526 if (y_a >= Y_UPPER_BOUND) {
527 if (x_abs.is_nan())
528 return x;
529 if (ybits.is_nan())
530 return y;
531
532 if (x_a == FloatBits::one().uintval())
533 return 1.0f;
534
535 bool is_overflow = (x_a < FloatBits::one().uintval()) == y_sign;
536 if (y_abs.is_inf() || x_a == FloatBits::inf().uintval())
537 return is_overflow ? FloatBits::inf().get_val() : 0.0f;
538
539 return is_overflow ? set_overflow() : set_underflow();
540 }
541
542 // y is finite and non-zero.
543
544 if (x_a == FloatBits::inf().uintval()) {
545 // pow(+-Inf, y):
546 // - y < 0: returns +-0.0 (negative if x = -Inf and y is odd integer).
547 // - y > 0: returns +-Inf (negative if x = -Inf and y is odd integer).
548 bool out_is_neg = x_sign && is_odd_integer(x: y);
549 Sign out_sign = out_is_neg ? Sign::NEG : Sign::POS;
550 return y_sign ? FloatBits::zero(sign: out_sign).get_val()
551 : FloatBits::inf(sign: out_sign).get_val();
552 }
553
554 if (x_a > FloatBits::inf().uintval()) {
555 // x is NaN.
556 return x;
557 }
558
559 if (x_a == 0) {
560 bool out_is_neg = x_sign && is_odd_integer(x: y);
561 if (y_sign) {
562 fputil::set_errno_if_required(EDOM);
563 fputil::raise_except_if_required(FE_DIVBYZERO);
564 return FloatBits::inf(sign: out_is_neg ? Sign::NEG : Sign::POS).get_val();
565 }
566 return out_is_neg ? -0.0f : 0.0f;
567 }
568
569 // Handle negative base x < 0:
570 // - If y is an integer: (-|x|)^y = (-1)^y * |x|^y.
571 // - If y is not an integer: (-|x|)^y is undefined in real numbers;
572 // raises FE_INVALID, sets errno to EDOM, and returns quiet NaN.
573 if (x_sign) {
574 if (is_integer(x: y)) {
575 x = -x;
576 if (is_odd_integer(x: y)) {
577 if constexpr (sizeof(SignType) == 8)
578 sign = 0x8000'0000'0000'0000ULL;
579 else
580 sign = 0x8000'0000U;
581 }
582 } else {
583 // pow( negative, non-integer ) = NaN
584 fputil::set_errno_if_required(EDOM);
585 fputil::raise_except_if_required(FE_INVALID);
586 return FloatBits::quiet_nan().get_val();
587 }
588 }
589
590 if (y_a <= Y_LOWER_BOUND) {
591 volatile float one = 1.0f;
592 volatile float eps = ((x_a < FloatBits::one().uintval()) == !y_sign)
593 ? -0x1.0p-50f
594 : 0x1.0p-50f;
595 return one + eps;
596 }
597
598 // Normalize denormal inputs.
599 if (x_a < FloatBits::min_normal().uintval()) {
600 int shift = cpp::countl_zero(value: x_a) - 8;
601 ex -= shift;
602 x = cpp::bit_cast<float>(from: x_a << shift);
603 }
604
605 return cpp::nullopt;
606}
607
608#ifndef LIBC_MATH_HAS_SMALL_TABLES
609// Check if x^y is an exact rounding boundary case (a 25-bit dyadic float,
610// which is either an exact 24-bit float or an exact halfway midpoint).
611//
612// Reference:
613// Lauter, C. and Lefevre, V., "Rounding Boundary Cases for Values of the
614// Exponential and Power Functions," IEEE Trans. Comput. 58(8):1063-1074.
615//
616// 1. x = 2^e: x^y = 2^(e * y) is in F_25 iff e * y is an integer.
617// 2. x = m * 2^e, y = n * 2^f with m, n odd integers:
618// x^y in F_25 implies:
619// - 0 <= y <= 15, n <= 15.
620// - f >= -3.
621// - e * y is an integer.
622// - When f < 0 (i.e. y = n / 2^(-f)): m^(2^f) is an integer, verified
623// by taking repeated square roots -f times.
624// - The resulting significand (m^(2^f))^n <= 2^25.
625//
626// Returns true and sets exact_m and exact_exp if x^y is an exact boundary.
627LIBC_INLINE bool is_exact_rounding_boundary(float x, float y, uint32_t &exact_m,
628 int &exact_exp) {
629 using FPBits = fputil::FPBits<float>;
630
631 FPBits xbits(x);
632 int x_e = 0;
633 uint32_t x_mant = xbits.get_mantissa();
634 if (LIBC_UNLIKELY(xbits.get_biased_exponent() == 0)) {
635 if (x_mant != 0) {
636 int shift = cpp::countl_zero(value: x_mant) - 8;
637 x_mant = (x_mant << shift) & FPBits::FRACTION_MASK;
638 x_e = -126 - shift;
639 }
640 } else {
641 x_e = xbits.get_exponent();
642 }
643
644 // Case 1: x = 2^e.
645 if (x_mant == 0) {
646 double e = static_cast<double>(x_e);
647 double ey = e * static_cast<double>(y);
648 if (fputil::nearest_integer(x: ey) == ey) {
649 exact_m = 1;
650 if (ey > 150.0)
651 exact_exp = 200;
652 else if (ey < -200.0)
653 exact_exp = -200;
654 else
655 exact_exp = static_cast<int>(ey);
656 return true;
657 }
658 return false;
659 }
660
661 // Case 2: x is not a power of 2.
662 if (y < 0.0f || y > 15.0f)
663 return false;
664
665 // Decompose y = n * 2^f with n odd.
666 FPBits ybits(y);
667 uint32_t y_mant = ybits.get_mantissa() | (1U << FPBits::FRACTION_LEN);
668 int y_exp = ybits.get_exponent() - static_cast<int>(FPBits::FRACTION_LEN);
669 int tz_y = cpp::countr_zero(value: y_mant);
670 uint32_t n = y_mant >> tz_y;
671 int f = y_exp + tz_y;
672
673 if (n > 15 || f < -3)
674 return false;
675
676 // Decompose x = m * 2^e with m odd.
677 uint32_t full_x_mant = x_mant | (1U << FPBits::FRACTION_LEN);
678 int tz_x = cpp::countr_zero(value: full_x_mant);
679 uint32_t m = full_x_mant >> tz_x;
680 int e = x_e - static_cast<int>(FPBits::FRACTION_LEN) + tz_x;
681
682 if (f < 0) {
683 // Non-integer exponent: y = n / 2^(-f) with -f in {1, 2, 3}.
684 double ey = static_cast<double>(e) * static_cast<double>(y);
685 if (fputil::nearest_integer(x: ey) != ey)
686 return false;
687
688 // Check if m^(2^f) is an integer by taking repeated square roots -f times.
689 int count = -f;
690 uint32_t cur = m;
691 for (int i = 0; i < count; ++i) {
692 uint32_t s =
693 static_cast<uint32_t>(fputil::sqrt<double>(x: static_cast<double>(cur)));
694 if (s * s != cur)
695 return false;
696 cur = s;
697 }
698
699 // Compute res = cur^n and check if res <= 2^25.
700 uint32_t res = 1;
701 for (uint32_t i = 0; i < n; ++i) {
702 if (res > (1U << 25) / cur)
703 return false;
704 res *= cur;
705 }
706 exact_m = res;
707 exact_exp = static_cast<int>(ey);
708 return true;
709 } else {
710 // Integer exponent: f >= 0, so y is an integer.
711 int y_int = static_cast<int>(y);
712 uint32_t res = 1;
713 for (int i = 0; i < y_int; ++i) {
714 if (res > (1U << 25) / m)
715 return false;
716 res *= m;
717 }
718 exact_m = res;
719 exact_exp = e * y_int;
720 return true;
721 }
722}
723#endif // !LIBC_MATH_HAS_SMALL_TABLES
724
725} // namespace powf_internal
726} // namespace math
727} // namespace LIBC_NAMESPACE_DECL
728
729#endif // LLVM_LIBC_SRC___SUPPORT_MATH_POWF_UTILS_H
730