2 * Configuration for math routines.
4 * Copyright (c) 2017-2020, Arm Limited.
5 * SPDX-License-Identifier: MIT
15 /* If defined to 1, return correct results for special cases in non-nearest
16 rounding modes (logf (1.0f) returns 0.0f with FE_DOWNWARD rather than -0.0f).
17 This may be set to 0 if there is no fenv support or if math functions only
18 get called in round to nearest mode. */
19 # define WANT_ROUNDING 1
22 /* If defined to 1, set errno in math functions according to ISO C. Many math
23 libraries do not set errno, so this is 0 by default. It may need to be
24 set to 1 if math.h has (math_errhandling & MATH_ERRNO) != 0. */
27 #ifndef WANT_ERRNO_UFLOW
28 /* Set errno to ERANGE if result underflows to 0 (in all rounding modes). */
29 # define WANT_ERRNO_UFLOW (WANT_ROUNDING && WANT_ERRNO)
32 /* Compiler can inline round as a single instruction. */
33 #ifndef HAVE_FAST_ROUND
35 # define HAVE_FAST_ROUND 1
37 # define HAVE_FAST_ROUND 0
41 /* Compiler can inline lround, but not (long)round(x). */
42 #ifndef HAVE_FAST_LROUND
43 # if __aarch64__ && (100*__GNUC__ + __GNUC_MINOR__) >= 408 && __NO_MATH_ERRNO__
44 # define HAVE_FAST_LROUND 1
46 # define HAVE_FAST_LROUND 0
50 /* Compiler can inline fma as a single instruction. */
52 # if defined FP_FAST_FMA || __aarch64__
53 # define HAVE_FAST_FMA 1
55 # define HAVE_FAST_FMA 0
59 /* Provide *_finite symbols and some of the glibc hidden symbols
60 so libmathlib can be used with binaries compiled against glibc
61 to interpose math functions with both static and dynamic linking. */
64 # define USE_GLIBC_ABI 1
66 # define USE_GLIBC_ABI 0
70 /* Optionally used extensions. */
72 # define HIDDEN __attribute__ ((__visibility__ ("hidden")))
73 # define NOINLINE __attribute__ ((noinline))
74 # define UNUSED __attribute__ ((unused))
75 # define likely(x) __builtin_expect (!!(x), 1)
76 # define unlikely(x) __builtin_expect (x, 0)
78 # define attribute_copy(f) __attribute__ ((copy (f)))
80 # define attribute_copy(f)
82 # define strong_alias(f, a) \
83 extern __typeof (f) a __attribute__ ((alias (#f))) attribute_copy (f);
84 # define hidden_alias(f, a) \
85 extern __typeof (f) a __attribute__ ((alias (#f), visibility ("hidden"))) \
91 # define likely(x) (x)
92 # define unlikely(x) (x)
96 /* When set, the roundtoint and converttoint functions are provided with
97 the semantics documented below. */
98 # define TOINT_INTRINSICS 1
100 /* Round x to nearest int in all rounding modes, ties have to be rounded
101 consistently with converttoint so the results match. If the result
102 would be outside of [-2^31, 2^31-1] then the semantics is unspecified. */
103 static inline double_t
104 roundtoint (double_t x)
109 /* Convert x to nearest int in all rounding modes, ties have to be rounded
110 consistently with roundtoint. If the result is not representible in an
111 int32_t then the semantics is unspecified. */
112 static inline int32_t
113 converttoint (double_t x)
115 # if HAVE_FAST_LROUND
118 return (long) round (x);
123 static inline uint32_t
145 static inline uint64_t
157 asdouble (uint64_t i)
167 #ifndef IEEE_754_2008_SNAN
168 # define IEEE_754_2008_SNAN 1
171 issignalingf_inline (float x)
173 uint32_t ix = asuint (x);
174 if (!IEEE_754_2008_SNAN)
175 return (ix & 0x7fc00000) == 0x7fc00000;
176 return 2 * (ix ^ 0x00400000) > 2u * 0x7fc00000;
180 issignaling_inline (double x)
182 uint64_t ix = asuint64 (x);
183 if (!IEEE_754_2008_SNAN)
184 return (ix & 0x7ff8000000000000) == 0x7ff8000000000000;
185 return 2 * (ix ^ 0x0008000000000000) > 2 * 0x7ff8000000000000ULL;
188 #if __aarch64__ && __GNUC__
189 /* Prevent the optimization of a floating-point expression. */
191 opt_barrier_float (float x)
193 __asm__ __volatile__ ("" : "+w" (x));
197 opt_barrier_double (double x)
199 __asm__ __volatile__ ("" : "+w" (x));
202 /* Force the evaluation of a floating-point expression for its side-effect. */
204 force_eval_float (float x)
206 __asm__ __volatile__ ("" : "+w" (x));
209 force_eval_double (double x)
211 __asm__ __volatile__ ("" : "+w" (x));
215 opt_barrier_float (float x)
217 volatile float y = x;
221 opt_barrier_double (double x)
223 volatile double y = x;
227 force_eval_float (float x)
229 volatile float y UNUSED = x;
232 force_eval_double (double x)
234 volatile double y UNUSED = x;
238 /* Evaluate an expression as the specified type, normally a type
239 cast should be enough, but compilers implement non-standard
240 excess-precision handling, so when FLT_EVAL_METHOD != 0 then
241 these functions may need to be customized. */
243 eval_as_float (float x)
248 eval_as_double (double x)
253 /* Error handling tail calls for special cases, with a sign argument.
254 The sign of the return value is set if the argument is non-zero. */
256 /* The result overflows. */
257 HIDDEN float __math_oflowf (uint32_t);
258 /* The result underflows to 0 in nearest rounding mode. */
259 HIDDEN float __math_uflowf (uint32_t);
260 /* The result underflows to 0 in some directed rounding mode only. */
261 HIDDEN float __math_may_uflowf (uint32_t);
262 /* Division by zero. */
263 HIDDEN float __math_divzerof (uint32_t);
264 /* The result overflows. */
265 HIDDEN double __math_oflow (uint32_t);
266 /* The result underflows to 0 in nearest rounding mode. */
267 HIDDEN double __math_uflow (uint32_t);
268 /* The result underflows to 0 in some directed rounding mode only. */
269 HIDDEN double __math_may_uflow (uint32_t);
270 /* Division by zero. */
271 HIDDEN double __math_divzero (uint32_t);
273 /* Error handling using input checking. */
275 /* Invalid input unless it is a quiet NaN. */
276 HIDDEN float __math_invalidf (float);
277 /* Invalid input unless it is a quiet NaN. */
278 HIDDEN double __math_invalid (double);
280 /* Error handling using output checking, only for errno setting. */
282 /* Check if the result overflowed to infinity. */
283 HIDDEN double __math_check_oflow (double);
284 /* Check if the result underflowed to 0. */
285 HIDDEN double __math_check_uflow (double);
287 /* Check if the result overflowed to infinity. */
289 check_oflow (double x)
291 return WANT_ERRNO ? __math_check_oflow (x) : x;
294 /* Check if the result underflowed to 0. */
296 check_uflow (double x)
298 return WANT_ERRNO ? __math_check_uflow (x) : x;
301 /* Check if the result overflowed to infinity. */
302 HIDDEN float __math_check_oflowf (float);
303 /* Check if the result underflowed to 0. */
304 HIDDEN float __math_check_uflowf (float);
306 /* Check if the result overflowed to infinity. */
308 check_oflowf (float x)
310 return WANT_ERRNO ? __math_check_oflowf (x) : x;
313 /* Check if the result underflowed to 0. */
315 check_uflowf (float x)
317 return WANT_ERRNO ? __math_check_uflowf (x) : x;
320 /* Shared between expf, exp2f and powf. */
321 #define EXP2F_TABLE_BITS 5
322 #define EXP2F_POLY_ORDER 3
323 extern const struct exp2f_data
325 uint64_t tab[1 << EXP2F_TABLE_BITS];
327 double poly[EXP2F_POLY_ORDER];
329 double invln2_scaled;
330 double poly_scaled[EXP2F_POLY_ORDER];
331 } __exp2f_data HIDDEN;
333 #define LOGF_TABLE_BITS 4
334 #define LOGF_POLY_ORDER 4
335 extern const struct logf_data
340 } tab[1 << LOGF_TABLE_BITS];
342 double poly[LOGF_POLY_ORDER - 1]; /* First order coefficient is 1. */
343 } __logf_data HIDDEN;
345 #define LOG2F_TABLE_BITS 4
346 #define LOG2F_POLY_ORDER 4
347 extern const struct log2f_data
352 } tab[1 << LOG2F_TABLE_BITS];
353 double poly[LOG2F_POLY_ORDER];
354 } __log2f_data HIDDEN;
356 #define POWF_LOG2_TABLE_BITS 4
357 #define POWF_LOG2_POLY_ORDER 5
359 # define POWF_SCALE_BITS EXP2F_TABLE_BITS
361 # define POWF_SCALE_BITS 0
363 #define POWF_SCALE ((double) (1 << POWF_SCALE_BITS))
364 extern const struct powf_log2_data
369 } tab[1 << POWF_LOG2_TABLE_BITS];
370 double poly[POWF_LOG2_POLY_ORDER];
371 } __powf_log2_data HIDDEN;
374 #define EXP_TABLE_BITS 7
375 #define EXP_POLY_ORDER 5
376 /* Use polynomial that is optimized for a wider input range. This may be
377 needed for good precision in non-nearest rounding and !TOINT_INTRINSICS. */
378 #define EXP_POLY_WIDE 0
379 /* Use close to nearest rounding toint when !TOINT_INTRINSICS. This may be
380 needed for good precision in non-nearest rouning and !EXP_POLY_WIDE. */
381 #define EXP_USE_TOINT_NARROW 0
382 #define EXP2_POLY_ORDER 5
383 #define EXP2_POLY_WIDE 0
384 extern const struct exp_data
390 double poly[4]; /* Last four coefficients. */
392 double exp2_poly[EXP2_POLY_ORDER];
393 uint64_t tab[2*(1 << EXP_TABLE_BITS)];
396 #define LOG_TABLE_BITS 7
397 #define LOG_POLY_ORDER 6
398 #define LOG_POLY1_ORDER 12
399 extern const struct log_data
403 double poly[LOG_POLY_ORDER - 1]; /* First coefficient is 1. */
404 double poly1[LOG_POLY1_ORDER - 1];
405 struct {double invc, logc;} tab[1 << LOG_TABLE_BITS];
407 struct {double chi, clo;} tab2[1 << LOG_TABLE_BITS];
411 #define LOG2_TABLE_BITS 6
412 #define LOG2_POLY_ORDER 7
413 #define LOG2_POLY1_ORDER 11
414 extern const struct log2_data
418 double poly[LOG2_POLY_ORDER - 1];
419 double poly1[LOG2_POLY1_ORDER - 1];
420 struct {double invc, logc;} tab[1 << LOG2_TABLE_BITS];
422 struct {double chi, clo;} tab2[1 << LOG2_TABLE_BITS];
424 } __log2_data HIDDEN;
426 #define POW_LOG_TABLE_BITS 7
427 #define POW_LOG_POLY_ORDER 8
428 extern const struct pow_log_data
432 double poly[POW_LOG_POLY_ORDER - 1]; /* First coefficient is 1. */
433 /* Note: the pad field is unused, but allows slightly faster indexing. */
434 struct {double invc, pad, logc, logctail;} tab[1 << POW_LOG_TABLE_BITS];
435 } __pow_log_data HIDDEN;
437 extern const struct erff_data
439 float erff_poly_A[6];
440 float erff_poly_B[7];
441 } __erff_data HIDDEN;
443 #define ERF_POLY_A_ORDER 19
444 #define ERF_POLY_A_NCOEFFS 10
445 #define ERFC_POLY_C_NCOEFFS 16
446 #define ERFC_POLY_D_NCOEFFS 18
447 #define ERFC_POLY_E_NCOEFFS 14
448 #define ERFC_POLY_F_NCOEFFS 17
449 extern const struct erf_data
451 double erf_poly_A[ERF_POLY_A_NCOEFFS];
452 double erf_ratio_N_A[5];
453 double erf_ratio_D_A[5];
454 double erf_ratio_N_B[7];
455 double erf_ratio_D_B[6];
456 double erfc_poly_C[ERFC_POLY_C_NCOEFFS];
457 double erfc_poly_D[ERFC_POLY_D_NCOEFFS];
458 double erfc_poly_E[ERFC_POLY_E_NCOEFFS];
459 double erfc_poly_F[ERFC_POLY_F_NCOEFFS];