source: trunk/libs/newlib/src/newlib/libm/common/fdlibm.h @ 620

Last change on this file since 620 was 444, checked in by satin@…, 6 years ago

add newlib,libalmos-mkh, restructure shared_syscalls.h and mini-libc
File size: 12.2 KB
Line 
1
2/* @(#)fdlibm.h 5.1 93/09/24 */
3/*
4 * ====================================================
5 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
6 *
7 * Developed at SunPro, a Sun Microsystems, Inc. business.
8 * Permission to use, copy, modify, and distribute this
9 * software is freely granted, provided that this notice
10 * is preserved.
11 * ====================================================
12 */
13
14/* REDHAT LOCAL: Include files.  */
15#include <math.h>
16#include <sys/types.h>
17#include <machine/ieeefp.h>
18
19/* REDHAT LOCAL: Default to XOPEN_MODE.  */
20#define _XOPEN_MODE
21
22/* Most routines need to check whether a float is finite, infinite, or not a
23   number, and many need to know whether the result of an operation will
24   overflow.  These conditions depend on whether the largest exponent is
25   used for NaNs & infinities, or whether it's used for finite numbers.  The
26   macros below wrap up that kind of information:
27
28   FLT_UWORD_IS_FINITE(X)
29        True if a positive float with bitmask X is finite.
30
31   FLT_UWORD_IS_NAN(X)
32        True if a positive float with bitmask X is not a number.
33
34   FLT_UWORD_IS_INFINITE(X)
35        True if a positive float with bitmask X is +infinity.
36
37   FLT_UWORD_MAX
38        The bitmask of FLT_MAX.
39
40   FLT_UWORD_HALF_MAX
41        The bitmask of FLT_MAX/2.
42
43   FLT_UWORD_EXP_MAX
44        The bitmask of the largest finite exponent (129 if the largest
45        exponent is used for finite numbers, 128 otherwise).
46
47   FLT_UWORD_LOG_MAX
48        The bitmask of log(FLT_MAX), rounded down.  This value is the largest
49        input that can be passed to exp() without producing overflow.
50
51   FLT_UWORD_LOG_2MAX
52        The bitmask of log(2*FLT_MAX), rounded down.  This value is the
53        largest input than can be passed to cosh() without producing
54        overflow.
55
56   FLT_LARGEST_EXP
57        The largest biased exponent that can be used for finite numbers
58        (255 if the largest exponent is used for finite numbers, 254
59        otherwise) */
60
61#ifdef _FLT_LARGEST_EXPONENT_IS_NORMAL
62#define FLT_UWORD_IS_FINITE(x) 1
63#define FLT_UWORD_IS_NAN(x) 0
64#define FLT_UWORD_IS_INFINITE(x) 0
65#define FLT_UWORD_MAX 0x7fffffff
66#define FLT_UWORD_EXP_MAX 0x43010000
67#define FLT_UWORD_LOG_MAX 0x42b2d4fc
68#define FLT_UWORD_LOG_2MAX 0x42b437e0
69#define HUGE ((float)0X1.FFFFFEP128)
70#else
71#define FLT_UWORD_IS_FINITE(x) ((x)<0x7f800000L)
72#define FLT_UWORD_IS_NAN(x) ((x)>0x7f800000L)
73#define FLT_UWORD_IS_INFINITE(x) ((x)==0x7f800000L)
74#define FLT_UWORD_MAX 0x7f7fffffL
75#define FLT_UWORD_EXP_MAX 0x43000000
76#define FLT_UWORD_LOG_MAX 0x42b17217
77#define FLT_UWORD_LOG_2MAX 0x42b2d4fc
78#define HUGE ((float)3.40282346638528860e+38)
79#endif
80#define FLT_UWORD_HALF_MAX (FLT_UWORD_MAX-(1L<<23))
81#define FLT_LARGEST_EXP (FLT_UWORD_MAX>>23)
82
83/* Many routines check for zero and subnormal numbers.  Such things depend
84   on whether the target supports denormals or not:
85
86   FLT_UWORD_IS_ZERO(X)
87        True if a positive float with bitmask X is +0.  Without denormals,
88        any float with a zero exponent is a +0 representation.  With
89        denormals, the only +0 representation is a 0 bitmask.
90
91   FLT_UWORD_IS_SUBNORMAL(X)
92        True if a non-zero positive float with bitmask X is subnormal.
93        (Routines should check for zeros first.)
94
95   FLT_UWORD_MIN
96        The bitmask of the smallest float above +0.  Call this number
97        REAL_FLT_MIN...
98
99   FLT_UWORD_EXP_MIN
100        The bitmask of the float representation of REAL_FLT_MIN's exponent.
101
102   FLT_UWORD_LOG_MIN
103        The bitmask of |log(REAL_FLT_MIN)|, rounding down.
104
105   FLT_SMALLEST_EXP
106        REAL_FLT_MIN's exponent - EXP_BIAS (1 if denormals are not supported,
107        -22 if they are).
108*/
109
110#ifdef _FLT_NO_DENORMALS
111#define FLT_UWORD_IS_ZERO(x) ((x)<0x00800000L)
112#define FLT_UWORD_IS_SUBNORMAL(x) 0
113#define FLT_UWORD_MIN 0x00800000
114#define FLT_UWORD_EXP_MIN 0x42fc0000
115#define FLT_UWORD_LOG_MIN 0x42aeac50
116#define FLT_SMALLEST_EXP 1
117#else
118#define FLT_UWORD_IS_ZERO(x) ((x)==0)
119#define FLT_UWORD_IS_SUBNORMAL(x) ((x)<0x00800000L)
120#define FLT_UWORD_MIN 0x00000001
121#define FLT_UWORD_EXP_MIN 0x43160000
122#define FLT_UWORD_LOG_MIN 0x42cff1b5
123#define FLT_SMALLEST_EXP -22
124#endif
125
126#ifdef __STDC__
127#undef __P
128#define __P(p)  p
129#else
130#define __P(p)  ()
131#endif
132
133/*
134 * set X_TLOSS = pi*2**52, which is possibly defined in <values.h>
135 * (one may replace the following line by "#include <values.h>")
136 */
137
138#define X_TLOSS         1.41484755040568800000e+16
139
140/* Functions that are not documented, and are not in <math.h>.  */
141
142#ifdef _SCALB_INT
143extern double scalb __P((double, int));
144#else
145extern double scalb __P((double, double));
146#endif
147extern double significand __P((double));
148
149extern long double __ieee754_hypotl __P((long double, long double));
150
151/* ieee style elementary functions */
152extern double __ieee754_sqrt __P((double));                     
153extern double __ieee754_acos __P((double));                     
154extern double __ieee754_acosh __P((double));                   
155extern double __ieee754_log __P((double));                     
156extern double __ieee754_atanh __P((double));                   
157extern double __ieee754_asin __P((double));                     
158extern double __ieee754_atan2 __P((double,double));                     
159extern double __ieee754_exp __P((double));
160extern double __ieee754_cosh __P((double));
161extern double __ieee754_fmod __P((double,double));
162extern double __ieee754_pow __P((double,double));
163extern double __ieee754_lgamma_r __P((double,int *));
164extern double __ieee754_gamma_r __P((double,int *));
165extern double __ieee754_log10 __P((double));
166extern double __ieee754_sinh __P((double));
167extern double __ieee754_hypot __P((double,double));
168extern double __ieee754_j0 __P((double));
169extern double __ieee754_j1 __P((double));
170extern double __ieee754_y0 __P((double));
171extern double __ieee754_y1 __P((double));
172extern double __ieee754_jn __P((int,double));
173extern double __ieee754_yn __P((int,double));
174extern double __ieee754_remainder __P((double,double));
175extern __int32_t __ieee754_rem_pio2 __P((double,double*));
176#ifdef _SCALB_INT
177extern double __ieee754_scalb __P((double,int));
178#else
179extern double __ieee754_scalb __P((double,double));
180#endif
181
182/* fdlibm kernel function */
183extern double __kernel_standard __P((double,double,int));
184extern double __kernel_sin __P((double,double,int));
185extern double __kernel_cos __P((double,double));
186extern double __kernel_tan __P((double,double,int));
187extern int    __kernel_rem_pio2 __P((double*,double*,int,int,int,const __int32_t*));
188
189/* Undocumented float functions.  */
190#ifdef _SCALB_INT
191extern float scalbf __P((float, int));
192#else
193extern float scalbf __P((float, float));
194#endif
195extern float significandf __P((float));
196
197/* ieee style elementary float functions */
198extern float __ieee754_sqrtf __P((float));                     
199extern float __ieee754_acosf __P((float));                     
200extern float __ieee754_acoshf __P((float));                     
201extern float __ieee754_logf __P((float));                       
202extern float __ieee754_atanhf __P((float));                     
203extern float __ieee754_asinf __P((float));                     
204extern float __ieee754_atan2f __P((float,float));                       
205extern float __ieee754_expf __P((float));
206extern float __ieee754_coshf __P((float));
207extern float __ieee754_fmodf __P((float,float));
208extern float __ieee754_powf __P((float,float));
209extern float __ieee754_lgammaf_r __P((float,int *));
210extern float __ieee754_gammaf_r __P((float,int *));
211extern float __ieee754_log10f __P((float));
212extern float __ieee754_sinhf __P((float));
213extern float __ieee754_hypotf __P((float,float));
214extern float __ieee754_j0f __P((float));
215extern float __ieee754_j1f __P((float));
216extern float __ieee754_y0f __P((float));
217extern float __ieee754_y1f __P((float));
218extern float __ieee754_jnf __P((int,float));
219extern float __ieee754_ynf __P((int,float));
220extern float __ieee754_remainderf __P((float,float));
221extern __int32_t __ieee754_rem_pio2f __P((float,float*));
222#ifdef _SCALB_INT
223extern float __ieee754_scalbf __P((float,int));
224#else
225extern float __ieee754_scalbf __P((float,float));
226#endif
227
228#if !__OBSOLETE_MATH
229/* The new math code does not provide separate wrapper function
230   for error handling, so the extern symbol is called directly.
231   This is valid as long as there are no namespace issues (the
232   extern symbol is reserved whenever the caller is reserved)
233   and there are no observable error handling side effects.  */
234# define __ieee754_expf(x) expf(x)
235# define __ieee754_logf(x) logf(x)
236# define __ieee754_powf(x,y) powf(x,y)
237#endif
238
239/* float versions of fdlibm kernel functions */
240extern float __kernel_sinf __P((float,float,int));
241extern float __kernel_cosf __P((float,float));
242extern float __kernel_tanf __P((float,float,int));
243extern int   __kernel_rem_pio2f __P((float*,float*,int,int,int,const __int32_t*));
244
245/* The original code used statements like
246        n0 = ((*(int*)&one)>>29)^1;             * index of high word *
247        ix0 = *(n0+(int*)&x);                   * high word of x *
248        ix1 = *((1-n0)+(int*)&x);               * low word of x *
249   to dig two 32 bit words out of the 64 bit IEEE floating point
250   value.  That is non-ANSI, and, moreover, the gcc instruction
251   scheduler gets it wrong.  We instead use the following macros.
252   Unlike the original code, we determine the endianness at compile
253   time, not at run time; I don't see much benefit to selecting
254   endianness at run time.  */
255
256#ifndef __IEEE_BIG_ENDIAN
257#ifndef __IEEE_LITTLE_ENDIAN
258 #error Must define endianness
259#endif
260#endif
261
262/* A union which permits us to convert between a double and two 32 bit
263   ints.  */
264
265#ifdef __IEEE_BIG_ENDIAN
266
267typedef union 
268{
269  double value;
270  struct 
271  {
272    __uint32_t msw;
273    __uint32_t lsw;
274  } parts;
275} ieee_double_shape_type;
276
277#endif
278
279#ifdef __IEEE_LITTLE_ENDIAN
280
281typedef union 
282{
283  double value;
284  struct 
285  {
286    __uint32_t lsw;
287    __uint32_t msw;
288  } parts;
289} ieee_double_shape_type;
290
291#endif
292
293/* Get two 32 bit ints from a double.  */
294
295#define EXTRACT_WORDS(ix0,ix1,d)                                \
296do {                                                            \
297  ieee_double_shape_type ew_u;                                  \
298  ew_u.value = (d);                                             \
299  (ix0) = ew_u.parts.msw;                                       \
300  (ix1) = ew_u.parts.lsw;                                       \
301} while (0)
302
303/* Get the more significant 32 bit int from a double.  */
304
305#define GET_HIGH_WORD(i,d)                                      \
306do {                                                            \
307  ieee_double_shape_type gh_u;                                  \
308  gh_u.value = (d);                                             \
309  (i) = gh_u.parts.msw;                                         \
310} while (0)
311
312/* Get the less significant 32 bit int from a double.  */
313
314#define GET_LOW_WORD(i,d)                                       \
315do {                                                            \
316  ieee_double_shape_type gl_u;                                  \
317  gl_u.value = (d);                                             \
318  (i) = gl_u.parts.lsw;                                         \
319} while (0)
320
321/* Set a double from two 32 bit ints.  */
322
323#define INSERT_WORDS(d,ix0,ix1)                                 \
324do {                                                            \
325  ieee_double_shape_type iw_u;                                  \
326  iw_u.parts.msw = (ix0);                                       \
327  iw_u.parts.lsw = (ix1);                                       \
328  (d) = iw_u.value;                                             \
329} while (0)
330
331/* Set the more significant 32 bits of a double from an int.  */
332
333#define SET_HIGH_WORD(d,v)                                      \
334do {                                                            \
335  ieee_double_shape_type sh_u;                                  \
336  sh_u.value = (d);                                             \
337  sh_u.parts.msw = (v);                                         \
338  (d) = sh_u.value;                                             \
339} while (0)
340
341/* Set the less significant 32 bits of a double from an int.  */
342
343#define SET_LOW_WORD(d,v)                                       \
344do {                                                            \
345  ieee_double_shape_type sl_u;                                  \
346  sl_u.value = (d);                                             \
347  sl_u.parts.lsw = (v);                                         \
348  (d) = sl_u.value;                                             \
349} while (0)
350
351/* A union which permits us to convert between a float and a 32 bit
352   int.  */
353
354typedef union
355{
356  float value;
357  __uint32_t word;
358} ieee_float_shape_type;
359
360/* Get a 32 bit int from a float.  */
361
362#define GET_FLOAT_WORD(i,d)                                     \
363do {                                                            \
364  ieee_float_shape_type gf_u;                                   \
365  gf_u.value = (d);                                             \
366  (i) = gf_u.word;                                              \
367} while (0)
368
369/* Set a float from a 32 bit int.  */
370
371#define SET_FLOAT_WORD(d,i)                                     \
372do {                                                            \
373  ieee_float_shape_type sf_u;                                   \
374  sf_u.word = (i);                                              \
375  (d) = sf_u.value;                                             \
376} while (0)
377
378/* Macros to avoid undefined behaviour that can arise if the amount
379   of a shift is exactly equal to the size of the shifted operand.  */
380
381#define SAFE_LEFT_SHIFT(op,amt)                                 \
382  (((amt) < 8 * sizeof(op)) ? ((op) << (amt)) : 0)
383
384#define SAFE_RIGHT_SHIFT(op,amt)                                \
385  (((amt) < 8 * sizeof(op)) ? ((op) >> (amt)) : 0)
386
387#ifdef  _COMPLEX_H
388
389/*
390 * Quoting from ISO/IEC 9899:TC2:
391 *
392 * 6.2.5.13 Types
393 * Each complex type has the same representation and alignment requirements as
394 * an array type containing exactly two elements of the corresponding real type;
395 * the first element is equal to the real part, and the second element to the
396 * imaginary part, of the complex number.
397 */
398typedef union {
399        float complex z;
400        float parts[2];
401} float_complex;
402
403typedef union {
404        double complex z;
405        double parts[2];
406} double_complex;
407
408typedef union {
409        long double complex z;
410        long double parts[2];
411} long_double_complex;
412
413#define REAL_PART(z)    ((z).parts[0])
414#define IMAG_PART(z)    ((z).parts[1])
415
416#endif  /* _COMPLEX_H */
417
Note: See TracBrowser for help on using the repository browser.