1 /* Configuration for double precision math routines.
2 Copyright (C) 2018-2022 Free Software Foundation, Inc.
3 This file is part of the GNU C Library.
4
5 The GNU C Library is free software; you can redistribute it and/or
6 modify it under the terms of the GNU Lesser General Public
7 License as published by the Free Software Foundation; either
8 version 2.1 of the License, or (at your option) any later version.
9
10 The GNU C Library is distributed in the hope that it will be useful,
11 but WITHOUT ANY WARRANTY; without even the implied warranty of
12 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13 Lesser General Public License for more details.
14
15 You should have received a copy of the GNU Lesser General Public
16 License along with the GNU C Library; if not, see
17 <https://www.gnu.org/licenses/>. */
18
19 #ifndef _MATH_CONFIG_H
20 #define _MATH_CONFIG_H
21
22 #include <math.h>
23 #include <math_private.h>
24 #include <nan-high-order-bit.h>
25 #include <stdint.h>
26
27 #ifndef WANT_ROUNDING
28 /* Correct special case results in non-nearest rounding modes. */
29 # define WANT_ROUNDING 1
30 #endif
31 #ifndef WANT_ERRNO
32 /* Set errno according to ISO C with (math_errhandling & MATH_ERRNO) != 0. */
33 # define WANT_ERRNO 1
34 #endif
35 #ifndef WANT_ERRNO_UFLOW
36 /* Set errno to ERANGE if result underflows to 0 (in all rounding modes). */
37 # define WANT_ERRNO_UFLOW (WANT_ROUNDING && WANT_ERRNO)
38 #endif
39
40 #ifndef TOINT_INTRINSICS
41 /* When set, the roundtoint and converttoint functions are provided with
42 the semantics documented below. */
43 # define TOINT_INTRINSICS 0
44 #endif
45
46 #if TOINT_INTRINSICS
47 /* Round x to nearest int in all rounding modes, ties have to be rounded
48 consistently with converttoint so the results match. If the result
49 would be outside of [-2^31, 2^31-1] then the semantics is unspecified. */
50 static inline double_t
51 roundtoint (double_t x);
52
53 /* Convert x to nearest int in all rounding modes, ties have to be rounded
54 consistently with roundtoint. If the result is not representible in an
55 int32_t then the semantics is unspecified. */
56 static inline int32_t
57 converttoint (double_t x);
58 #endif
59
60 static inline uint64_t
asuint64(double f)61 asuint64 (double f)
62 {
63 union
64 {
65 double f;
66 uint64_t i;
67 } u = {f};
68 return u.i;
69 }
70
71 static inline double
asdouble(uint64_t i)72 asdouble (uint64_t i)
73 {
74 union
75 {
76 uint64_t i;
77 double f;
78 } u = {i};
79 return u.f;
80 }
81
82 static inline int
issignaling_inline(double x)83 issignaling_inline (double x)
84 {
85 uint64_t ix = asuint64 (x);
86 if (HIGH_ORDER_BIT_IS_SET_FOR_SNAN)
87 return (ix & 0x7ff8000000000000) == 0x7ff8000000000000;
88 return 2 * (ix ^ 0x0008000000000000) > 2 * 0x7ff8000000000000ULL;
89 }
90
91 #define NOINLINE __attribute__ ((noinline))
92
93 /* Error handling tail calls for special cases, with a sign argument.
94 The sign of the return value is set if the argument is non-zero. */
95
96 /* The result overflows. */
97 attribute_hidden double __math_oflow (uint32_t);
98 /* The result underflows to 0 in nearest rounding mode. */
99 attribute_hidden double __math_uflow (uint32_t);
100 /* The result underflows to 0 in some directed rounding mode only. */
101 attribute_hidden double __math_may_uflow (uint32_t);
102 /* Division by zero. */
103 attribute_hidden double __math_divzero (uint32_t);
104
105 /* Error handling using input checking. */
106
107 /* Invalid input unless it is a quiet NaN. */
108 attribute_hidden double __math_invalid (double);
109
110 /* Error handling using output checking, only for errno setting. */
111
112 /* Check if the result overflowed to infinity. */
113 attribute_hidden double __math_check_oflow (double);
114 /* Check if the result underflowed to 0. */
115 attribute_hidden double __math_check_uflow (double);
116
117 /* Check if the result overflowed to infinity. */
118 static inline double
check_oflow(double x)119 check_oflow (double x)
120 {
121 return WANT_ERRNO ? __math_check_oflow (x) : x;
122 }
123
124 /* Check if the result underflowed to 0. */
125 static inline double
check_uflow(double x)126 check_uflow (double x)
127 {
128 return WANT_ERRNO ? __math_check_uflow (x) : x;
129 }
130
131 #define EXP_TABLE_BITS 7
132 #define EXP_POLY_ORDER 5
133 #define EXP2_POLY_ORDER 5
134 extern const struct exp_data
135 {
136 double invln2N;
137 double shift;
138 double negln2hiN;
139 double negln2loN;
140 double poly[4]; /* Last four coefficients. */
141 double exp2_shift;
142 double exp2_poly[EXP2_POLY_ORDER];
143 uint64_t tab[2*(1 << EXP_TABLE_BITS)];
144 } __exp_data attribute_hidden;
145
146 #define LOG_TABLE_BITS 7
147 #define LOG_POLY_ORDER 6
148 #define LOG_POLY1_ORDER 12
149 extern const struct log_data
150 {
151 double ln2hi;
152 double ln2lo;
153 double poly[LOG_POLY_ORDER - 1]; /* First coefficient is 1. */
154 double poly1[LOG_POLY1_ORDER - 1];
155 /* See e_log_data.c for details. */
156 struct {double invc, logc;} tab[1 << LOG_TABLE_BITS];
157 #ifndef __FP_FAST_FMA
158 struct {double chi, clo;} tab2[1 << LOG_TABLE_BITS];
159 #endif
160 } __log_data attribute_hidden;
161
162 #define LOG2_TABLE_BITS 6
163 #define LOG2_POLY_ORDER 7
164 #define LOG2_POLY1_ORDER 11
165 extern const struct log2_data
166 {
167 double invln2hi;
168 double invln2lo;
169 double poly[LOG2_POLY_ORDER - 1];
170 double poly1[LOG2_POLY1_ORDER - 1];
171 /* See e_log2_data.c for details. */
172 struct {double invc, logc;} tab[1 << LOG2_TABLE_BITS];
173 #ifndef __FP_FAST_FMA
174 struct {double chi, clo;} tab2[1 << LOG2_TABLE_BITS];
175 #endif
176 } __log2_data attribute_hidden;
177
178 #define POW_LOG_TABLE_BITS 7
179 #define POW_LOG_POLY_ORDER 8
180 extern const struct pow_log_data
181 {
182 double ln2hi;
183 double ln2lo;
184 double poly[POW_LOG_POLY_ORDER - 1]; /* First coefficient is 1. */
185 /* Note: the pad field is unused, but allows slightly faster indexing. */
186 /* See e_pow_log_data.c for details. */
187 struct {double invc, pad, logc, logctail;} tab[1 << POW_LOG_TABLE_BITS];
188 } __pow_log_data attribute_hidden;
189
190 #endif
191