1 /* e_logf.c -- float version of e_log.c.
2 * Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com.
3 * adapted for log2 by Ulrich Drepper <drepper@cygnus.com>
7 * ====================================================
8 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
10 * Developed at SunPro, a Sun Microsystems, Inc. business.
11 * Permission to use, copy, modify, and distribute this
12 * software is freely granted, provided that this notice
14 * ====================================================
19 #include "math_private.h"
26 two25
= 3.355443200e+07, /* 0x4c000000 */
27 Lg1
= 6.6666668653e-01, /* 3F2AAAAB */
28 Lg2
= 4.0000000596e-01, /* 3ECCCCCD */
29 Lg3
= 2.8571429849e-01, /* 3E924925 */
30 Lg4
= 2.2222198546e-01, /* 3E638E29 */
31 Lg5
= 1.8183572590e-01, /* 3E3A3325 */
32 Lg6
= 1.5313838422e-01, /* 3E1CD04F */
33 Lg7
= 1.4798198640e-01; /* 3E178897 */
36 static const float zero
= 0.0;
38 static float zero
= 0.0;
42 float __log2f(float x
)
48 float hfsq
,f
,s
,z
,R
,w
,t1
,t2
,dk
;
54 if (ix
< 0x00800000) { /* x < 2**-126 */
55 if ((ix
&0x7fffffff)==0)
56 return -two25
/zero
; /* log(+-0)=-inf */
57 if (ix
<0) return (x
-x
)/zero
; /* log(-#) = NaN */
58 k
-= 25; x
*= two25
; /* subnormal number, scale up x */
61 if (ix
>= 0x7f800000) return x
+x
;
64 i
= (ix
+(0x95f64<<3))&0x800000;
65 SET_FLOAT_WORD(x
,ix
|(i
^0x3f800000)); /* normalize x or x/2 */
69 if((0x007fffff&(15+ix
))<16) { /* |f| < 2**-20 */
70 if(f
==zero
) return dk
;
71 R
= f
*f
*((float)0.5-(float)0.33333333333333333*f
);
79 t1
= w
*(Lg2
+w
*(Lg4
+w
*Lg6
));
80 t2
= z
*(Lg1
+w
*(Lg3
+w
*(Lg5
+w
*Lg7
)));
85 return dk
-((hfsq
-(s
*(hfsq
+R
)))-f
);
87 return dk
-((s
*(f
-R
))-f
);
90 weak_alias (__log2f
, log2f
)