Chromium Code Reviews
chromiumcodereview-hr@appspot.gserviceaccount.com (chromiumcodereview-hr) | Please choose your nickname with Settings | Help | Chromium Project | Gerrit Changes | Sign out
(354)

Side by Side Diff: fusl/src/math/logl.c

Issue 1714623002: [fusl] clang-format fusl (Closed) Base URL: git@github.com:domokit/mojo.git@master
Patch Set: headers too Created 4 years, 10 months ago
Use n/p to move between diff chunks; N/P to move between comments. Draft comments are only viewable by you.
Jump to:
View unified diff | Download patch
OLDNEW
1 /* origin: OpenBSD /usr/src/lib/libm/src/ld80/e_logl.c */ 1 /* origin: OpenBSD /usr/src/lib/libm/src/ld80/e_logl.c */
2 /* 2 /*
3 * Copyright (c) 2008 Stephen L. Moshier <steve@moshier.net> 3 * Copyright (c) 2008 Stephen L. Moshier <steve@moshier.net>
4 * 4 *
5 * Permission to use, copy, modify, and distribute this software for any 5 * Permission to use, copy, modify, and distribute this software for any
6 * purpose with or without fee is hereby granted, provided that the above 6 * purpose with or without fee is hereby granted, provided that the above
7 * copyright notice and this permission notice appear in all copies. 7 * copyright notice and this permission notice appear in all copies.
8 * 8 *
9 * THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES 9 * THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES
10 * WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF 10 * WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF
(...skipping 37 matching lines...) Expand 10 before | Expand all | Expand 10 after
48 * IEEE exp(+-10000) 100000 5.39e-20 2.34e-20 48 * IEEE exp(+-10000) 100000 5.39e-20 2.34e-20
49 * 49 *
50 * In the tests over the interval exp(+-10000), the logarithms 50 * In the tests over the interval exp(+-10000), the logarithms
51 * of the random arguments were uniformly distributed over 51 * of the random arguments were uniformly distributed over
52 * [-10000, +10000]. 52 * [-10000, +10000].
53 */ 53 */
54 54
55 #include "libm.h" 55 #include "libm.h"
56 56
57 #if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024 57 #if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
58 long double logl(long double x) 58 long double logl(long double x) {
59 { 59 return log(x);
60 » return log(x);
61 } 60 }
62 #elif LDBL_MANT_DIG == 64 && LDBL_MAX_EXP == 16384 61 #elif LDBL_MANT_DIG == 64 && LDBL_MAX_EXP == 16384
63 /* Coefficients for log(1+x) = x - x**2/2 + x**3 P(x)/Q(x) 62 /* Coefficients for log(1+x) = x - x**2/2 + x**3 P(x)/Q(x)
64 * 1/sqrt(2) <= x < sqrt(2) 63 * 1/sqrt(2) <= x < sqrt(2)
65 * Theoretical peak relative error = 2.32e-20 64 * Theoretical peak relative error = 2.32e-20
66 */ 65 */
67 static const long double P[] = { 66 static const long double P[] = {
68 4.5270000862445199635215E-5L, 67 4.5270000862445199635215E-5L, 4.9854102823193375972212E-1L,
69 4.9854102823193375972212E-1L, 68 6.5787325942061044846969E0L, 2.9911919328553073277375E1L,
70 6.5787325942061044846969E0L, 69 6.0949667980987787057556E1L, 5.7112963590585538103336E1L,
71 2.9911919328553073277375E1L, 70 2.0039553499201281259648E1L,
72 6.0949667980987787057556E1L,
73 5.7112963590585538103336E1L,
74 2.0039553499201281259648E1L,
75 }; 71 };
76 static const long double Q[] = { 72 static const long double Q[] = {
77 /* 1.0000000000000000000000E0,*/ 73 /* 1.0000000000000000000000E0,*/
78 1.5062909083469192043167E1L, 74 1.5062909083469192043167E1L, 8.3047565967967209469434E1L,
79 8.3047565967967209469434E1L, 75 2.2176239823732856465394E2L, 3.0909872225312059774938E2L,
80 2.2176239823732856465394E2L, 76 2.1642788614495947685003E2L, 6.0118660497603843919306E1L,
81 3.0909872225312059774938E2L,
82 2.1642788614495947685003E2L,
83 6.0118660497603843919306E1L,
84 }; 77 };
85 78
86 /* Coefficients for log(x) = z + z^3 P(z^2)/Q(z^2), 79 /* Coefficients for log(x) = z + z^3 P(z^2)/Q(z^2),
87 * where z = 2(x-1)/(x+1) 80 * where z = 2(x-1)/(x+1)
88 * 1/sqrt(2) <= x < sqrt(2) 81 * 1/sqrt(2) <= x < sqrt(2)
89 * Theoretical peak relative error = 6.16e-22 82 * Theoretical peak relative error = 6.16e-22
90 */ 83 */
91 static const long double R[4] = { 84 static const long double R[4] = {
92 1.9757429581415468984296E-3L, 85 1.9757429581415468984296E-3L, -7.1990767473014147232598E-1L,
93 -7.1990767473014147232598E-1L, 86 1.0777257190312272158094E1L, -3.5717684488096787370998E1L,
94 1.0777257190312272158094E1L,
95 -3.5717684488096787370998E1L,
96 }; 87 };
97 static const long double S[4] = { 88 static const long double S[4] = {
98 /* 1.00000000000000000000E0L,*/ 89 /* 1.00000000000000000000E0L,*/
99 -2.6201045551331104417768E1L, 90 -2.6201045551331104417768E1L, 1.9361891836232102174846E2L,
100 1.9361891836232102174846E2L, 91 -4.2861221385716144629696E2L,
101 -4.2861221385716144629696E2L,
102 }; 92 };
103 static const long double C1 = 6.9314575195312500000000E-1L; 93 static const long double C1 = 6.9314575195312500000000E-1L;
104 static const long double C2 = 1.4286068203094172321215E-6L; 94 static const long double C2 = 1.4286068203094172321215E-6L;
105 95
106 #define SQRTH 0.70710678118654752440L 96 #define SQRTH 0.70710678118654752440L
107 97
108 long double logl(long double x) 98 long double logl(long double x) {
109 { 99 long double y, z;
110 » long double y, z; 100 int e;
111 » int e;
112 101
113 » if (isnan(x)) 102 if (isnan(x))
114 » » return x; 103 return x;
115 » if (x == INFINITY) 104 if (x == INFINITY)
116 » » return x; 105 return x;
117 » if (x <= 0.0) { 106 if (x <= 0.0) {
118 » » if (x == 0.0) 107 if (x == 0.0)
119 » » » return -1/(x*x); /* -inf with divbyzero */ 108 return -1 / (x * x); /* -inf with divbyzero */
120 » » return 0/0.0f; /* nan with invalid */ 109 return 0 / 0.0f; /* nan with invalid */
121 » } 110 }
122 111
123 » /* separate mantissa from exponent */ 112 /* separate mantissa from exponent */
124 » /* Note, frexp is used so that denormal numbers 113 /* Note, frexp is used so that denormal numbers
125 » * will be handled properly. 114 * will be handled properly.
126 » */ 115 */
127 » x = frexpl(x, &e); 116 x = frexpl(x, &e);
128 117
129 » /* logarithm using log(x) = z + z**3 P(z)/Q(z), 118 /* logarithm using log(x) = z + z**3 P(z)/Q(z),
130 » * where z = 2(x-1)/(x+1) 119 * where z = 2(x-1)/(x+1)
131 » */ 120 */
132 » if (e > 2 || e < -2) { 121 if (e > 2 || e < -2) {
133 » » if (x < SQRTH) { /* 2(2x-1)/(2x+1) */ 122 if (x < SQRTH) { /* 2(2x-1)/(2x+1) */
134 » » » e -= 1; 123 e -= 1;
135 » » » z = x - 0.5; 124 z = x - 0.5;
136 » » » y = 0.5 * z + 0.5; 125 y = 0.5 * z + 0.5;
137 » » } else { /* 2 (x-1)/(x+1) */ 126 } else { /* 2 (x-1)/(x+1) */
138 » » » z = x - 0.5; 127 z = x - 0.5;
139 » » » z -= 0.5; 128 z -= 0.5;
140 » » » y = 0.5 * x + 0.5; 129 y = 0.5 * x + 0.5;
141 » » } 130 }
142 » » x = z / y; 131 x = z / y;
143 » » z = x*x; 132 z = x * x;
144 » » z = x * (z * __polevll(z, R, 3) / __p1evll(z, S, 3)); 133 z = x * (z * __polevll(z, R, 3) / __p1evll(z, S, 3));
145 » » z = z + e * C2; 134 z = z + e * C2;
146 » » z = z + x; 135 z = z + x;
147 » » z = z + e * C1; 136 z = z + e * C1;
148 » » return z; 137 return z;
149 » } 138 }
150 139
151 » /* logarithm using log(1+x) = x - .5x**2 + x**3 P(x)/Q(x) */ 140 /* logarithm using log(1+x) = x - .5x**2 + x**3 P(x)/Q(x) */
152 » if (x < SQRTH) { 141 if (x < SQRTH) {
153 » » e -= 1; 142 e -= 1;
154 » » x = 2.0*x - 1.0; 143 x = 2.0 * x - 1.0;
155 » } else { 144 } else {
156 » » x = x - 1.0; 145 x = x - 1.0;
157 » } 146 }
158 » z = x*x; 147 z = x * x;
159 » y = x * (z * __polevll(x, P, 6) / __p1evll(x, Q, 6)); 148 y = x * (z * __polevll(x, P, 6) / __p1evll(x, Q, 6));
160 » y = y + e * C2; 149 y = y + e * C2;
161 » z = y - 0.5*z; 150 z = y - 0.5 * z;
162 » /* Note, the sum of above terms does not exceed x/4, 151 /* Note, the sum of above terms does not exceed x/4,
163 » * so it contributes at most about 1/4 lsb to the error. 152 * so it contributes at most about 1/4 lsb to the error.
164 » */ 153 */
165 » z = z + x; 154 z = z + x;
166 » z = z + e * C1; /* This sum has an error of 1/2 lsb. */ 155 z = z + e * C1; /* This sum has an error of 1/2 lsb. */
167 » return z; 156 return z;
168 } 157 }
169 #elif LDBL_MANT_DIG == 113 && LDBL_MAX_EXP == 16384 158 #elif LDBL_MANT_DIG == 113 && LDBL_MAX_EXP == 16384
170 // TODO: broken implementation to make things compile 159 // TODO: broken implementation to make things compile
171 long double logl(long double x) 160 long double logl(long double x) {
172 { 161 return log(x);
173 » return log(x);
174 } 162 }
175 #endif 163 #endif
OLDNEW

Powered by Google App Engine
This is Rietveld 408576698