1 /* origin: FreeBSD /usr/src/lib/msun/ld80/s_exp2l.c and /usr/src/lib/msun/ld128/s_exp2l.c */
3 * Copyright (c) 2005-2008 David Schultz <das@FreeBSD.ORG>
6 * Redistribution and use in source and binary forms, with or without
7 * modification, are permitted provided that the following conditions
9 * 1. Redistributions of source code must retain the above copyright
10 * notice, this list of conditions and the following disclaimer.
11 * 2. Redistributions in binary form must reproduce the above copyright
12 * notice, this list of conditions and the following disclaimer in the
13 * documentation and/or other materials provided with the distribution.
15 * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND
16 * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
17 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
18 * ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
19 * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
20 * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
21 * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
22 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
23 * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
24 * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
30 #if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
31 long double exp2l(long double x
)
35 #elif LDBL_MANT_DIG == 64 && LDBL_MAX_EXP == 16384
37 #define TBLSIZE (1 << TBLBITS)
40 redux
= 0x1.8p63
/ TBLSIZE
,
41 P1
= 0x1.62e42fefa39efp
-1,
42 P2
= 0x1.ebfbdff82c58fp
-3,
43 P3
= 0x1.c6b08d7049fap
-5,
44 P4
= 0x1.3b2ab6fba4da5p
-7,
45 P5
= 0x1.5d8804780a736p
-10,
46 P6
= 0x1.430918835e33dp
-13;
48 static const double tbl
[TBLSIZE
* 2] = {
49 0x1.6a09e667f3bcdp
-1, -0x1.bdd3413b2648p
-55,
50 0x1.6c012750bdabfp
-1, -0x1.2895667ff0cp
-57,
51 0x1.6dfb23c651a2fp
-1, -0x1.bbe3a683c88p
-58,
52 0x1.6ff7df9519484p
-1, -0x1.83c0f25860fp
-56,
53 0x1.71f75e8ec5f74p
-1, -0x1.16e4786887bp
-56,
54 0x1.73f9a48a58174p
-1, -0x1.0a8d96c65d5p
-55,
55 0x1.75feb564267c9p
-1, -0x1.0245957316ep
-55,
56 0x1.780694fde5d3fp
-1, 0x1.866b80a0216p
-55,
57 0x1.7a11473eb0187p
-1, -0x1.41577ee0499p
-56,
58 0x1.7c1ed0130c132p
-1, 0x1.f124cd1164ep
-55,
59 0x1.7e2f336cf4e62p
-1, 0x1.05d02ba157ap
-57,
60 0x1.80427543e1a12p
-1, -0x1.27c86626d97p
-55,
61 0x1.82589994cce13p
-1, -0x1.d4c1dd41533p
-55,
62 0x1.8471a4623c7adp
-1, -0x1.8d684a341cep
-56,
63 0x1.868d99b4492edp
-1, -0x1.fc6f89bd4f68p
-55,
64 0x1.88ac7d98a6699p
-1, 0x1.994c2f37cb5p
-55,
65 0x1.8ace5422aa0dbp
-1, 0x1.6e9f156864bp
-55,
66 0x1.8cf3216b5448cp
-1, -0x1.0d55e32e9e4p
-57,
67 0x1.8f1ae99157736p
-1, 0x1.5cc13a2e397p
-56,
68 0x1.9145b0b91ffc6p
-1, -0x1.dd6792e5825p
-55,
69 0x1.93737b0cdc5e5p
-1, -0x1.75fc781b58p
-58,
70 0x1.95a44cbc8520fp
-1, -0x1.64b7c96a5fp
-57,
71 0x1.97d829fde4e5p
-1, -0x1.d185b7c1b86p
-55,
72 0x1.9a0f170ca07bap
-1, -0x1.173bd91cee6p
-55,
73 0x1.9c49182a3f09p
-1, 0x1.c7c46b071f2p
-57,
74 0x1.9e86319e32323p
-1, 0x1.824ca78e64cp
-57,
75 0x1.a0c667b5de565p
-1, -0x1.359495d1cd5p
-55,
76 0x1.a309bec4a2d33p
-1, 0x1.6305c7ddc368p
-55,
77 0x1.a5503b23e255dp
-1, -0x1.d2f6edb8d42p
-55,
78 0x1.a799e1330b358p
-1, 0x1.bcb7ecac564p
-55,
79 0x1.a9e6b5579fdbfp
-1, 0x1.0fac90ef7fdp
-55,
80 0x1.ac36bbfd3f37ap
-1, -0x1.f9234cae76dp
-56,
81 0x1.ae89f995ad3adp
-1, 0x1.7a1cd345dcc8p
-55,
82 0x1.b0e07298db666p
-1, -0x1.bdef54c80e4p
-55,
83 0x1.b33a2b84f15fbp
-1, -0x1.2805e3084d8p
-58,
84 0x1.b59728de5593ap
-1, -0x1.c71dfbbba6ep
-55,
85 0x1.b7f76f2fb5e47p
-1, -0x1.5584f7e54acp
-57,
86 0x1.ba5b030a1064ap
-1, -0x1.efcd30e5429p
-55,
87 0x1.bcc1e904bc1d2p
-1, 0x1.23dd07a2d9fp
-56,
88 0x1.bf2c25bd71e09p
-1, -0x1.efdca3f6b9c8p
-55,
89 0x1.c199bdd85529cp
-1, 0x1.11065895049p
-56,
90 0x1.c40ab5fffd07ap
-1, 0x1.b4537e083c6p
-55,
91 0x1.c67f12e57d14bp
-1, 0x1.2884dff483c8p
-55,
92 0x1.c8f6d9406e7b5p
-1, 0x1.1acbc48805cp
-57,
93 0x1.cb720dcef9069p
-1, 0x1.503cbd1e94ap
-57,
94 0x1.cdf0b555dc3fap
-1, -0x1.dd83b53829dp
-56,
95 0x1.d072d4a07897cp
-1, -0x1.cbc3743797a8p
-55,
96 0x1.d2f87080d89f2p
-1, -0x1.d487b719d858p
-55,
97 0x1.d5818dcfba487p
-1, 0x1.2ed02d75b37p
-56,
98 0x1.d80e316c98398p
-1, -0x1.11ec18bedep
-55,
99 0x1.da9e603db3285p
-1, 0x1.c2300696db5p
-55,
100 0x1.dd321f301b46p
-1, 0x1.2da5778f019p
-55,
101 0x1.dfc97337b9b5fp
-1, -0x1.1a5cd4f184b8p
-55,
102 0x1.e264614f5a129p
-1, -0x1.7b627817a148p
-55,
103 0x1.e502ee78b3ff6p
-1, 0x1.39e8980a9cdp
-56,
104 0x1.e7a51fbc74c83p
-1, 0x1.2d522ca0c8ep
-55,
105 0x1.ea4afa2a490dap
-1, -0x1.e9c23179c288p
-55,
106 0x1.ecf482d8e67f1p
-1, -0x1.c93f3b411ad8p
-55,
107 0x1.efa1bee615a27p
-1, 0x1.dc7f486a4b68p
-55,
108 0x1.f252b376bba97p
-1, 0x1.3a1a5bf0d8e8p
-55,
109 0x1.f50765b6e454p
-1, 0x1.9d3e12dd8a18p
-55,
110 0x1.f7bfdad9cbe14p
-1, -0x1.dbb12d00635p
-55,
111 0x1.fa7c1819e90d8p
-1, 0x1.74853f3a593p
-56,
112 0x1.fd3c22b8f71f1p
-1, 0x1.2eb74966578p
-58,
114 0x1.0163da9fb3335p
+0, 0x1.b61299ab8cd8p
-54,
115 0x1.02c9a3e778061p
+0, -0x1.19083535b08p
-56,
116 0x1.04315e86e7f85p
+0, -0x1.0a31c1977c98p
-54,
117 0x1.059b0d3158574p
+0, 0x1.d73e2a475b4p
-55,
118 0x1.0706b29ddf6dep
+0, -0x1.c91dfe2b13cp
-55,
119 0x1.0874518759bc8p
+0, 0x1.186be4bb284p
-57,
120 0x1.09e3ecac6f383p
+0, 0x1.14878183161p
-54,
121 0x1.0b5586cf9890fp
+0, 0x1.8a62e4adc61p
-54,
122 0x1.0cc922b7247f7p
+0, 0x1.01edc16e24f8p
-54,
123 0x1.0e3ec32d3d1a2p
+0, 0x1.03a1727c58p
-59,
124 0x1.0fb66affed31bp
+0, -0x1.b9bedc44ebcp
-57,
125 0x1.11301d0125b51p
+0, -0x1.6c51039449bp
-54,
126 0x1.12abdc06c31ccp
+0, -0x1.1b514b36ca8p
-58,
127 0x1.1429aaea92dep
+0, -0x1.32fbf9af1368p
-54,
128 0x1.15a98c8a58e51p
+0, 0x1.2406ab9eeabp
-55,
129 0x1.172b83c7d517bp
+0, -0x1.19041b9d78ap
-55,
130 0x1.18af9388c8deap
+0, -0x1.11023d1970f8p
-54,
131 0x1.1a35beb6fcb75p
+0, 0x1.e5b4c7b4969p
-55,
132 0x1.1bbe084045cd4p
+0, -0x1.95386352ef6p
-54,
133 0x1.1d4873168b9aap
+0, 0x1.e016e00a264p
-54,
134 0x1.1ed5022fcd91dp
+0, -0x1.1df98027bb78p
-54,
135 0x1.2063b88628cd6p
+0, 0x1.dc775814a85p
-55,
136 0x1.21f49917ddc96p
+0, 0x1.2a97e9494a6p
-55,
137 0x1.2387a6e756238p
+0, 0x1.9b07eb6c7058p
-54,
138 0x1.251ce4fb2a63fp
+0, 0x1.ac155bef4f5p
-55,
139 0x1.26b4565e27cddp
+0, 0x1.2bd339940eap
-55,
140 0x1.284dfe1f56381p
+0, -0x1.a4c3a8c3f0d8p
-54,
141 0x1.29e9df51fdee1p
+0, 0x1.612e8afad12p
-55,
142 0x1.2b87fd0dad99p
+0, -0x1.10adcd6382p
-59,
143 0x1.2d285a6e4030bp
+0, 0x1.0024754db42p
-54,
144 0x1.2ecafa93e2f56p
+0, 0x1.1ca0f45d524p
-56,
145 0x1.306fe0a31b715p
+0, 0x1.6f46ad23183p
-55,
146 0x1.32170fc4cd831p
+0, 0x1.a9ce78e1804p
-55,
147 0x1.33c08b26416ffp
+0, 0x1.327218436598p
-54,
148 0x1.356c55f929ff1p
+0, -0x1.b5cee5c4e46p
-55,
149 0x1.371a7373aa9cbp
+0, -0x1.63aeabf42ebp
-54,
150 0x1.38cae6d05d866p
+0, -0x1.e958d3c99048p
-54,
151 0x1.3a7db34e59ff7p
+0, -0x1.5e436d661f6p
-56,
152 0x1.3c32dc313a8e5p
+0, -0x1.efff8375d2ap
-54,
153 0x1.3dea64c123422p
+0, 0x1.ada0911f09fp
-55,
154 0x1.3fa4504ac801cp
+0, -0x1.7d023f956fap
-54,
155 0x1.4160a21f72e2ap
+0, -0x1.ef3691c309p
-58,
156 0x1.431f5d950a897p
+0, -0x1.1c7dde35f7ap
-55,
157 0x1.44e086061892dp
+0, 0x1.89b7a04ef8p
-59,
158 0x1.46a41ed1d0057p
+0, 0x1.c944bd1648a8p
-54,
159 0x1.486a2b5c13cdp
+0, 0x1.3c1a3b69062p
-56,
160 0x1.4a32af0d7d3dep
+0, 0x1.9cb62f3d1be8p
-54,
161 0x1.4bfdad5362a27p
+0, 0x1.d4397afec42p
-56,
162 0x1.4dcb299fddd0dp
+0, 0x1.8ecdbbc6a78p
-54,
163 0x1.4f9b2769d2ca7p
+0, -0x1.4b309d25958p
-54,
164 0x1.516daa2cf6642p
+0, -0x1.f768569bd94p
-55,
165 0x1.5342b569d4f82p
+0, -0x1.07abe1db13dp
-55,
166 0x1.551a4ca5d920fp
+0, -0x1.d689cefede6p
-55,
167 0x1.56f4736b527dap
+0, 0x1.9bb2c011d938p
-54,
168 0x1.58d12d497c7fdp
+0, 0x1.295e15b9a1ep
-55,
169 0x1.5ab07dd485429p
+0, 0x1.6324c0546478p
-54,
170 0x1.5c9268a5946b7p
+0, 0x1.c4b1b81698p
-60,
171 0x1.5e76f15ad2148p
+0, 0x1.ba6f93080e68p
-54,
172 0x1.605e1b976dc09p
+0, -0x1.3e2429b56de8p
-54,
173 0x1.6247eb03a5585p
+0, -0x1.383c17e40b48p
-54,
174 0x1.6434634ccc32p
+0, -0x1.c483c759d89p
-55,
175 0x1.6623882552225p
+0, -0x1.bb60987591cp
-54,
176 0x1.68155d44ca973p
+0, 0x1.038ae44f74p
-57,
180 * exp2l(x): compute the base 2 exponential of x
182 * Accuracy: Peak error < 0.511 ulp.
184 * Method: (equally-spaced tables)
187 * x = 2**k + y, for integer k and |y| <= 1/2.
188 * Thus we have exp2l(x) = 2**k * exp2(y).
191 * y = i/TBLSIZE + z for integer i near y * TBLSIZE.
192 * Thus we have exp2(y) = exp2(i/TBLSIZE) * exp2(z),
193 * with |z| <= 2**-(TBLBITS+1).
195 * We compute exp2(i/TBLSIZE) via table lookup and exp2(z) via a
196 * degree-6 minimax polynomial with maximum error under 2**-69.
197 * The table entries each have 104 bits of accuracy, encoded as
198 * a pair of double precision values.
200 long double exp2l(long double x
)
202 union ldshape u
= {x
};
203 int e
= u
.i
.se
& 0x7fff;
206 union {uint32_t u
; int32_t i
;} k
;
208 /* Filter out exceptional cases. */
209 if (e
>= 0x3fff + 13) { /* |x| >= 8192 or x is NaN */
210 if (u
.i
.se
>= 0x3fff + 14 && u
.i
.se
>> 15 == 0)
212 return x
* 0x1p
16383L;
213 if (e
== 0x7fff) /* -inf or -nan */
216 if (x
<= -16446 || x
- 0x1p
63 + 0x1p
63 != x
)
218 FORCE_EVAL((float)(-0x1p
-149/x
));
222 } else if (e
< 0x3fff - 64) {
227 * Reduce x, computing z, i0, and k. The low bits of x + redux
228 * contain the 16-bit integer part of the exponent (k) followed by
229 * TBLBITS fractional bits (i0). We use bit tricks to extract these
230 * as integers, then set z to the remainder.
232 * Example: Suppose x is 0xabc.123456p0 and TBLBITS is 8.
233 * Then the low-order word of x + redux is 0x000abc12,
234 * We split this into k = 0xabc and i0 = 0x12 (adjusted to
235 * index into the table), then we compute z = 0x0.003456p0.
238 i0
= u
.i
.m
+ TBLSIZE
/ 2;
239 k
.u
= i0
/ TBLSIZE
* TBLSIZE
;
245 /* Compute r = exp2l(y) = exp2lt[i0] * p(z). */
246 long double t_hi
= tbl
[2*i0
];
247 long double t_lo
= tbl
[2*i0
+ 1];
248 /* XXX This gives > 1 ulp errors outside of FE_TONEAREST mode */
249 r
= t_lo
+ (t_hi
+ t_lo
) * z
* (P1
+ z
* (P2
+ z
* (P3
+ z
* (P4
250 + z
* (P5
+ z
* P6
))))) + t_hi
;
252 return scalbnl(r
, k
.i
);
254 #elif LDBL_MANT_DIG == 113 && LDBL_MAX_EXP == 16384
256 #define TBLSIZE (1 << TBLBITS)
258 static const long double
259 P1
= 0x1.62e42fefa39ef35793c7673007e6p
-1L,
260 P2
= 0x1.ebfbdff82c58ea86f16b06ec9736p
-3L,
261 P3
= 0x1.c6b08d704a0bf8b33a762bad3459p
-5L,
262 P4
= 0x1.3b2ab6fba4e7729ccbbe0b4f3fc2p
-7L,
263 P5
= 0x1.5d87fe78a67311071dee13fd11d9p
-10L,
264 P6
= 0x1.430912f86c7876f4b663b23c5fe5p
-13L;
267 P7
= 0x1.ffcbfc588b041p
-17,
268 P8
= 0x1.62c0223a5c7c7p
-20,
269 P9
= 0x1.b52541ff59713p
-24,
270 P10
= 0x1.e4cf56a391e22p
-28,
271 redux
= 0x1.8p112
/ TBLSIZE
;
273 static const long double tbl
[TBLSIZE
] = {
274 0x1.6a09e667f3bcc908b2fb1366dfeap
-1L,
275 0x1.6c012750bdabeed76a99800f4edep
-1L,
276 0x1.6dfb23c651a2ef220e2cbe1bc0d4p
-1L,
277 0x1.6ff7df9519483cf87e1b4f3e1e98p
-1L,
278 0x1.71f75e8ec5f73dd2370f2ef0b148p
-1L,
279 0x1.73f9a48a58173bd5c9a4e68ab074p
-1L,
280 0x1.75feb564267c8bf6e9aa33a489a8p
-1L,
281 0x1.780694fde5d3f619ae02808592a4p
-1L,
282 0x1.7a11473eb0186d7d51023f6ccb1ap
-1L,
283 0x1.7c1ed0130c1327c49334459378dep
-1L,
284 0x1.7e2f336cf4e62105d02ba1579756p
-1L,
285 0x1.80427543e1a11b60de67649a3842p
-1L,
286 0x1.82589994cce128acf88afab34928p
-1L,
287 0x1.8471a4623c7acce52f6b97c6444cp
-1L,
288 0x1.868d99b4492ec80e41d90ac2556ap
-1L,
289 0x1.88ac7d98a669966530bcdf2d4cc0p
-1L,
290 0x1.8ace5422aa0db5ba7c55a192c648p
-1L,
291 0x1.8cf3216b5448bef2aa1cd161c57ap
-1L,
292 0x1.8f1ae991577362b982745c72eddap
-1L,
293 0x1.9145b0b91ffc588a61b469f6b6a0p
-1L,
294 0x1.93737b0cdc5e4f4501c3f2540ae8p
-1L,
295 0x1.95a44cbc8520ee9b483695a0e7fep
-1L,
296 0x1.97d829fde4e4f8b9e920f91e8eb6p
-1L,
297 0x1.9a0f170ca07b9ba3109b8c467844p
-1L,
298 0x1.9c49182a3f0901c7c46b071f28dep
-1L,
299 0x1.9e86319e323231824ca78e64c462p
-1L,
300 0x1.a0c667b5de564b29ada8b8cabbacp
-1L,
301 0x1.a309bec4a2d3358c171f770db1f4p
-1L,
302 0x1.a5503b23e255c8b424491caf88ccp
-1L,
303 0x1.a799e1330b3586f2dfb2b158f31ep
-1L,
304 0x1.a9e6b5579fdbf43eb243bdff53a2p
-1L,
305 0x1.ac36bbfd3f379c0db966a3126988p
-1L,
306 0x1.ae89f995ad3ad5e8734d17731c80p
-1L,
307 0x1.b0e07298db66590842acdfc6fb4ep
-1L,
308 0x1.b33a2b84f15faf6bfd0e7bd941b0p
-1L,
309 0x1.b59728de559398e3881111648738p
-1L,
310 0x1.b7f76f2fb5e46eaa7b081ab53ff6p
-1L,
311 0x1.ba5b030a10649840cb3c6af5b74cp
-1L,
312 0x1.bcc1e904bc1d2247ba0f45b3d06cp
-1L,
313 0x1.bf2c25bd71e088408d7025190cd0p
-1L,
314 0x1.c199bdd85529c2220cb12a0916bap
-1L,
315 0x1.c40ab5fffd07a6d14df820f17deap
-1L,
316 0x1.c67f12e57d14b4a2137fd20f2a26p
-1L,
317 0x1.c8f6d9406e7b511acbc48805c3f6p
-1L,
318 0x1.cb720dcef90691503cbd1e949d0ap
-1L,
319 0x1.cdf0b555dc3f9c44f8958fac4f12p
-1L,
320 0x1.d072d4a07897b8d0f22f21a13792p
-1L,
321 0x1.d2f87080d89f18ade123989ea50ep
-1L,
322 0x1.d5818dcfba48725da05aeb66dff8p
-1L,
323 0x1.d80e316c98397bb84f9d048807a0p
-1L,
324 0x1.da9e603db3285708c01a5b6d480cp
-1L,
325 0x1.dd321f301b4604b695de3c0630c0p
-1L,
326 0x1.dfc97337b9b5eb968cac39ed284cp
-1L,
327 0x1.e264614f5a128a12761fa17adc74p
-1L,
328 0x1.e502ee78b3ff6273d130153992d0p
-1L,
329 0x1.e7a51fbc74c834b548b2832378a4p
-1L,
330 0x1.ea4afa2a490d9858f73a18f5dab4p
-1L,
331 0x1.ecf482d8e67f08db0312fb949d50p
-1L,
332 0x1.efa1bee615a27771fd21a92dabb6p
-1L,
333 0x1.f252b376bba974e8696fc3638f24p
-1L,
334 0x1.f50765b6e4540674f84b762861a6p
-1L,
335 0x1.f7bfdad9cbe138913b4bfe72bd78p
-1L,
336 0x1.fa7c1819e90d82e90a7e74b26360p
-1L,
337 0x1.fd3c22b8f71f10975ba4b32bd006p
-1L,
338 0x1.0000000000000000000000000000p
+0L,
339 0x1.0163da9fb33356d84a66ae336e98p
+0L,
340 0x1.02c9a3e778060ee6f7caca4f7a18p
+0L,
341 0x1.04315e86e7f84bd738f9a20da442p
+0L,
342 0x1.059b0d31585743ae7c548eb68c6ap
+0L,
343 0x1.0706b29ddf6ddc6dc403a9d87b1ep
+0L,
344 0x1.0874518759bc808c35f25d942856p
+0L,
345 0x1.09e3ecac6f3834521e060c584d5cp
+0L,
346 0x1.0b5586cf9890f6298b92b7184200p
+0L,
347 0x1.0cc922b7247f7407b705b893dbdep
+0L,
348 0x1.0e3ec32d3d1a2020742e4f8af794p
+0L,
349 0x1.0fb66affed31af232091dd8a169ep
+0L,
350 0x1.11301d0125b50a4ebbf1aed9321cp
+0L,
351 0x1.12abdc06c31cbfb92bad324d6f84p
+0L,
352 0x1.1429aaea92ddfb34101943b2588ep
+0L,
353 0x1.15a98c8a58e512480d573dd562aep
+0L,
354 0x1.172b83c7d517adcdf7c8c50eb162p
+0L,
355 0x1.18af9388c8de9bbbf70b9a3c269cp
+0L,
356 0x1.1a35beb6fcb753cb698f692d2038p
+0L,
357 0x1.1bbe084045cd39ab1e72b442810ep
+0L,
358 0x1.1d4873168b9aa7805b8028990be8p
+0L,
359 0x1.1ed5022fcd91cb8819ff61121fbep
+0L,
360 0x1.2063b88628cd63b8eeb0295093f6p
+0L,
361 0x1.21f49917ddc962552fd29294bc20p
+0L,
362 0x1.2387a6e75623866c1fadb1c159c0p
+0L,
363 0x1.251ce4fb2a63f3582ab7de9e9562p
+0L,
364 0x1.26b4565e27cdd257a673281d3068p
+0L,
365 0x1.284dfe1f5638096cf15cf03c9fa0p
+0L,
366 0x1.29e9df51fdee12c25d15f5a25022p
+0L,
367 0x1.2b87fd0dad98ffddea46538fca24p
+0L,
368 0x1.2d285a6e4030b40091d536d0733ep
+0L,
369 0x1.2ecafa93e2f5611ca0f45d5239a4p
+0L,
370 0x1.306fe0a31b7152de8d5a463063bep
+0L,
371 0x1.32170fc4cd8313539cf1c3009330p
+0L,
372 0x1.33c08b26416ff4c9c8610d96680ep
+0L,
373 0x1.356c55f929ff0c94623476373be4p
+0L,
374 0x1.371a7373aa9caa7145502f45452ap
+0L,
375 0x1.38cae6d05d86585a9cb0d9bed530p
+0L,
376 0x1.3a7db34e59ff6ea1bc9299e0a1fep
+0L,
377 0x1.3c32dc313a8e484001f228b58cf0p
+0L,
378 0x1.3dea64c12342235b41223e13d7eep
+0L,
379 0x1.3fa4504ac801ba0bf701aa417b9cp
+0L,
380 0x1.4160a21f72e29f84325b8f3dbacap
+0L,
381 0x1.431f5d950a896dc704439410b628p
+0L,
382 0x1.44e086061892d03136f409df0724p
+0L,
383 0x1.46a41ed1d005772512f459229f0ap
+0L,
384 0x1.486a2b5c13cd013c1a3b69062f26p
+0L,
385 0x1.4a32af0d7d3de672d8bcf46f99b4p
+0L,
386 0x1.4bfdad5362a271d4397afec42e36p
+0L,
387 0x1.4dcb299fddd0d63b36ef1a9e19dep
+0L,
388 0x1.4f9b2769d2ca6ad33d8b69aa0b8cp
+0L,
389 0x1.516daa2cf6641c112f52c84d6066p
+0L,
390 0x1.5342b569d4f81df0a83c49d86bf4p
+0L,
391 0x1.551a4ca5d920ec52ec620243540cp
+0L,
392 0x1.56f4736b527da66ecb004764e61ep
+0L,
393 0x1.58d12d497c7fd252bc2b7343d554p
+0L,
394 0x1.5ab07dd48542958c93015191e9a8p
+0L,
395 0x1.5c9268a5946b701c4b1b81697ed4p
+0L,
396 0x1.5e76f15ad21486e9be4c20399d12p
+0L,
397 0x1.605e1b976dc08b076f592a487066p
+0L,
398 0x1.6247eb03a5584b1f0fa06fd2d9eap
+0L,
399 0x1.6434634ccc31fc76f8714c4ee122p
+0L,
400 0x1.66238825522249127d9e29b92ea2p
+0L,
401 0x1.68155d44ca973081c57227b9f69ep
+0L,
404 static const float eps
[TBLSIZE
] = {
536 * exp2l(x): compute the base 2 exponential of x
538 * Accuracy: Peak error < 0.502 ulp.
540 * Method: (accurate tables)
543 * x = 2**k + y, for integer k and |y| <= 1/2.
544 * Thus we have exp2(x) = 2**k * exp2(y).
547 * y = i/TBLSIZE + z - eps[i] for integer i near y * TBLSIZE.
548 * Thus we have exp2(y) = exp2(i/TBLSIZE) * exp2(z - eps[i]),
549 * with |z - eps[i]| <= 2**-8 + 2**-98 for the table used.
551 * We compute exp2(i/TBLSIZE) via table lookup and exp2(z - eps[i]) via
552 * a degree-10 minimax polynomial with maximum error under 2**-120.
553 * The values in exp2t[] and eps[] are chosen such that
554 * exp2t[i] = exp2(i/TBLSIZE + eps[i]), and eps[i] is a small offset such
555 * that exp2t[i] is accurate to 2**-122.
557 * Note that the range of i is +-TBLSIZE/2, so we actually index the tables
558 * by i0 = i + TBLSIZE/2.
560 * This method is due to Gal, with many details due to Gal and Bachelis:
562 * Gal, S. and Bachelis, B. An Accurate Elementary Mathematical Library
563 * for the IEEE Floating Point Standard. TOMS 17(1), 26-46 (1991).
568 union ldshape u
= {x
};
569 int e
= u
.i
.se
& 0x7fff;
572 union {uint32_t u
; int32_t i
;} k
;
574 /* Filter out exceptional cases. */
575 if (e
>= 0x3fff + 14) { /* |x| >= 16384 or x is NaN */
576 if (u
.i
.se
>= 0x3fff + 15 && u
.i
.se
>> 15 == 0)
578 return x
* 0x1p
16383L;
579 if (e
== 0x7fff) /* -inf or -nan */
582 if (x
<= -16495 || x
- 0x1p
112 + 0x1p
112 != x
)
584 FORCE_EVAL((float)(-0x1p
-149/x
));
588 } else if (e
< 0x3fff - 114) {
593 * Reduce x, computing z, i0, and k. The low bits of x + redux
594 * contain the 16-bit integer part of the exponent (k) followed by
595 * TBLBITS fractional bits (i0). We use bit tricks to extract these
596 * as integers, then set z to the remainder.
598 * Example: Suppose x is 0xabc.123456p0 and TBLBITS is 8.
599 * Then the low-order word of x + redux is 0x000abc12,
600 * We split this into k = 0xabc and i0 = 0x12 (adjusted to
601 * index into the table), then we compute z = 0x0.003456p0.
604 i0
= u
.i2
.lo
+ TBLSIZE
/ 2;
605 k
.u
= i0
/ TBLSIZE
* TBLSIZE
;
611 /* Compute r = exp2(y) = exp2t[i0] * p(z - eps[i]). */
614 r
= t
+ t
* z
* (P1
+ z
* (P2
+ z
* (P3
+ z
* (P4
+ z
* (P5
+ z
* (P6
615 + z
* (P7
+ z
* (P8
+ z
* (P9
+ z
* P10
)))))))));
617 return scalbnl(r
, k
.i
);