1 /*
2 * CDDL HEADER START
3 *
4 * The contents of this file are subject to the terms of the
5 * Common Development and Distribution License (the "License").
6 * You may not use this file except in compliance with the License.
7 *
8 * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE
9 * or http://www.opensolaris.org/os/licensing.
10 * See the License for the specific language governing permissions
11 * and limitations under the License.
12 *
13 * When distributing Covered Code, include this CDDL HEADER in each
14 * file and include the License file at usr/src/OPENSOLARIS.LICENSE.
15 * If applicable, add the following below this CDDL HEADER, with the
16 * fields enclosed by brackets "[]" replaced with your own identifying
17 * information: Portions Copyright [yyyy] [name of copyright owner]
18 *
19 * CDDL HEADER END
20 */
21 /*
22 * Copyright 2011 Nexenta Systems, Inc. All rights reserved.
23 */
24 /*
25 * Copyright 2006 Sun Microsystems, Inc. All rights reserved.
26 * Use is subject to license terms.
27 */
28
29 /* INDENT OFF */
30 /*
31 * exp10(x)
32 * Code by K.C. Ng for SUN 4.0 libm.
33 * Method :
34 * n = nint(x*(log10/log2));
35 * exp10(x) = 10**x = exp(x*ln(10)) = exp(n*ln2+(x*ln10-n*ln2))
36 * = 2**n*exp(ln10*(x-n*log2/log10)))
37 * If x is an integer < 23 then use repeat multiplication. For
38 * 10**22 is the largest representable integer.
39 */
40 /* INDENT ON */
41
42 #include "libm.h"
43
44 static const double C[] = {
45 3.3219280948736234787, /* log(10)/log(2) */
46 2.3025850929940456840, /* log(10) */
47 3.0102999565860955045E-1, /* log(2)/log(10) high */
48 5.3716447674669983622E-12, /* log(2)/log(10) low */
49 0.0,
50 0.5,
51 1.0,
52 10.0,
53 1.0e300,
54 1.0e-300,
55 };
56
57 #define lg10 C[0]
58 #define ln10 C[1]
59 #define logt2hi C[2]
60 #define logt2lo C[3]
61 #define zero C[4]
62 #define half C[5]
63 #define one C[6]
64 #define ten C[7]
65 #define huge C[8]
66 #define tiny C[9]
67
68 double
69 exp10(double x) {
70 double t, pt;
71 int ix, hx, k;
72
73 ix = ((int *)&x)[HIWORD];
74 hx = ix & ~0x80000000;
75
76 if (hx >= 0x4074a000) { /* |x| >= 330 or x is nan */
77 if (hx >= 0x7ff00000) { /* x is inf or nan */
78 if (ix == 0xfff00000 && ((int *)&x)[LOWORD] == 0)
79 return (zero);
80 return (x * x);
81 }
82 t = (ix < 0)? tiny : huge;
83 return (t * t);
84 }
85
86 if (hx < 0x3c000000)
87 return (one + x);
88
89 k = (int)x;
90 if (0 <= k && k < 23 && (double)k == x) {
91 /* x is a small positive integer */
92 t = one;
93 pt = ten;
94 if (k & 1)
95 t = ten;
96 k >>= 1;
97 while (k) {
98 pt *= pt;
99 if (k & 1)
100 t *= pt;
101 k >>= 1;
102 }
103 return (t);
104 }
105 t = x * lg10;
106 k = (int)((ix < 0)? t - half : t + half);
107 return (scalbn(exp(ln10 * ((x - k * logt2hi) - k * logt2lo)), k));
108 }
|
1 /*
2 * CDDL HEADER START
3 *
4 * The contents of this file are subject to the terms of the
5 * Common Development and Distribution License (the "License").
6 * You may not use this file except in compliance with the License.
7 *
8 * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE
9 * or http://www.opensolaris.org/os/licensing.
10 * See the License for the specific language governing permissions
11 * and limitations under the License.
12 *
13 * When distributing Covered Code, include this CDDL HEADER in each
14 * file and include the License file at usr/src/OPENSOLARIS.LICENSE.
15 * If applicable, add the following below this CDDL HEADER, with the
16 * fields enclosed by brackets "[]" replaced with your own identifying
17 * information: Portions Copyright [yyyy] [name of copyright owner]
18 *
19 * CDDL HEADER END
20 */
21
22 /*
23 * Copyright 2011 Nexenta Systems, Inc. All rights reserved.
24 */
25
26 /*
27 * Copyright 2006 Sun Microsystems, Inc. All rights reserved.
28 * Use is subject to license terms.
29 */
30
31
32 /*
33 * exp10(x)
34 * Code by K.C. Ng for SUN 4.0 libm.
35 * Method :
36 * n = nint(x*(log10/log2));
37 * exp10(x) = 10**x = exp(x*ln(10)) = exp(n*ln2+(x*ln10-n*ln2))
38 * = 2**n*exp(ln10*(x-n*log2/log10)))
39 * If x is an integer < 23 then use repeat multiplication. For
40 * 10**22 is the largest representable integer.
41 */
42
43 #include "libm.h"
44
45 static const double C[] = {
46 3.3219280948736234787, /* log(10)/log(2) */
47 2.3025850929940456840, /* log(10) */
48 3.0102999565860955045E-1, /* log(2)/log(10) high */
49 5.3716447674669983622E-12, /* log(2)/log(10) low */
50 0.0,
51 0.5,
52 1.0,
53 10.0,
54 1.0e300,
55 1.0e-300,
56 };
57
58 #define lg10 C[0]
59 #define ln10 C[1]
60 #define logt2hi C[2]
61 #define logt2lo C[3]
62 #define zero C[4]
63 #define half C[5]
64 #define one C[6]
65 #define ten C[7]
66 #define huge C[8]
67 #define tiny C[9]
68
69 double
70 exp10(double x)
71 {
72 double t, pt;
73 int ix, hx, k;
74
75 ix = ((int *)&x)[HIWORD];
76 hx = ix & ~0x80000000;
77
78 if (hx >= 0x4074a000) { /* |x| >= 330 or x is nan */
79 if (hx >= 0x7ff00000) { /* x is inf or nan */
80 if (ix == 0xfff00000 && ((int *)&x)[LOWORD] == 0)
81 return (zero);
82
83 return (x * x);
84 }
85
86 t = (ix < 0) ? tiny : huge;
87 return (t * t);
88 }
89
90 if (hx < 0x3c000000)
91 return (one + x);
92
93 k = (int)x;
94
95 if (0 <= k && k < 23 && (double)k == x) {
96 /* x is a small positive integer */
97 t = one;
98 pt = ten;
99
100 if (k & 1)
101 t = ten;
102
103 k >>= 1;
104
105 while (k) {
106 pt *= pt;
107
108 if (k & 1)
109 t *= pt;
110
111 k >>= 1;
112 }
113
114 return (t);
115 }
116
117 t = x * lg10;
118 k = (int)((ix < 0) ? t - half : t + half);
119 return (scalbn(exp(ln10 * ((x - k * logt2hi) - k * logt2lo)), k));
120 }
|