1 /*
2 * Double-precision SVE pow(x, y) function.
3 *
4 * Copyright (c) 2022-2025, Arm Limited.
5 * SPDX-License-Identifier: MIT OR Apache-2.0 WITH LLVM-exception
6 */
7
8 #include "sv_math.h"
9 #include "test_sig.h"
10 #include "test_defs.h"
11
12 /* This version share a similar algorithm as AOR scalar pow.
13
14 The core computation consists in computing pow(x, y) as
15
16 exp (y * log (x)).
17
18 The algorithms for exp and log are very similar to scalar exp and log.
19 The log relies on table lookup for 3 variables and an order 8 polynomial.
20 It returns a high and a low contribution that are then passed to the exp,
21 to minimise the loss of accuracy in both routines.
22 The exp is based on 8-bit table lookup for scale and order-4 polynomial.
23 The SVE algorithm drops the tail in the exp computation at the price of
24 a lower accuracy, slightly above 1ULP.
25 The SVE algorithm also drops the special treatement of small (< 2^-65) and
26 large (> 2^63) finite values of |y|, as they only affect non-round to
27 nearest modes.
28
29 Maximum measured error is 1.04 ULPs:
30 SV_NAME_D2 (pow) (0x1.3d2d45bc848acp+63, -0x1.a48a38b40cd43p-12)
31 got 0x1.f7116284221fcp-1
32 want 0x1.f7116284221fdp-1. */
33
34 /* Data is defined in v_pow_log_data.c. */
35 #define N_LOG (1 << V_POW_LOG_TABLE_BITS)
36 #define Off 0x3fe6955500000000
37
38 /* Data is defined in v_pow_exp_data.c. */
39 #define N_EXP (1 << V_POW_EXP_TABLE_BITS)
40 #define SignBias (0x800 << V_POW_EXP_TABLE_BITS)
41 #define SmallExp 0x3c9 /* top12(0x1p-54). */
42 #define BigExp 0x408 /* top12(512.). */
43 #define ThresExp 0x03f /* BigExp - SmallExp. */
44 #define HugeExp 0x409 /* top12(1024.). */
45
46 /* Constants associated with pow. */
47 #define SmallBoundX 0x1p-126
48 #define SmallPowX 0x001 /* top12(0x1p-126). */
49 #define BigPowX 0x7ff /* top12(INFINITY). */
50 #define ThresPowX 0x7fe /* BigPowX - SmallPowX. */
51 #define SmallPowY 0x3be /* top12(0x1.e7b6p-65). */
52 #define BigPowY 0x43e /* top12(0x1.749p62). */
53 #define ThresPowY 0x080 /* BigPowY - SmallPowY. */
54
55 static const struct data
56 {
57 double log_c0, log_c2, log_c4, log_c6, ln2_hi, ln2_lo;
58 double log_c1, log_c3, log_c5, off;
59 double n_over_ln2, exp_c2, ln2_over_n_hi, ln2_over_n_lo;
60 double exp_c0, exp_c1;
61 } data = {
62 .log_c0 = -0x1p-1,
63 .log_c1 = -0x1.555555555556p-1,
64 .log_c2 = 0x1.0000000000006p-1,
65 .log_c3 = 0x1.999999959554ep-1,
66 .log_c4 = -0x1.555555529a47ap-1,
67 .log_c5 = -0x1.2495b9b4845e9p0,
68 .log_c6 = 0x1.0002b8b263fc3p0,
69 .off = Off,
70 .exp_c0 = 0x1.fffffffffffd4p-2,
71 .exp_c1 = 0x1.5555571d6ef9p-3,
72 .exp_c2 = 0x1.5555576a5adcep-5,
73 .ln2_hi = 0x1.62e42fefa3800p-1,
74 .ln2_lo = 0x1.ef35793c76730p-45,
75 .n_over_ln2 = 0x1.71547652b82fep0 * N_EXP,
76 .ln2_over_n_hi = 0x1.62e42fefc0000p-9,
77 .ln2_over_n_lo = -0x1.c610ca86c3899p-45,
78 };
79
80 /* Check if x is an integer. */
81 static inline svbool_t
sv_isint(svbool_t pg,svfloat64_t x)82 sv_isint (svbool_t pg, svfloat64_t x)
83 {
84 return svcmpeq (pg, svrintz_z (pg, x), x);
85 }
86
87 /* Check if x is real not integer valued. */
88 static inline svbool_t
sv_isnotint(svbool_t pg,svfloat64_t x)89 sv_isnotint (svbool_t pg, svfloat64_t x)
90 {
91 return svcmpne (pg, svrintz_z (pg, x), x);
92 }
93
94 /* Check if x is an odd integer. */
95 static inline svbool_t
sv_isodd(svbool_t pg,svfloat64_t x)96 sv_isodd (svbool_t pg, svfloat64_t x)
97 {
98 svfloat64_t y = svmul_x (svptrue_b64 (), x, 0.5);
99 return sv_isnotint (pg, y);
100 }
101
102 /* Returns 0 if not int, 1 if odd int, 2 if even int. The argument is
103 the bit representation of a non-zero finite floating-point value. */
104 static inline int
checkint(uint64_t iy)105 checkint (uint64_t iy)
106 {
107 int e = iy >> 52 & 0x7ff;
108 if (e < 0x3ff)
109 return 0;
110 if (e > 0x3ff + 52)
111 return 2;
112 if (iy & ((1ULL << (0x3ff + 52 - e)) - 1))
113 return 0;
114 if (iy & (1ULL << (0x3ff + 52 - e)))
115 return 1;
116 return 2;
117 }
118
119 /* Top 12 bits (sign and exponent of each double float lane). */
120 static inline svuint64_t
sv_top12(svfloat64_t x)121 sv_top12 (svfloat64_t x)
122 {
123 return svlsr_x (svptrue_b64 (), svreinterpret_u64 (x), 52);
124 }
125
126 /* Returns 1 if input is the bit representation of 0, infinity or nan. */
127 static inline int
zeroinfnan(uint64_t i)128 zeroinfnan (uint64_t i)
129 {
130 return 2 * i - 1 >= 2 * asuint64 (INFINITY) - 1;
131 }
132
133 /* Returns 1 if input is the bit representation of 0, infinity or nan. */
134 static inline svbool_t
sv_zeroinfnan(svbool_t pg,svuint64_t i)135 sv_zeroinfnan (svbool_t pg, svuint64_t i)
136 {
137 return svcmpge (pg, svsub_x (pg, svadd_x (pg, i, i), 1),
138 2 * asuint64 (INFINITY) - 1);
139 }
140
141 /* Handle cases that may overflow or underflow when computing the result that
142 is scale*(1+TMP) without intermediate rounding. The bit representation of
143 scale is in SBITS, however it has a computed exponent that may have
144 overflown into the sign bit so that needs to be adjusted before using it as
145 a double. (int32_t)KI is the k used in the argument reduction and exponent
146 adjustment of scale, positive k here means the result may overflow and
147 negative k means the result may underflow. */
148 static inline double
specialcase(double tmp,uint64_t sbits,uint64_t ki)149 specialcase (double tmp, uint64_t sbits, uint64_t ki)
150 {
151 double scale;
152 if ((ki & 0x80000000) == 0)
153 {
154 /* k > 0, the exponent of scale might have overflowed by <= 460. */
155 sbits -= 1009ull << 52;
156 scale = asdouble (sbits);
157 return 0x1p1009 * (scale + scale * tmp);
158 }
159 /* k < 0, need special care in the subnormal range. */
160 sbits += 1022ull << 52;
161 /* Note: sbits is signed scale. */
162 scale = asdouble (sbits);
163 double y = scale + scale * tmp;
164 return 0x1p-1022 * y;
165 }
166
167 /* Scalar fallback for special cases of SVE pow's exp. */
168 static inline svfloat64_t
sv_call_specialcase(svfloat64_t x1,svuint64_t u1,svuint64_t u2,svfloat64_t y,svbool_t cmp)169 sv_call_specialcase (svfloat64_t x1, svuint64_t u1, svuint64_t u2,
170 svfloat64_t y, svbool_t cmp)
171 {
172 svbool_t p = svpfirst (cmp, svpfalse ());
173 while (svptest_any (cmp, p))
174 {
175 double sx1 = svclastb (p, 0, x1);
176 uint64_t su1 = svclastb (p, 0, u1);
177 uint64_t su2 = svclastb (p, 0, u2);
178 double elem = specialcase (sx1, su1, su2);
179 svfloat64_t y2 = sv_f64 (elem);
180 y = svsel (p, y2, y);
181 p = svpnext_b64 (cmp, p);
182 }
183 return y;
184 }
185
186 /* Compute y+TAIL = log(x) where the rounded result is y and TAIL has about
187 additional 15 bits precision. IX is the bit representation of x, but
188 normalized in the subnormal range using the sign bit for the exponent. */
189 static inline svfloat64_t
sv_log_inline(svbool_t pg,svuint64_t ix,svfloat64_t * tail,const struct data * d)190 sv_log_inline (svbool_t pg, svuint64_t ix, svfloat64_t *tail,
191 const struct data *d)
192 {
193 /* x = 2^k z; where z is in range [Off,2*Off) and exact.
194 The range is split into N subintervals.
195 The ith subinterval contains z and c is near its center. */
196 svuint64_t tmp = svsub_x (pg, ix, d->off);
197 svuint64_t i = svand_x (pg, svlsr_x (pg, tmp, 52 - V_POW_LOG_TABLE_BITS),
198 sv_u64 (N_LOG - 1));
199 svint64_t k = svasr_x (pg, svreinterpret_s64 (tmp), 52);
200 svuint64_t iz = svsub_x (pg, ix, svlsl_x (pg, svreinterpret_u64 (k), 52));
201 svfloat64_t z = svreinterpret_f64 (iz);
202 svfloat64_t kd = svcvt_f64_x (pg, k);
203
204 /* log(x) = k*Ln2 + log(c) + log1p(z/c-1). */
205 /* SVE lookup requires 3 separate lookup tables, as opposed to scalar version
206 that uses array of structures. We also do the lookup earlier in the code
207 to make sure it finishes as early as possible. */
208 svfloat64_t invc = svld1_gather_index (pg, __v_pow_log_data.invc, i);
209 svfloat64_t logc = svld1_gather_index (pg, __v_pow_log_data.logc, i);
210 svfloat64_t logctail = svld1_gather_index (pg, __v_pow_log_data.logctail, i);
211
212 /* Note: 1/c is j/N or j/N/2 where j is an integer in [N,2N) and
213 |z/c - 1| < 1/N, so r = z/c - 1 is exactly representible. */
214 svfloat64_t r = svmad_x (pg, z, invc, -1.0);
215 /* k*Ln2 + log(c) + r. */
216
217 svfloat64_t ln2_hilo = svld1rq_f64 (svptrue_b64 (), &d->ln2_hi);
218 svfloat64_t t1 = svmla_lane_f64 (logc, kd, ln2_hilo, 0);
219 svfloat64_t t2 = svadd_x (pg, t1, r);
220 svfloat64_t lo1 = svmla_lane_f64 (logctail, kd, ln2_hilo, 1);
221 svfloat64_t lo2 = svadd_x (pg, svsub_x (pg, t1, t2), r);
222
223 /* Evaluation is optimized assuming superscalar pipelined execution. */
224
225 svfloat64_t log_c02 = svld1rq_f64 (svptrue_b64 (), &d->log_c0);
226 svfloat64_t ar = svmul_lane_f64 (r, log_c02, 0);
227 svfloat64_t ar2 = svmul_x (svptrue_b64 (), r, ar);
228 svfloat64_t ar3 = svmul_x (svptrue_b64 (), r, ar2);
229 /* k*Ln2 + log(c) + r + A[0]*r*r. */
230 svfloat64_t hi = svadd_x (pg, t2, ar2);
231 svfloat64_t lo3 = svmls_x (pg, ar2, ar, r);
232 svfloat64_t lo4 = svadd_x (pg, svsub_x (pg, t2, hi), ar2);
233 /* p = log1p(r) - r - A[0]*r*r. */
234 /* p = (ar3 * (A[1] + r * A[2] + ar2 * (A[3] + r * A[4] + ar2 * (A[5] + r *
235 A[6])))). */
236
237 svfloat64_t log_c46 = svld1rq_f64 (svptrue_b64 (), &d->log_c4);
238 svfloat64_t a56 = svmla_lane_f64 (sv_f64 (d->log_c5), r, log_c46, 1);
239 svfloat64_t a34 = svmla_lane_f64 (sv_f64 (d->log_c3), r, log_c46, 0);
240 svfloat64_t a12 = svmla_lane_f64 (sv_f64 (d->log_c1), r, log_c02, 1);
241 svfloat64_t p = svmla_x (pg, a34, ar2, a56);
242 p = svmla_x (pg, a12, ar2, p);
243 p = svmul_x (svptrue_b64 (), ar3, p);
244 svfloat64_t lo = svadd_x (
245 pg, svadd_x (pg, svsub_x (pg, svadd_x (pg, lo1, lo2), lo3), lo4), p);
246 svfloat64_t y = svadd_x (pg, hi, lo);
247 *tail = svadd_x (pg, svsub_x (pg, hi, y), lo);
248 return y;
249 }
250
251 static inline svfloat64_t
sv_exp_core(svbool_t pg,svfloat64_t x,svfloat64_t xtail,svuint64_t sign_bias,svfloat64_t * tmp,svuint64_t * sbits,svuint64_t * ki,const struct data * d)252 sv_exp_core (svbool_t pg, svfloat64_t x, svfloat64_t xtail,
253 svuint64_t sign_bias, svfloat64_t *tmp, svuint64_t *sbits,
254 svuint64_t *ki, const struct data *d)
255 {
256 /* exp(x) = 2^(k/N) * exp(r), with exp(r) in [2^(-1/2N),2^(1/2N)]. */
257 /* x = ln2/N*k + r, with int k and r in [-ln2/2N, ln2/2N]. */
258 svfloat64_t n_over_ln2_and_c2 = svld1rq_f64 (svptrue_b64 (), &d->n_over_ln2);
259 svfloat64_t z = svmul_lane_f64 (x, n_over_ln2_and_c2, 0);
260 /* z - kd is in [-1, 1] in non-nearest rounding modes. */
261 svfloat64_t kd = svrinta_x (pg, z);
262 *ki = svreinterpret_u64 (svcvt_s64_x (pg, kd));
263
264 svfloat64_t ln2_over_n_hilo
265 = svld1rq_f64 (svptrue_b64 (), &d->ln2_over_n_hi);
266 svfloat64_t r = x;
267 r = svmls_lane_f64 (r, kd, ln2_over_n_hilo, 0);
268 r = svmls_lane_f64 (r, kd, ln2_over_n_hilo, 1);
269 /* The code assumes 2^-200 < |xtail| < 2^-8/N. */
270 r = svadd_x (pg, r, xtail);
271 /* 2^(k/N) ~= scale. */
272 svuint64_t idx = svand_x (pg, *ki, N_EXP - 1);
273 svuint64_t top
274 = svlsl_x (pg, svadd_x (pg, *ki, sign_bias), 52 - V_POW_EXP_TABLE_BITS);
275 /* This is only a valid scale when -1023*N < k < 1024*N. */
276 *sbits = svld1_gather_index (pg, __v_pow_exp_data.sbits, idx);
277 *sbits = svadd_x (pg, *sbits, top);
278 /* exp(x) = 2^(k/N) * exp(r) ~= scale + scale * (exp(r) - 1). */
279 svfloat64_t r2 = svmul_x (svptrue_b64 (), r, r);
280 *tmp = svmla_lane_f64 (sv_f64 (d->exp_c1), r, n_over_ln2_and_c2, 1);
281 *tmp = svmla_x (pg, sv_f64 (d->exp_c0), r, *tmp);
282 *tmp = svmla_x (pg, r, r2, *tmp);
283 svfloat64_t scale = svreinterpret_f64 (*sbits);
284 /* Note: tmp == 0 or |tmp| > 2^-200 and scale > 2^-739, so there
285 is no spurious underflow here even without fma. */
286 z = svmla_x (pg, scale, scale, *tmp);
287 return z;
288 }
289
290 /* Computes sign*exp(x+xtail) where |xtail| < 2^-8/N and |xtail| <= |x|.
291 The sign_bias argument is SignBias or 0 and sets the sign to -1 or 1. */
292 static inline svfloat64_t
sv_exp_inline(svbool_t pg,svfloat64_t x,svfloat64_t xtail,svuint64_t sign_bias,const struct data * d)293 sv_exp_inline (svbool_t pg, svfloat64_t x, svfloat64_t xtail,
294 svuint64_t sign_bias, const struct data *d)
295 {
296 /* 3 types of special cases: tiny (uflow and spurious uflow), huge (oflow)
297 and other cases of large values of x (scale * (1 + TMP) oflow). */
298 svuint64_t abstop = svand_x (pg, sv_top12 (x), 0x7ff);
299 /* |x| is large (|x| >= 512) or tiny (|x| <= 0x1p-54). */
300 svbool_t uoflow = svcmpge (pg, svsub_x (pg, abstop, SmallExp), ThresExp);
301
302 svfloat64_t tmp;
303 svuint64_t sbits, ki;
304 if (unlikely (svptest_any (pg, uoflow)))
305 {
306 svfloat64_t z
307 = sv_exp_core (pg, x, xtail, sign_bias, &tmp, &sbits, &ki, d);
308
309 /* |x| is tiny (|x| <= 0x1p-54). */
310 svbool_t uflow
311 = svcmpge (pg, svsub_x (pg, abstop, SmallExp), 0x80000000);
312 uflow = svand_z (pg, uoflow, uflow);
313 /* |x| is huge (|x| >= 1024). */
314 svbool_t oflow = svcmpge (pg, abstop, HugeExp);
315 oflow = svand_z (pg, uoflow, svbic_z (pg, oflow, uflow));
316
317 /* For large |x| values (512 < |x| < 1024) scale * (1 + TMP) can overflow
318 or underflow. */
319 svbool_t special = svbic_z (pg, uoflow, svorr_z (pg, uflow, oflow));
320
321 /* Update result with special and large cases. */
322 z = sv_call_specialcase (tmp, sbits, ki, z, special);
323
324 /* Handle underflow and overflow. */
325 svbool_t x_is_neg = svcmplt (pg, x, 0);
326 svuint64_t sign_mask
327 = svlsl_x (pg, sign_bias, 52 - V_POW_EXP_TABLE_BITS);
328 svfloat64_t res_uoflow
329 = svsel (x_is_neg, sv_f64 (0.0), sv_f64 (INFINITY));
330 res_uoflow = svreinterpret_f64 (
331 svorr_x (pg, svreinterpret_u64 (res_uoflow), sign_mask));
332 /* Avoid spurious underflow for tiny x. */
333 svfloat64_t res_spurious_uflow
334 = svreinterpret_f64 (svorr_x (pg, sign_mask, 0x3ff0000000000000));
335
336 z = svsel (oflow, res_uoflow, z);
337 z = svsel (uflow, res_spurious_uflow, z);
338 return z;
339 }
340
341 return sv_exp_core (pg, x, xtail, sign_bias, &tmp, &sbits, &ki, d);
342 }
343
344 static inline double
pow_sc(double x,double y)345 pow_sc (double x, double y)
346 {
347 uint64_t ix = asuint64 (x);
348 uint64_t iy = asuint64 (y);
349 /* Special cases: |x| or |y| is 0, inf or nan. */
350 if (unlikely (zeroinfnan (iy)))
351 {
352 if (2 * iy == 0)
353 return issignaling_inline (x) ? x + y : 1.0;
354 if (ix == asuint64 (1.0))
355 return issignaling_inline (y) ? x + y : 1.0;
356 if (2 * ix > 2 * asuint64 (INFINITY) || 2 * iy > 2 * asuint64 (INFINITY))
357 return x + y;
358 if (2 * ix == 2 * asuint64 (1.0))
359 return 1.0;
360 if ((2 * ix < 2 * asuint64 (1.0)) == !(iy >> 63))
361 return 0.0; /* |x|<1 && y==inf or |x|>1 && y==-inf. */
362 return y * y;
363 }
364 if (unlikely (zeroinfnan (ix)))
365 {
366 double_t x2 = x * x;
367 if (ix >> 63 && checkint (iy) == 1)
368 x2 = -x2;
369 return (iy >> 63) ? 1 / x2 : x2;
370 }
371 return x;
372 }
373
SV_NAME_D2(pow)374 svfloat64_t SV_NAME_D2 (pow) (svfloat64_t x, svfloat64_t y, const svbool_t pg)
375 {
376 const struct data *d = ptr_barrier (&data);
377
378 /* This preamble handles special case conditions used in the final scalar
379 fallbacks. It also updates ix and sign_bias, that are used in the core
380 computation too, i.e., exp( y * log (x) ). */
381 svuint64_t vix0 = svreinterpret_u64 (x);
382 svuint64_t viy0 = svreinterpret_u64 (y);
383
384 /* Negative x cases. */
385 svbool_t xisneg = svcmplt (pg, x, 0);
386
387 /* Set sign_bias and ix depending on sign of x and nature of y. */
388 svbool_t yint_or_xpos = pg;
389 svuint64_t sign_bias = sv_u64 (0);
390 svuint64_t vix = vix0;
391 if (unlikely (svptest_any (pg, xisneg)))
392 {
393 /* Determine nature of y. */
394 yint_or_xpos = sv_isint (xisneg, y);
395 svbool_t yisodd_xisneg = sv_isodd (xisneg, y);
396 /* ix set to abs(ix) if y is integer. */
397 vix = svand_m (yint_or_xpos, vix0, 0x7fffffffffffffff);
398 /* Set to SignBias if x is negative and y is odd. */
399 sign_bias = svsel (yisodd_xisneg, sv_u64 (SignBias), sv_u64 (0));
400 }
401
402 /* Small cases of x: |x| < 0x1p-126. */
403 svbool_t xsmall = svaclt (yint_or_xpos, x, SmallBoundX);
404 if (unlikely (svptest_any (yint_or_xpos, xsmall)))
405 {
406 /* Normalize subnormal x so exponent becomes negative. */
407 svuint64_t vtopx = svlsr_x (svptrue_b64 (), vix, 52);
408 svbool_t topx_is_null = svcmpeq (xsmall, vtopx, 0);
409
410 svuint64_t vix_norm = svreinterpret_u64 (svmul_m (xsmall, x, 0x1p52));
411 vix_norm = svand_m (xsmall, vix_norm, 0x7fffffffffffffff);
412 vix_norm = svsub_m (xsmall, vix_norm, 52ULL << 52);
413 vix = svsel (topx_is_null, vix_norm, vix);
414 }
415
416 /* y_hi = log(ix, &y_lo). */
417 svfloat64_t vlo;
418 svfloat64_t vhi = sv_log_inline (yint_or_xpos, vix, &vlo, d);
419
420 /* z = exp(y_hi, y_lo, sign_bias). */
421 svfloat64_t vehi = svmul_x (svptrue_b64 (), y, vhi);
422 svfloat64_t vemi = svmls_x (yint_or_xpos, vehi, y, vhi);
423 svfloat64_t velo = svnmls_x (yint_or_xpos, vemi, y, vlo);
424 svfloat64_t vz = sv_exp_inline (yint_or_xpos, vehi, velo, sign_bias, d);
425
426 /* Cases of finite y and finite negative x. */
427 vz = svsel (yint_or_xpos, vz, sv_f64 (__builtin_nan ("")));
428
429 /* Special cases of x or y: zero, inf and nan. */
430 svbool_t xspecial = sv_zeroinfnan (svptrue_b64 (), vix0);
431 svbool_t yspecial = sv_zeroinfnan (svptrue_b64 (), viy0);
432 svbool_t special = svorr_z (svptrue_b64 (), xspecial, yspecial);
433
434 /* Cases of zero/inf/nan x or y. */
435 if (unlikely (svptest_any (svptrue_b64 (), special)))
436 vz = sv_call2_f64 (pow_sc, x, y, vz, special);
437
438 return vz;
439 }
440
441 TEST_SIG (SV, D, 2, pow)
442 TEST_ULP (SV_NAME_D2 (pow), 0.55)
443 TEST_DISABLE_FENV (SV_NAME_D2 (pow))
444 /* Wide intervals spanning the whole domain but shared between x and y. */
445 #define SV_POW_INTERVAL2(xlo, xhi, ylo, yhi, n) \
446 TEST_INTERVAL2 (SV_NAME_D2 (pow), xlo, xhi, ylo, yhi, n) \
447 TEST_INTERVAL2 (SV_NAME_D2 (pow), xlo, xhi, -ylo, -yhi, n) \
448 TEST_INTERVAL2 (SV_NAME_D2 (pow), -xlo, -xhi, ylo, yhi, n) \
449 TEST_INTERVAL2 (SV_NAME_D2 (pow), -xlo, -xhi, -ylo, -yhi, n)
450 #define EXPAND(str) str##000000000
451 #define SHL52(str) EXPAND (str)
452 SV_POW_INTERVAL2 (0, SHL52 (SmallPowX), 0, inf, 40000)
453 SV_POW_INTERVAL2 (SHL52 (SmallPowX), SHL52 (BigPowX), 0, inf, 40000)
454 SV_POW_INTERVAL2 (SHL52 (BigPowX), inf, 0, inf, 40000)
455 SV_POW_INTERVAL2 (0, inf, 0, SHL52 (SmallPowY), 40000)
456 SV_POW_INTERVAL2 (0, inf, SHL52 (SmallPowY), SHL52 (BigPowY), 40000)
457 SV_POW_INTERVAL2 (0, inf, SHL52 (BigPowY), inf, 40000)
458 SV_POW_INTERVAL2 (0, inf, 0, inf, 1000)
459 /* x~1 or y~1. */
460 SV_POW_INTERVAL2 (0x1p-1, 0x1p1, 0x1p-10, 0x1p10, 10000)
461 SV_POW_INTERVAL2 (0x1.ep-1, 0x1.1p0, 0x1p8, 0x1p16, 10000)
462 SV_POW_INTERVAL2 (0x1p-500, 0x1p500, 0x1p-1, 0x1p1, 10000)
463 /* around estimated argmaxs of ULP error. */
464 SV_POW_INTERVAL2 (0x1p-300, 0x1p-200, 0x1p-20, 0x1p-10, 10000)
465 SV_POW_INTERVAL2 (0x1p50, 0x1p100, 0x1p-20, 0x1p-10, 10000)
466 /* x is negative, y is odd or even integer, or y is real not integer. */
467 TEST_INTERVAL2 (SV_NAME_D2 (pow), -0.0, -10.0, 3.0, 3.0, 10000)
468 TEST_INTERVAL2 (SV_NAME_D2 (pow), -0.0, -10.0, 4.0, 4.0, 10000)
469 TEST_INTERVAL2 (SV_NAME_D2 (pow), -0.0, -10.0, 0.0, 10.0, 10000)
470 TEST_INTERVAL2 (SV_NAME_D2 (pow), 0.0, 10.0, -0.0, -10.0, 10000)
471 /* |x| is inf, y is odd or even integer, or y is real not integer. */
472 SV_POW_INTERVAL2 (inf, inf, 0.5, 0.5, 1)
473 SV_POW_INTERVAL2 (inf, inf, 1.0, 1.0, 1)
474 SV_POW_INTERVAL2 (inf, inf, 2.0, 2.0, 1)
475 SV_POW_INTERVAL2 (inf, inf, 3.0, 3.0, 1)
476 /* 0.0^y. */
477 SV_POW_INTERVAL2 (0.0, 0.0, 0.0, 0x1p120, 1000)
478 /* 1.0^y. */
479 TEST_INTERVAL2 (SV_NAME_D2 (pow), 1.0, 1.0, 0.0, 0x1p-50, 1000)
480 TEST_INTERVAL2 (SV_NAME_D2 (pow), 1.0, 1.0, 0x1p-50, 1.0, 1000)
481 TEST_INTERVAL2 (SV_NAME_D2 (pow), 1.0, 1.0, 1.0, 0x1p100, 1000)
482 TEST_INTERVAL2 (SV_NAME_D2 (pow), 1.0, 1.0, -1.0, -0x1p120, 1000)
483 CLOSE_SVE_ATTR
484