OLD | NEW |
1 // Copyright 2010 the V8 project authors. All rights reserved. | 1 // Copyright 2010 the V8 project authors. All rights reserved. |
2 // Redistribution and use in source and binary forms, with or without | 2 // Redistribution and use in source and binary forms, with or without |
3 // modification, are permitted provided that the following conditions are | 3 // modification, are permitted provided that the following conditions are |
4 // met: | 4 // met: |
5 // | 5 // |
6 // * Redistributions of source code must retain the above copyright | 6 // * Redistributions of source code must retain the above copyright |
7 // notice, this list of conditions and the following disclaimer. | 7 // notice, this list of conditions and the following disclaimer. |
8 // * Redistributions in binary form must reproduce the above | 8 // * Redistributions in binary form must reproduce the above |
9 // copyright notice, this list of conditions and the following | 9 // copyright notice, this list of conditions and the following |
10 // disclaimer in the documentation and/or other materials provided | 10 // disclaimer in the documentation and/or other materials provided |
(...skipping 13 matching lines...) Expand all Loading... |
24 // THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT | 24 // THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT |
25 // (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE | 25 // (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE |
26 // OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. | 26 // OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. |
27 | 27 |
28 #include <stdarg.h> | 28 #include <stdarg.h> |
29 #include <limits.h> | 29 #include <limits.h> |
30 | 30 |
31 #include "v8.h" | 31 #include "v8.h" |
32 | 32 |
33 #include "strtod.h" | 33 #include "strtod.h" |
34 #include "cached-powers.h" | 34 // #include "cached-powers.h" |
35 #include "double.h" | |
36 | 35 |
37 namespace v8 { | 36 namespace v8 { |
38 namespace internal { | 37 namespace internal { |
39 | 38 |
40 // 2^53 = 9007199254740992. | 39 // 2^53 = 9007199254740992. |
41 // Any integer with at most 15 decimal digits will hence fit into a double | 40 // Any integer with at most 15 decimal digits will hence fit into a double |
42 // (which has a 53bit significand) without loss of precision. | 41 // (which has a 53bit significand) without loss of precision. |
43 static const int kMaxExactDoubleIntegerDecimalDigits = 15; | 42 static const int kMaxExactDoubleIntegerDecimalDigits = 15; |
44 // 2^64 = 18446744073709551616 > 10^19 | 43 // 2^64 = 18446744073709551616 |
| 44 // Any integer with at most 19 digits will hence fit into a 64bit datatype. |
45 static const int kMaxUint64DecimalDigits = 19; | 45 static const int kMaxUint64DecimalDigits = 19; |
46 | |
47 // Max double: 1.7976931348623157 x 10^308 | 46 // Max double: 1.7976931348623157 x 10^308 |
48 // Min non-zero double: 4.9406564584124654 x 10^-324 | 47 // Min non-zero double: 4.9406564584124654 x 10^-324 |
49 // Any x >= 10^309 is interpreted as +infinity. | 48 // Any x >= 10^309 is interpreted as +infinity. |
50 // Any x <= 10^-324 is interpreted as 0. | 49 // Any x <= 10^-324 is interpreted as 0. |
51 // Note that 2.5e-324 (despite being smaller than the min double) will be read | 50 // Note that 2.5e-324 (despite being smaller than the min double) will be read |
52 // as non-zero (equal to the min non-zero double). | 51 // as non-zero (equal to the min non-zero double). |
53 static const int kMaxDecimalPower = 309; | 52 static const int kMaxDecimalPower = 309; |
54 static const int kMinDecimalPower = -324; | 53 static const int kMinDecimalPower = -324; |
55 | 54 |
56 // 2^64 = 18446744073709551616 | |
57 static const uint64_t kMaxUint64 = V8_2PART_UINT64_C(0xFFFFFFFF, FFFFFFFF); | |
58 | |
59 | |
60 static const double exact_powers_of_ten[] = { | 55 static const double exact_powers_of_ten[] = { |
61 1.0, // 10^0 | 56 1.0, // 10^0 |
62 10.0, | 57 10.0, |
63 100.0, | 58 100.0, |
64 1000.0, | 59 1000.0, |
65 10000.0, | 60 10000.0, |
66 100000.0, | 61 100000.0, |
67 1000000.0, | 62 1000000.0, |
68 10000000.0, | 63 10000000.0, |
69 100000000.0, | 64 100000000.0, |
(...skipping 65 matching lines...) Expand 10 before | Expand all | Expand 10 after Loading... |
135 static Vector<const char> TrimTrailingZeros(Vector<const char> buffer) { | 130 static Vector<const char> TrimTrailingZeros(Vector<const char> buffer) { |
136 for (int i = buffer.length() - 1; i >= 0; --i) { | 131 for (int i = buffer.length() - 1; i >= 0; --i) { |
137 if (buffer[i] != '0') { | 132 if (buffer[i] != '0') { |
138 return Vector<const char>(buffer.start(), i + 1); | 133 return Vector<const char>(buffer.start(), i + 1); |
139 } | 134 } |
140 } | 135 } |
141 return Vector<const char>(buffer.start(), 0); | 136 return Vector<const char>(buffer.start(), 0); |
142 } | 137 } |
143 | 138 |
144 | 139 |
145 // Reads digits from the buffer and converts them to a uint64. | 140 uint64_t ReadUint64(Vector<const char> buffer) { |
146 // Reads in as many digits as fit into a uint64. | 141 ASSERT(buffer.length() <= kMaxUint64DecimalDigits); |
147 // When the string starts with "1844674407370955161" no further digit is read. | |
148 // Since 2^64 = 18446744073709551616 it would still be possible read another | |
149 // digit if it was less or equal than 6, but this would complicate the code. | |
150 static uint64_t ReadUint64(Vector<const char> buffer, | |
151 int* number_of_read_digits) { | |
152 uint64_t result = 0; | 142 uint64_t result = 0; |
153 int i = 0; | 143 for (int i = 0; i < buffer.length(); ++i) { |
154 while (i < buffer.length() && result <= (kMaxUint64 / 10 - 1)) { | 144 int digit = buffer[i] - '0'; |
155 int digit = buffer[i++] - '0'; | |
156 ASSERT(0 <= digit && digit <= 9); | 145 ASSERT(0 <= digit && digit <= 9); |
157 result = 10 * result + digit; | 146 result = 10 * result + digit; |
158 } | 147 } |
159 *number_of_read_digits = i; | |
160 return result; | 148 return result; |
161 } | 149 } |
162 | 150 |
163 | 151 |
164 // Reads a DiyFp from the buffer. | |
165 // The returned DiyFp is not necessarily normalized. | |
166 // If remaining_decimals is zero then the returned DiyFp is accurate. | |
167 // Otherwise it has been rounded and has error of at most 1/2 ulp. | |
168 static void ReadDiyFp(Vector<const char> buffer, | |
169 DiyFp* result, | |
170 int* remaining_decimals) { | |
171 int read_digits; | |
172 uint64_t significand = ReadUint64(buffer, &read_digits); | |
173 if (buffer.length() == read_digits) { | |
174 *result = DiyFp(significand, 0); | |
175 *remaining_decimals = 0; | |
176 } else { | |
177 // Round the significand. | |
178 if (buffer[read_digits] >= '5') { | |
179 significand++; | |
180 } | |
181 // Compute the binary exponent. | |
182 int exponent = 0; | |
183 *result = DiyFp(significand, exponent); | |
184 *remaining_decimals = buffer.length() - read_digits; | |
185 } | |
186 } | |
187 | |
188 | |
189 static bool DoubleStrtod(Vector<const char> trimmed, | 152 static bool DoubleStrtod(Vector<const char> trimmed, |
190 int exponent, | 153 int exponent, |
191 double* result) { | 154 double* result) { |
192 #if (defined(V8_TARGET_ARCH_IA32) || defined(USE_SIMULATOR)) && !defined(WIN32) | 155 #if (defined(V8_TARGET_ARCH_IA32) || defined(USE_SIMULATOR)) && !defined(WIN32) |
193 // On x86 the floating-point stack can be 64 or 80 bits wide. If it is | 156 // On x86 the floating-point stack can be 64 or 80 bits wide. If it is |
194 // 80 bits wide (as is the case on Linux) then double-rounding occurs and the | 157 // 80 bits wide (as is the case on Linux) then double-rounding occurs and the |
195 // result is not accurate. | 158 // result is not accurate. |
196 // We know that Windows32 uses 64 bits and is therefore accurate. | 159 // We know that Windows32 uses 64 bits and is therefore accurate. |
197 // Note that the ARM simulator is compiled for 32bits. It therefore exhibits | 160 // Note that the ARM simulator is compiled for 32bits. It therefore exhibits |
198 // the same problem. | 161 // the same problem. |
199 return false; | 162 return false; |
200 #endif | 163 #endif |
201 if (trimmed.length() <= kMaxExactDoubleIntegerDecimalDigits) { | 164 if (trimmed.length() <= kMaxExactDoubleIntegerDecimalDigits) { |
202 int read_digits; | |
203 // The trimmed input fits into a double. | 165 // The trimmed input fits into a double. |
204 // If the 10^exponent (resp. 10^-exponent) fits into a double too then we | 166 // If the 10^exponent (resp. 10^-exponent) fits into a double too then we |
205 // can compute the result-double simply by multiplying (resp. dividing) the | 167 // can compute the result-double simply by multiplying (resp. dividing) the |
206 // two numbers. | 168 // two numbers. |
207 // This is possible because IEEE guarantees that floating-point operations | 169 // This is possible because IEEE guarantees that floating-point operations |
208 // return the best possible approximation. | 170 // return the best possible approximation. |
209 if (exponent < 0 && -exponent < kExactPowersOfTenSize) { | 171 if (exponent < 0 && -exponent < kExactPowersOfTenSize) { |
210 // 10^-exponent fits into a double. | 172 // 10^-exponent fits into a double. |
211 *result = static_cast<double>(ReadUint64(trimmed, &read_digits)); | 173 *result = static_cast<double>(ReadUint64(trimmed)); |
212 ASSERT(read_digits == trimmed.length()); | |
213 *result /= exact_powers_of_ten[-exponent]; | 174 *result /= exact_powers_of_ten[-exponent]; |
214 return true; | 175 return true; |
215 } | 176 } |
216 if (0 <= exponent && exponent < kExactPowersOfTenSize) { | 177 if (0 <= exponent && exponent < kExactPowersOfTenSize) { |
217 // 10^exponent fits into a double. | 178 // 10^exponent fits into a double. |
218 *result = static_cast<double>(ReadUint64(trimmed, &read_digits)); | 179 *result = static_cast<double>(ReadUint64(trimmed)); |
219 ASSERT(read_digits == trimmed.length()); | |
220 *result *= exact_powers_of_ten[exponent]; | 180 *result *= exact_powers_of_ten[exponent]; |
221 return true; | 181 return true; |
222 } | 182 } |
223 int remaining_digits = | 183 int remaining_digits = |
224 kMaxExactDoubleIntegerDecimalDigits - trimmed.length(); | 184 kMaxExactDoubleIntegerDecimalDigits - trimmed.length(); |
225 if ((0 <= exponent) && | 185 if ((0 <= exponent) && |
226 (exponent - remaining_digits < kExactPowersOfTenSize)) { | 186 (exponent - remaining_digits < kExactPowersOfTenSize)) { |
227 // The trimmed string was short and we can multiply it with | 187 // The trimmed string was short and we can multiply it with |
228 // 10^remaining_digits. As a result the remaining exponent now fits | 188 // 10^remaining_digits. As a result the remaining exponent now fits |
229 // into a double too. | 189 // into a double too. |
230 *result = static_cast<double>(ReadUint64(trimmed, &read_digits)); | 190 *result = static_cast<double>(ReadUint64(trimmed)); |
231 ASSERT(read_digits == trimmed.length()); | |
232 *result *= exact_powers_of_ten[remaining_digits]; | 191 *result *= exact_powers_of_ten[remaining_digits]; |
233 *result *= exact_powers_of_ten[exponent - remaining_digits]; | 192 *result *= exact_powers_of_ten[exponent - remaining_digits]; |
234 return true; | 193 return true; |
235 } | 194 } |
236 } | 195 } |
237 return false; | 196 return false; |
238 } | 197 } |
239 | 198 |
240 | 199 |
241 // Returns 10^exponent as an exact DiyFp. | |
242 // The given exponent must be in the range [1; kDecimalExponentDistance[. | |
243 static DiyFp AdjustmentPowerOfTen(int exponent) { | |
244 ASSERT(0 < exponent); | |
245 ASSERT(exponent < PowersOfTenCache::kDecimalExponentDistance); | |
246 // Simply hardcode the remaining powers for the given decimal exponent | |
247 // distance. | |
248 ASSERT(PowersOfTenCache::kDecimalExponentDistance == 8); | |
249 switch (exponent) { | |
250 case 1: return DiyFp(V8_2PART_UINT64_C(0xa0000000, 00000000), -60); | |
251 case 2: return DiyFp(V8_2PART_UINT64_C(0xc8000000, 00000000), -57); | |
252 case 3: return DiyFp(V8_2PART_UINT64_C(0xfa000000, 00000000), -54); | |
253 case 4: return DiyFp(V8_2PART_UINT64_C(0x9c400000, 00000000), -50); | |
254 case 5: return DiyFp(V8_2PART_UINT64_C(0xc3500000, 00000000), -47); | |
255 case 6: return DiyFp(V8_2PART_UINT64_C(0xf4240000, 00000000), -44); | |
256 case 7: return DiyFp(V8_2PART_UINT64_C(0x98968000, 00000000), -40); | |
257 default: | |
258 UNREACHABLE(); | |
259 return DiyFp(0, 0); | |
260 } | |
261 } | |
262 | |
263 | |
264 // If the function returns true then the result is the correct double. | |
265 // Otherwise it is either the correct double or the double that is just below | |
266 // the correct double. | |
267 static bool DiyFpStrtod(Vector<const char> buffer, | |
268 int exponent, | |
269 double* result) { | |
270 DiyFp input; | |
271 int remaining_decimals; | |
272 ReadDiyFp(buffer, &input, &remaining_decimals); | |
273 // Since we may have dropped some digits the input is not accurate. | |
274 // If remaining_decimals is different than 0 than the error is at most | |
275 // .5 ulp (unit in the last place). | |
276 // We don't want to deal with fractions and therefore keep a common | |
277 // denominator. | |
278 const int kDenominatorLog = 3; | |
279 const int kDenominator = 1 << kDenominatorLog; | |
280 // Move the remaining decimals into the exponent. | |
281 exponent += remaining_decimals; | |
282 int error = (remaining_decimals == 0 ? 0 : kDenominator / 2); | |
283 | |
284 int old_e = input.e(); | |
285 input.Normalize(); | |
286 error <<= old_e - input.e(); | |
287 | |
288 ASSERT(exponent <= PowersOfTenCache::kMaxDecimalExponent); | |
289 if (exponent < PowersOfTenCache::kMinDecimalExponent) { | |
290 *result = 0.0; | |
291 return true; | |
292 } | |
293 DiyFp cached_power; | |
294 int cached_decimal_exponent; | |
295 PowersOfTenCache::GetCachedPowerForDecimalExponent(exponent, | |
296 &cached_power, | |
297 &cached_decimal_exponent); | |
298 | |
299 if (cached_decimal_exponent != exponent) { | |
300 int adjustment_exponent = exponent - cached_decimal_exponent; | |
301 DiyFp adjustment_power = AdjustmentPowerOfTen(adjustment_exponent); | |
302 input.Multiply(adjustment_power); | |
303 if (kMaxUint64DecimalDigits - buffer.length() >= adjustment_exponent) { | |
304 // The product of input with the adjustment power fits into a 64 bit | |
305 // integer. | |
306 ASSERT(DiyFp::kSignificandSize == 64); | |
307 } else { | |
308 // The adjustment power is exact. There is hence only an error of 0.5. | |
309 error += kDenominator / 2; | |
310 } | |
311 } | |
312 | |
313 input.Multiply(cached_power); | |
314 // The error introduced by a multiplication of a*b equals | |
315 // error_a + error_b + error_a*error_b/2^64 + 0.5 | |
316 // Substituting a with 'input' and b with 'cached_power' we have | |
317 // error_b = 0.5 (all cached powers have an error of less than 0.5 ulp), | |
318 // error_ab = 0 or 1 / kDenominator > error_a*error_b/ 2^64 | |
319 int error_b = kDenominator / 2; | |
320 int error_ab = (error == 0 ? 0 : 1); // We round up to 1. | |
321 int fixed_error = kDenominator / 2; | |
322 error += error_b + error_ab + fixed_error; | |
323 | |
324 old_e = input.e(); | |
325 input.Normalize(); | |
326 error <<= old_e - input.e(); | |
327 | |
328 // See if the double's significand changes if we add/subtract the error. | |
329 int order_of_magnitude = DiyFp::kSignificandSize + input.e(); | |
330 int effective_significand_size = | |
331 Double::SignificandSizeForOrderOfMagnitude(order_of_magnitude); | |
332 int precision_digits_count = | |
333 DiyFp::kSignificandSize - effective_significand_size; | |
334 if (precision_digits_count + kDenominatorLog >= DiyFp::kSignificandSize) { | |
335 // This can only happen for very small denormals. In this case the | |
336 // half-way multiplied by the denominator exceeds the range of an uint64. | |
337 // Simply shift everything to the right. | |
338 int shift_amount = (precision_digits_count + kDenominatorLog) - | |
339 DiyFp::kSignificandSize + 1; | |
340 input.set_f(input.f() >> shift_amount); | |
341 input.set_e(input.e() + shift_amount); | |
342 // We add 1 for the lost precision of error, and kDenominator for | |
343 // the lost precision of input.f(). | |
344 error = (error >> shift_amount) + 1 + kDenominator; | |
345 precision_digits_count -= shift_amount; | |
346 } | |
347 // We use uint64_ts now. This only works if the DiyFp uses uint64_ts too. | |
348 ASSERT(DiyFp::kSignificandSize == 64); | |
349 ASSERT(precision_digits_count < 64); | |
350 uint64_t one64 = 1; | |
351 uint64_t precision_bits_mask = (one64 << precision_digits_count) - 1; | |
352 uint64_t precision_bits = input.f() & precision_bits_mask; | |
353 uint64_t half_way = one64 << (precision_digits_count - 1); | |
354 precision_bits *= kDenominator; | |
355 half_way *= kDenominator; | |
356 // If the last_bits are too close to the half-way case than we are too | |
357 // inaccurate and round down. In this case we return false so that we can | |
358 // fall back to a more precise algorithm. | |
359 uint64_t significand = input.f(); | |
360 if (precision_bits >= half_way + error) { | |
361 significand = (significand >> precision_digits_count) + 1; | |
362 exponent = input.e() + precision_digits_count; | |
363 } else { | |
364 significand = (significand >> precision_digits_count); | |
365 exponent = input.e() + precision_digits_count; | |
366 } | |
367 Double d = Double(significand, exponent); | |
368 *result = d.value(); | |
369 if (half_way - error < precision_bits && precision_bits < half_way + error) { | |
370 // Too imprecise. The caller will have to fall back to a slower version. | |
371 // However the returned number is guaranteed to be either the correct | |
372 // double, or the next-lower double. | |
373 return false; | |
374 } else { | |
375 return true; | |
376 } | |
377 } | |
378 | |
379 | |
380 double Strtod(Vector<const char> buffer, int exponent) { | 200 double Strtod(Vector<const char> buffer, int exponent) { |
381 Vector<const char> left_trimmed = TrimLeadingZeros(buffer); | 201 Vector<const char> left_trimmed = TrimLeadingZeros(buffer); |
382 Vector<const char> trimmed = TrimTrailingZeros(left_trimmed); | 202 Vector<const char> trimmed = TrimTrailingZeros(left_trimmed); |
383 exponent += left_trimmed.length() - trimmed.length(); | 203 exponent += left_trimmed.length() - trimmed.length(); |
384 if (trimmed.length() == 0) return 0.0; | 204 if (trimmed.length() == 0) return 0.0; |
385 if (exponent + trimmed.length() - 1 >= kMaxDecimalPower) return V8_INFINITY; | 205 if (exponent + trimmed.length() - 1 >= kMaxDecimalPower) return V8_INFINITY; |
386 if (exponent + trimmed.length() <= kMinDecimalPower) return 0.0; | 206 if (exponent + trimmed.length() <= kMinDecimalPower) return 0.0; |
387 | |
388 double result; | 207 double result; |
389 if (DoubleStrtod(trimmed, exponent, &result) || | 208 if (DoubleStrtod(trimmed, exponent, &result)) { |
390 DiyFpStrtod(trimmed, exponent, &result)) { | |
391 return result; | 209 return result; |
392 } | 210 } |
393 return old_strtod(trimmed, exponent); | 211 return old_strtod(trimmed, exponent); |
394 } | 212 } |
395 | 213 |
396 } } // namespace v8::internal | 214 } } // namespace v8::internal |
OLD | NEW |