Home | History | Annotate | Line # | Download | only in src
s_logb.c revision 1.4
      1  1.1  jtc /* @(#)s_logb.c 5.1 93/09/24 */
      2  1.1  jtc /*
      3  1.1  jtc  * ====================================================
      4  1.1  jtc  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
      5  1.1  jtc  *
      6  1.1  jtc  * Developed at SunPro, a Sun Microsystems, Inc. business.
      7  1.1  jtc  * Permission to use, copy, modify, and distribute this
      8  1.1  jtc  * software is freely granted, provided that this notice
      9  1.1  jtc  * is preserved.
     10  1.1  jtc  * ====================================================
     11  1.1  jtc  */
     12  1.3  jtc 
     13  1.3  jtc #ifndef lint
     14  1.4  jtc static char rcsid[] = "$Id: s_logb.c,v 1.4 1994/03/03 17:04:41 jtc Exp $";
     15  1.3  jtc #endif
     16  1.1  jtc 
     17  1.1  jtc /*
     18  1.1  jtc  * double logb(x)
     19  1.1  jtc  * IEEE 754 logb. Included to pass IEEE test suite. Not recommend.
     20  1.1  jtc  * Use ilogb instead.
     21  1.1  jtc  */
     22  1.1  jtc 
     23  1.2  jtc #include <math.h>
     24  1.4  jtc #include <machine/endian.h>
     25  1.1  jtc 
     26  1.4  jtc #if BYTE_ORDER == LITTLE_ENDIAN
     27  1.4  jtc #define n0	1
     28  1.1  jtc #else
     29  1.4  jtc #define n0	0
     30  1.1  jtc #endif
     31  1.1  jtc 
     32  1.1  jtc #ifdef __STDC__
     33  1.1  jtc 	double logb(double x)
     34  1.1  jtc #else
     35  1.1  jtc 	double logb(x)
     36  1.1  jtc 	double x;
     37  1.1  jtc #endif
     38  1.1  jtc {
     39  1.4  jtc 	int lx,ix;
     40  1.1  jtc 	ix = (*(n0+(int*)&x))&0x7fffffff;	/* high |x| */
     41  1.1  jtc 	lx = *(1-n0+(int*)&x);			/* low x */
     42  1.1  jtc 	if((ix|lx)==0) return -1.0/fabs(x);
     43  1.1  jtc 	if(ix>=0x7ff00000) return x*x;
     44  1.1  jtc 	if((ix>>=20)==0) 			/* IEEE 754 logb */
     45  1.1  jtc 		return -1022.0;
     46  1.1  jtc 	else
     47  1.1  jtc 		return (double) (ix-1023);
     48  1.1  jtc }
     49