overby.me
this repo has no description
1/*
2 * Copyright (c) 2002 Apple Computer, Inc. All rights reserved.
3 *
4 * @APPLE_LICENSE_HEADER_START@
5 *
6 * The contents of this file constitute Original Code as defined in and
7 * are subject to the Apple Public Source License Version 1.1 (the
8 * "License"). You may not use this file except in compliance with the
9 * License. Please obtain a copy of the License at
10 * http://www.apple.com/publicsource and read it before using this file.
11 *
12 * This Original Code and all software distributed under the License are
13 * distributed on an "AS IS" basis, WITHOUT WARRANTY OF ANY KIND, EITHER
14 * EXPRESS OR IMPLIED, AND APPLE HEREBY DISCLAIMS ALL SUCH WARRANTIES,
15 * INCLUDING WITHOUT LIMITATION, ANY WARRANTIES OF MERCHANTABILITY,
16 * FITNESS FOR A PARTICULAR PURPOSE OR NON-INFRINGEMENT. Please see the
17 * License for the specific language governing rights and limitations
18 * under the License.
19 *
20 * @APPLE_LICENSE_HEADER_END@
21 */
22/*******************************************************************************
23* *
24* File: fpmacros.c *
25* *
26* Contains: C source code for PowerPC implementations of the inquiry *
27* macros, as defined in C99. *
28* *
29* Copyright � 1992-2001 Apple Computer, Inc. All rights reserved. *
30* *
31* Written by A. Sazegari and Jon Okada, started on November 29 1992. *
32* Modified and ported by Robert A. Murley (ram) for Mac OS X. *
33* *
34* fpmacros is a new file that brings all of C99 macros together. *
35* *
36* A MathLib v4 file. *
37* *
38* Macros __isnormald, __isfinited, __isnand and __inf were previously *
39* in fp.c. *
40* Macros __fpclassifyf, __isnormalf, __isfinitef, __isnanf and *
41* __signbitf were previously in fpfloatfunc.c. *
42* Macro __fpclassifyd was in classify.c. *
43* Macro __signbitd was in sign.c. *
44* *
45* Change History (most recent first): *
46* *
47* 06 Nov 01 ram commented out warning about Intel architectures. *
48* changed i386 stubs to call abort(). *
49* 02 Nov 01 ram added stubs for i386 routines. *
50* 08 Oct 01 ram removed <CoreServices/CoreServices.h>. *
51* changed compiler errors to warnings. *
52* 24 Sep 01 ram corrected mantissa mask in fpclassifyf and isnanf. *
53* 18 Sep 01 ali added <CoreServices/CoreServices.h> to get to <fp.h>. *
54* 10 Sep 01 ali added macros to detect PowerPC and correct compiler. *
55* 09 Sep 01 ali added more comments. *
56* 05 Sep 01 ram added __inf routine. *
57* added #ifdef __ppc__. *
58* 07 Jul 01 ram first created from fpfloatfunc.c, fp.c, *
59* classify.c and sign.c in MathLib v3 Mac OS9. *
60* replaced DblInHex typedef with hexdouble. *
61* *
62* W A R N I N G: *
63* These routines require a 64-bit double precision IEEE-754 model. *
64* They are written for PowerPC only and are expecting the compiler *
65* to generate the correct sequence of multiply-add fused instructions. *
66* *
67* These routines are not intended for 32-bit Intel architectures. *
68* *
69* A version of gcc higher than 932 is required. *
70* *
71* GCC compiler options: *
72* optimization level 3 (-O3) *
73* -fschedule-insns -finline-functions -funroll-all-loops *
74* *
75*******************************************************************************/
76#include "math.h"
77#include "float.h"
78#include "fp_private.h"
79#include "fenv.h"
80
81#if !defined(BUILDING_FOR_CARBONCORE_LEGACY)
82
83/******************************************************************************
84* No other functions are called by these routines outside of fpmacros.c. *
85******************************************************************************/
86
87int __math_errhandling ( void )
88{
89 return (MATH_ERREXCEPT); // return the bitmask indicating the error discipline(s) in use.
90}
91
92/**************************************************************************
93 Function __fpclassifyf
94 Returns the classification code of the argument float x, as defined in
95 C99.
96**************************************************************************/
97
98int __fpclassifyf ( float x )
99{
100 uint32_t iexp;
101 hexsingle z;
102
103 z.fval = x;
104 __NOOP;
105 __NOOP;
106 __NOOP;
107 iexp = z.lval & 0x7f800000; // isolate float exponent
108
109 if (iexp == 0x7f800000) { // NaN or INF case
110 if ((z.lval & 0x007fffff) == 0)
111 return FP_INFINITE;
112 else if ((z.lval & fQuietNan) != 0)
113 return FP_QNAN;
114 else
115 return FP_SNAN;
116 }
117
118 if (iexp != 0) // normal float
119 return FP_NORMAL;
120
121 if ((z.lval & 0x007fffff) == 0)
122 return FP_ZERO; // zero
123 else
124 return FP_SUBNORMAL; //must be subnormal
125}
126
127
128/*************************************************************************
129 Function __fpclassifyd
130 Returns the classification code of the argument double x, as
131 defined in C99.
132*************************************************************************/
133
134int __fpclassifyd ( double arg )
135{
136 uint32_t exponent;
137 hexdouble x;
138
139 x.d = arg;
140 __NOOP;
141 __NOOP;
142 __NOOP;
143
144 exponent = x.i.hi & 0x7ff00000;
145 if ( exponent == 0x7ff00000 )
146 {
147 if ( ( ( x.i.hi & 0x000fffff ) | x.i.lo ) == 0 )
148 return FP_INFINITE;
149 else
150 return ( x.i.hi & dQuietNan ) ? FP_QNAN : FP_SNAN;
151 }
152 else if ( exponent != 0)
153 return FP_NORMAL;
154 else
155 {
156 if ( ( ( x.i.hi & 0x000fffff ) | x.i.lo ) == 0 )
157 return FP_ZERO;
158 else
159 return FP_SUBNORMAL;
160 }
161}
162
163static const double twoTo52 = 0x1.0p+52; // 2^52
164static const double mtwoTo52 = -0x1.0p+52;
165
166static inline double fastfloor( double f ) __attribute__((always_inline));
167static inline double fastfloor( double f )
168{
169 double b, c, d, e, g, h, t;
170
171 c = __fsel( f, mtwoTo52, twoTo52 ); b = fabs( f );
172 d = (f - c) + c; e = b - twoTo52;
173 g = f - d;
174 h = __fsel( g, 0.0, 1.0 );
175 t = d - h;
176 return __fsel( e, f, t );
177}
178
179// Called from fpclassify macro iff sizeof(double) != sizeof(long double)
180int __fpclassify ( long double arg )
181{
182 uint32_t exponent;
183 int e;
184 hexdouble x, xx;
185 hexdbldbl u;
186
187 u.dd.head = 0.0;
188 u.dd.tail = 0.0;
189
190 u.ld = arg;
191 __NOOP;
192 __NOOP;
193 __NOOP;
194
195 x.d = u.dd.head;
196 xx.d = fabs( u.dd.tail );
197 __NOOP;
198 __NOOP;
199 __NOOP;
200
201 exponent = x.i.hi & 0x7ff00000;
202 e = (exponent >> 20) - 1023; // unbiased exponent
203
204 if ( exponent == 0x7ff00000 )
205 {
206 if ( ( ( x.i.hi & 0x000fffff ) | x.i.lo ) == 0 )
207 return FP_INFINITE;
208 else
209 return ( x.i.hi & dQuietNan ) ? FP_QNAN : FP_SNAN;
210 }
211 else if ( exponent == 0 )
212 {
213 if ( ( ( x.i.hi & 0x000fffff ) | x.i.lo ) == 0 )
214 return FP_ZERO;
215 else
216 return FP_SUBNORMAL;
217 }
218 else if ( e < (LDBL_MIN_EXP - 1) )
219 {
220 return FP_SUBNORMAL;
221 }
222 else if ( xx.d == 0.0 )
223 {
224 return FP_NORMAL;
225 }
226 else /* xx.d > 0 */
227 {
228 double t, u;
229 int n = LDBL_MANT_DIG - e; // scaling arg by 2^n puts LDBL_MANT_DIG bits left of the binary point
230
231 t = ldexp( xx.d, n ); // t holds the tails contribution to the last LDBL_MANT_DIG bits
232 if ( t == 0.0 )
233 return FP_SUPERNORMAL; // contribution of xx.d entirely beyond LDBL_MANT_DIG places to the right
234
235 u = fastfloor( t ); // chop off any "SUPERNORMAL" bits
236
237 if ( t == u )
238 return FP_NORMAL; // no "SUPERNORMAL" bits detected
239 else
240 return FP_SUPERNORMAL;
241 }
242}
243
244int __issupernormal ( long double x )
245{
246 return (FP_SUPERNORMAL == __fpclassify( x ));
247}
248
249/***********************************************************************
250 Function __isnormalf
251 Returns nonzero if and only if x is a normalized float number and
252 zero otherwise.
253***********************************************************************/
254
255int __isnormalf ( float x )
256{
257 uint32_t iexp;
258 hexsingle z;
259
260 z.fval = x;
261 __NOOP;
262 __NOOP;
263 __NOOP;
264 iexp = z.lval & 0x7f800000; /* isolate float exponent */
265 return ((iexp != 0x7f800000) && (iexp != 0));
266}
267
268
269/***********************************************************************
270 Function __isnormald
271 Returns nonzero if and only if x is a normalized double number and
272 zero otherwise.
273***********************************************************************/
274
275int __isnormald ( double x )
276{
277 uint32_t iexp;
278 hexdouble z;
279
280 z.d = x;
281 __NOOP;
282 __NOOP;
283 __NOOP;
284 iexp = z.i.hi & 0x7ff00000; /* isolate float exponent */
285 return ((iexp != 0x7ff00000) && (iexp != 0));
286}
287
288int __isnormal ( long double x )
289{
290 return (FP_NORMAL == __fpclassify( x ));
291}
292
293/***********************************************************************
294 Function __isfinitef
295 Returns nonzero if and only if x is a finite (normal, subnormal,
296 or zero) float number and zero otherwise.
297***********************************************************************/
298
299int __isfinitef ( float x )
300{
301 hexsingle z;
302
303 z.fval = x;
304 __NOOP;
305 __NOOP;
306 __NOOP;
307 return ((z.lval & 0x7f800000) != 0x7f800000);
308}
309
310
311/***********************************************************************
312 Function __isfinited
313 Returns nonzero if and only if x is a finite (normal, subnormal,
314 or zero) double number and zero otherwise.
315***********************************************************************/
316int __isfinited ( double x )
317{
318 hexdouble z;
319
320 z.d = x;
321 __NOOP;
322 __NOOP;
323 __NOOP;
324 return ((z.i.hi & 0x7ff00000) != 0x7ff00000);
325}
326
327int __isfinite ( long double x )
328{
329 hexdbldbl u;
330
331 u.ld = x;
332 __NOOP;
333 __NOOP;
334 __NOOP;
335 return ((u.hh.hexhead.i.hi & 0x7ff00000) != 0x7ff00000);
336}
337
338/***********************************************************************
339 Function __isinff
340 Returns nonzero if and only if x is an infinite float number and zero
341 otherwise.
342***********************************************************************/
343
344int __isinff ( float x )
345{
346 hexsingle z;
347
348 z.fval = x;
349 __NOOP;
350 __NOOP;
351 __NOOP;
352 return (((z.lval&0x7f800000) == 0x7f800000) && ((z.lval&0x007fffff) == 0));
353}
354
355
356/***********************************************************************
357 Function __isinfd
358 Returns nonzero if and only if x is an infinite double number and zero
359 otherwise.
360***********************************************************************/
361int __isinfd ( double x )
362{
363 hexdouble z;
364
365 z.d = x;
366 __NOOP;
367 __NOOP;
368 __NOOP;
369 return (((z.i.hi&0x7ff00000) == 0x7ff00000) && (((z.i.hi&0x000fffff) | z.i.lo) == 0));
370}
371
372int __isinf ( long double x )
373{
374 hexdbldbl u;
375
376 u.ld = x;
377 __NOOP;
378 __NOOP;
379 __NOOP;
380 return (((u.hh.hexhead.i.hi&0x7ff00000) == 0x7ff00000) && (((u.hh.hexhead.i.hi&0x000fffff) | u.hh.hexhead.i.lo) == 0));
381}
382
383// Retain the following for legacy ABI compatability as it shipped in 10.2 and 10.3.
384int __isinfl ( long double x )
385{
386 return __isinf( x );
387}
388
389/***********************************************************************
390 Function __isnanf
391 Returns nonzero if and only if x is a float NaN and zero otherwise.
392***********************************************************************/
393
394int __isnanf ( float x )
395{
396 hexsingle z;
397
398 z.fval = x;
399 __NOOP;
400 __NOOP;
401 __NOOP;
402 return (((z.lval&0x7f800000) == 0x7f800000) && ((z.lval&0x007fffff) != 0));
403}
404
405
406/***********************************************************************
407 Function __isnand
408 Returns nonzero if and only if x is a double NaN and zero otherwise.
409***********************************************************************/
410int __isnand ( double x )
411{
412 hexdouble z;
413
414 z.d = x;
415 __NOOP;
416 __NOOP;
417 __NOOP;
418 return (((z.i.hi&0x7ff00000) == 0x7ff00000) && (((z.i.hi&0x000fffff) | z.i.lo) != 0));
419}
420
421int __isnan ( long double x )
422{
423 hexdbldbl u;
424
425 u.ld = x;
426 __NOOP;
427 __NOOP;
428 __NOOP;
429 return (((u.hh.hexhead.i.hi&0x7ff00000) == 0x7ff00000) && (((u.hh.hexhead.i.hi&0x000fffff) | u.hh.hexhead.i.lo) != 0));
430}
431
432/***********************************************************************
433 Function __signbitf
434 Returns nonzero if and only if the sign bit of the float number x is
435 set and zero otherwise.
436***********************************************************************/
437
438int __signbitf ( float x )
439{
440 hexsingle z;
441
442 z.fval = x;
443 __NOOP;
444 __NOOP;
445 __NOOP;
446 return (((int32_t)z.lval) < 0);
447}
448
449
450/***********************************************************************
451 Function __signbitd
452 Returns nonzero if and only if the sign bit of the double number x is
453 set and zero otherwise.
454***********************************************************************/
455
456int __signbitd ( double arg )
457{
458 hexdouble z;
459
460 z.d = arg;
461 __NOOP;
462 __NOOP;
463 __NOOP;
464 return (((int32_t)z.i.hi) < 0);
465}
466
467int __signbitl ( long double arg )
468{
469 hexdbldbl u;
470
471 u.ld = arg;
472 __NOOP;
473 __NOOP;
474 __NOOP;
475 return (((int32_t)u.hh.hexhead.i.hi) < 0);
476}
477
478// Retain the following for legacy ABI compatability as it shipped in 10.2 and 10.3.
479int __signbit ( long double arg )
480{
481 return __signbitl( arg );
482}
483
484float __nan ( void )
485{
486 static const hexsingle aQuietNAN = { 0x7fc00000 };
487 return ( aQuietNAN.fval );
488}
489
490/***********************************************************************
491 Function __inf
492 Returns the value of positive infinity for a double.
493***********************************************************************/
494
495float __inff ( void )
496{
497 static const hexsingle PosINF = { 0x7f800000 };
498 return ( PosINF.fval );
499}
500
501#else /* BUILDING_FOR_CARBONCORE_LEGACY */
502
503double __inf ( void )
504{
505 static const hexdouble PosINF = HEXDOUBLE(0x7ff00000, 0x00000000);
506 return ( PosINF.d );
507}
508
509#endif /* BUILDING_FOR_CARBONCORE_LEGACY */