1
1
use super :: scalbnf;
2
2
3
- const HALF : [ f32 ; 2 ] = [ 0.5 , -0.5 ] ;
4
- const LN2_HI : f32 = 6.9314575195e-01 ; /* 0x3f317200 */
5
- const LN2_LO : f32 = 1.4286067653e-06 ; /* 0x35bfbe8e */
6
- const INV_LN2 : f32 = 1.4426950216e+00 ; /* 0x3fb8aa3b */
3
+ const HALF : [ f32 ; 2 ] = [ 0.5 , -0.5 ] ;
4
+ const LN2_HI : f32 = 6.9314575195e-01 ; /* 0x3f317200 */
5
+ const LN2_LO : f32 = 1.4286067653e-06 ; /* 0x35bfbe8e */
6
+ const INV_LN2 : f32 = 1.4426950216e+00 ; /* 0x3fb8aa3b */
7
7
/*
8
8
* Domain [-0.34568, 0.34568], range ~[-4.278e-9, 4.447e-9]:
9
9
* |x*(exp(x)+1)/(exp(x)-1) - p(x)| < 2**-27.74
10
10
*/
11
- const P1 : f32 = 1.6666625440e-1 ; /* 0xaaaa8f.0p-26 */
12
- const P2 : f32 = -2.7667332906e-3 ; /* -0xb55215.0p-32 */
11
+ const P1 : f32 = 1.6666625440e-1 ; /* 0xaaaa8f.0p-26 */
12
+ const P2 : f32 = -2.7667332906e-3 ; /* -0xb55215.0p-32 */
13
13
14
14
#[ inline]
15
15
pub fn expf ( mut x : f32 ) -> f32 {
16
- let x1p127 = f32:: from_bits ( 0x7f000000 ) ; // 0x1p127f === 2 ^ 127
17
- let x1p_126 = f32:: from_bits ( 0x800000 ) ; // 0x1p-126f === 2 ^ -126 /*original 0x1p-149f ??????????? */
18
-
16
+ let x1p127 = f32:: from_bits ( 0x7f000000 ) ; // 0x1p127f === 2 ^ 127
17
+ let x1p_126 = f32:: from_bits ( 0x800000 ) ; // 0x1p-126f === 2 ^ -126 /*original 0x1p-149f ??????????? */
18
+
19
19
let mut hx = x. to_bits ( ) ;
20
- let sign = ( hx >> 31 ) as i32 ; /* sign bit of x */
21
- let signb : bool = sign != 0 ;
22
- hx &= 0x7fffffff ; /* high word of |x| */
23
-
20
+ let sign = ( hx >> 31 ) as i32 ; /* sign bit of x */
21
+ let signb: bool = sign != 0 ;
22
+ hx &= 0x7fffffff ; /* high word of |x| */
23
+
24
24
/* special cases */
25
- if hx >= 0x42aeac50 { /* if |x| >= -87.33655f or NaN */
26
- if hx > 0x7f800000 { /* NaN */
25
+ if hx >= 0x42aeac50 {
26
+ /* if |x| >= -87.33655f or NaN */
27
+ if hx > 0x7f800000 {
28
+ /* NaN */
27
29
return x;
28
30
}
29
- if ( hx >= 0x42b17218 ) && ( !signb) { /* x >= 88.722839f */
31
+ if ( hx >= 0x42b17218 ) && ( !signb) {
32
+ /* x >= 88.722839f */
30
33
/* overflow */
31
34
x *= x1p127;
32
35
return x;
33
36
}
34
37
if signb {
35
38
/* underflow */
36
- force_eval ! ( -x1p_126/x) ;
37
- if hx >= 0x42cff1b5 { /* x <= -103.972084f */
38
- return 0.
39
+ force_eval ! ( -x1p_126 / x) ;
40
+ if hx >= 0x42cff1b5 {
41
+ /* x <= -103.972084f */
42
+ return 0. ;
39
43
}
40
44
}
41
45
}
42
-
46
+
43
47
/* argument reduction */
44
- let k : i32 ;
45
- let hi : f32 ;
46
- let lo : f32 ;
47
- if hx > 0x3eb17218 { /* if |x| > 0.5 ln2 */
48
- if hx > 0x3f851592 { /* if |x| > 1.5 ln2 */
49
- k = ( INV_LN2 * x + HALF [ sign as usize ] ) as i32 ;
48
+ let k: i32 ;
49
+ let hi: f32 ;
50
+ let lo: f32 ;
51
+ if hx > 0x3eb17218 {
52
+ /* if |x| > 0.5 ln2 */
53
+ if hx > 0x3f851592 {
54
+ /* if |x| > 1.5 ln2 */
55
+ k = ( INV_LN2 * x + HALF [ sign as usize ] ) as i32 ;
50
56
} else {
51
57
k = 1 - sign - sign;
52
58
}
53
59
let kf = k as f32 ;
54
- hi = x - kf* LN2_HI ; /* k*ln2hi is exact here */
55
- lo = kf* LN2_LO ;
60
+ hi = x - kf * LN2_HI ; /* k*ln2hi is exact here */
61
+ lo = kf * LN2_LO ;
56
62
x = hi - lo;
57
- } else if hx > 0x39000000 { /* |x| > 2**-14 */
63
+ } else if hx > 0x39000000 {
64
+ /* |x| > 2**-14 */
58
65
k = 0 ;
59
66
hi = x;
60
67
lo = 0. ;
@@ -63,11 +70,11 @@ pub fn expf(mut x: f32) -> f32 {
63
70
force_eval ! ( x1p127 + x) ;
64
71
return 1. + x;
65
72
}
66
-
73
+
67
74
/* x is now in primary range */
68
- let xx = x* x;
69
- let c = x - xx* ( P1 +xx * P2 ) ;
70
- let y = 1. + ( x* c/ ( 2. - c) - lo + hi) ;
75
+ let xx = x * x;
76
+ let c = x - xx * ( P1 + xx * P2 ) ;
77
+ let y = 1. + ( x * c / ( 2. - c) - lo + hi) ;
71
78
if k == 0 {
72
79
y
73
80
} else {
0 commit comments