1 /* origin: FreeBSD /usr/src/lib/msun/src/math_private.h */
2 /*
3  * ====================================================
4  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
5  *
6  * Developed at SunPro, a Sun Microsystems, Inc. business.
7  * Permission to use, copy, modify, and distribute this
8  * software is freely granted, provided that this notice
9  * is preserved.
10  * ====================================================
11  */
12 
13 #pragma once
14 
15 #include <complex.h>
16 #include <endian.h>
17 #include <float.h>
18 #include <math.h>
19 #include <stdint.h>
20 
21 // Clang does not support the C99 rounding mode pragma. Support seems
22 // unlikely to be coming soon, but for reference the clang/llvm bug
23 // tracking this fact may be found at:
24 // https://llvm.org/bugs/show_bug.cgi?id=8100
25 #define PRAGMA_STDC_FENV_ACCESS_ON
26 
27 // GCC complains about the STDC CX_LIMITED_RANGE pragma.
28 #define PRAGMA_STDC_CX_LIMITED_RANGE_ON
29 
30 #if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
31 #elif LDBL_MANT_DIG == 64 && LDBL_MAX_EXP == 16384 && __BYTE_ORDER == __LITTLE_ENDIAN
32 union ldshape {
33     long double f;
34     struct {
35         uint64_t m;
36         uint16_t se;
37     } i;
38 };
39 #elif LDBL_MANT_DIG == 113 && LDBL_MAX_EXP == 16384 && __BYTE_ORDER == __LITTLE_ENDIAN
40 union ldshape {
41     long double f;
42     struct {
43         uint64_t lo;
44         uint32_t mid;
45         uint16_t top;
46         uint16_t se;
47     } i;
48     struct {
49         uint64_t lo;
50         uint64_t hi;
51     } i2;
52 };
53 #elif LDBL_MANT_DIG == 113 && LDBL_MAX_EXP == 16384 && __BYTE_ORDER == __BIG_ENDIAN
54 union ldshape {
55     long double f;
56     struct {
57         uint16_t se;
58         uint16_t top;
59         uint32_t mid;
60         uint64_t lo;
61     } i;
62     struct {
63         uint64_t hi;
64         uint64_t lo;
65     } i2;
66 };
67 #else
68 #error Unsupported long double representation
69 #endif
70 
71 #define FORCE_EVAL(x)                             \
72     do {                                          \
73         if (sizeof(x) == sizeof(float)) {         \
74             volatile float __x;                   \
75             __x = (x);                            \
76             (void)__x;                            \
77         } else if (sizeof(x) == sizeof(double)) { \
78             volatile double __x;                  \
79             __x = (x);                            \
80             (void)__x;                            \
81         } else {                                  \
82             volatile long double __x;             \
83             __x = (x);                            \
84             (void)__x;                            \
85         }                                         \
86     } while (0)
87 
88 /* Get two 32 bit ints from a double.  */
89 #define EXTRACT_WORDS(hi, lo, d) \
90     do {                         \
91         union {                  \
92             double f;            \
93             uint64_t i;          \
94         } __u;                   \
95         __u.f = (d);             \
96         (hi) = __u.i >> 32;      \
97         (lo) = (uint32_t)__u.i;  \
98     } while (0)
99 
100 /* Get the more significant 32 bit int from a double.  */
101 #define GET_HIGH_WORD(hi, d) \
102     do {                     \
103         union {              \
104             double f;        \
105             uint64_t i;      \
106         } __u;               \
107         __u.f = (d);         \
108         (hi) = __u.i >> 32;  \
109     } while (0)
110 
111 /* Get the less significant 32 bit int from a double.  */
112 #define GET_LOW_WORD(lo, d)     \
113     do {                        \
114         union {                 \
115             double f;           \
116             uint64_t i;         \
117         } __u;                  \
118         __u.f = (d);            \
119         (lo) = (uint32_t)__u.i; \
120     } while (0)
121 
122 /* Set a double from two 32 bit ints.  */
123 #define INSERT_WORDS(d, hi, lo)                          \
124     do {                                                 \
125         union {                                          \
126             double f;                                    \
127             uint64_t i;                                  \
128         } __u;                                           \
129         __u.i = ((uint64_t)(hi) << 32) | (uint32_t)(lo); \
130         (d) = __u.f;                                     \
131     } while (0)
132 
133 /* Set the more significant 32 bits of a double from an int.  */
134 #define SET_HIGH_WORD(d, hi)           \
135     do {                               \
136         union {                        \
137             double f;                  \
138             uint64_t i;                \
139         } __u;                         \
140         __u.f = (d);                   \
141         __u.i &= 0xffffffff;           \
142         __u.i |= (uint64_t)(hi) << 32; \
143         (d) = __u.f;                   \
144     } while (0)
145 
146 /* Set the less significant 32 bits of a double from an int.  */
147 #define SET_LOW_WORD(d, lo)             \
148     do {                                \
149         union {                         \
150             double f;                   \
151             uint64_t i;                 \
152         } __u;                          \
153         __u.f = (d);                    \
154         __u.i &= 0xffffffff00000000ull; \
155         __u.i |= (uint32_t)(lo);        \
156         (d) = __u.f;                    \
157     } while (0)
158 
159 /* Get a 32 bit int from a float.  */
160 #define GET_FLOAT_WORD(w, d) \
161     do {                     \
162         union {              \
163             float f;         \
164             uint32_t i;      \
165         } __u;               \
166         __u.f = (d);         \
167         (w) = __u.i;         \
168     } while (0)
169 
170 /* Set a float from a 32 bit int.  */
171 #define SET_FLOAT_WORD(d, w) \
172     do {                     \
173         union {              \
174             float f;         \
175             uint32_t i;      \
176         } __u;               \
177         __u.i = (w);         \
178         (d) = __u.f;         \
179     } while (0)
180 
181 #undef __CMPLX
182 #undef CMPLX
183 #undef CMPLXF
184 #undef CMPLXL
185 
186 #define __CMPLX(x, y, t)    \
187     ((union {               \
188          _Complex t __z;    \
189          t __xy[2];         \
190      }){.__xy = {(x), (y)}} \
191          .__z)
192 
193 #define CMPLX(x, y) __CMPLX(x, y, double)
194 #define CMPLXF(x, y) __CMPLX(x, y, float)
195 #define CMPLXL(x, y) __CMPLX(x, y, long double)
196 
197 /* fdlibm kernel functions */
198 
199 int __rem_pio2_large(double*, double*, int, int, int);
200 
201 int __rem_pio2(double, double*);
202 double __sin(double, double, int);
203 double __cos(double, double);
204 double __tan(double, double, int);
205 double __expo2(double);
206 double complex __ldexp_cexp(double complex, int);
207 
208 int __rem_pio2f(float, double*);
209 float __sindf(double);
210 float __cosdf(double);
211 float __tandf(double, int);
212 float __expo2f(float);
213 float complex __ldexp_cexpf(float complex, int);
214 
215 int __rem_pio2l(long double, long double*);
216 long double __sinl(long double, long double, int);
217 long double __cosl(long double, long double);
218 long double __tanl(long double, long double, int);
219 
220 /* polynomial evaluation */
221 long double __polevll(long double, const long double*, int);
222 long double __p1evll(long double, const long double*, int);
223