1 /* origin: FreeBSD /usr/src/lib/msun/src/math_private.h */
3 * ====================================================
4 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
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
10 * ====================================================
27 * The original fdlibm code used statements like:
28 * n0 = ((*(int*)&one)>>29)^1; * index of high word *
29 * ix0 = *(n0+(int*)&x); * high word of x *
30 * ix1 = *((1-n0)+(int*)&x); * low word of x *
31 * to dig two 32 bit words out of the 64 bit IEEE floating point
32 * value. That is non-ANSI, and, moreover, the gcc instruction
33 * scheduler gets it wrong. We instead use the following macros.
34 * Unlike the original code, we determine the endianness at compile
35 * time, not at run time; I don't see much benefit to selecting
36 * endianness at run time.
49 /* Get two 32 bit ints from a double. */
50 #define EXTRACT_WORDS(hi,lo,d) \
54 (hi) = __u.bits >> 32; \
55 (lo) = (uint32_t)__u.bits; \
58 /* Get a 64 bit int from a double. */
59 #define EXTRACT_WORD64(i,d) \
66 /* Get the more significant 32 bit int from a double. */
67 #define GET_HIGH_WORD(i,d) \
71 (i) = __u.bits >> 32; \
74 /* Get the less significant 32 bit int from a double. */
75 #define GET_LOW_WORD(i,d) \
79 (i) = (uint32_t)__u.bits; \
82 /* Set a double from two 32 bit ints. */
83 #define INSERT_WORDS(d,hi,lo) \
86 __u.bits = ((uint64_t)(hi) << 32) | (uint32_t)(lo); \
90 /* Set a double from a 64 bit int. */
91 #define INSERT_WORD64(d,i) \
98 /* Set the more significant 32 bits of a double from an int. */
99 #define SET_HIGH_WORD(d,hi) \
103 __u.bits &= 0xffffffff; \
104 __u.bits |= (uint64_t)(hi) << 32; \
108 /* Set the less significant 32 bits of a double from an int. */
109 #define SET_LOW_WORD(d,lo) \
113 __u.bits &= 0xffffffff00000000ull; \
114 __u.bits |= (uint32_t)(lo); \
118 /* Get a 32 bit int from a float. */
119 #define GET_FLOAT_WORD(i,d) \
126 /* Set a float from a 32 bit int. */
127 #define SET_FLOAT_WORD(d,i) \
134 /* fdlibm kernel functions */
135 int __rem_pio2_large(double*,double*,int,int,int);
137 int __rem_pio2(double,double*);
138 double __sin(double,double,int);
139 double __cos(double,double);
140 double __tan(double,double,int);
141 double __ldexp_exp(double,int);
143 double complex __ldexp_cexp(double complex,int);
146 int __rem_pio2f(float,double*);
147 float __sindf(double);
148 float __cosdf(double);
149 float __tandf(double,int);
150 float __ldexp_expf(float,int);
152 float complex __ldexp_cexpf(float complex,int);
155 /* long double precision kernel functions */
156 long double __sinl(long double, long double, int);
157 long double __cosl(long double, long double);
158 long double __tanl(long double, long double, int);
160 /* polynomial evaluation */
161 long double __polevll(long double, long double *, int);
162 long double __p1evll(long double, long double *, int);
166 * Common routine to process the arguments to nan(), nanf(), and nanl().
168 void _scan_nan(uint32_t *__words, int __num_words, const char *__s);
170 // TODO: not needed when -fexcess-precision=standard is supported (>=gcc4.5)
172 * Attempt to get strict C99 semantics for assignment with non-C99 compilers.
175 #define STRICT_ASSIGN(type, lval, rval) do { \
176 volatile type __v = (rval); \
180 #define STRICT_ASSIGN(type, lval, rval) ((lval) = (type)(rval))
186 * C99 specifies that complex numbers have the same representation as
187 * an array of two elements, where the first element is the real part
188 * and the second element is the imaginary part.
199 long double complex f;
201 } long_double_complex;
202 #define REALPART(z) ((z).a[0])
203 #define IMAGPART(z) ((z).a[1])
206 * Inline functions that can be used to construct complex values.
208 * The C99 standard intends x+I*y to be used for this, but x+I*y is
209 * currently unusable in general since gcc introduces many overflow,
210 * underflow, sign and efficiency bugs by rewriting I*y as
211 * (0.0+I)*(y+0.0*I) and laboriously computing the full complex product.
212 * In particular, I*Inf is corrupted to NaN+I*Inf, and I*-0 is corrupted
215 static inline float complex
216 cpackf(float x, float y)
225 static inline double complex
226 cpack(double x, double y)
235 static inline long double complex
236 cpackl(long double x, long double y)
238 long_double_complex z;