1 | /**************************************************************************** |
2 | ** |
3 | ** Copyright (C) 2014 Digia Plc and/or its subsidiary(-ies). |
4 | ** Contact: http://www.qt-project.org/legal |
5 | ** |
6 | ** This file is part of the QtCore module of the Qt Toolkit. |
7 | ** |
8 | ** $QT_BEGIN_LICENSE:LGPL$ |
9 | ** Commercial License Usage |
10 | ** Licensees holding valid commercial Qt licenses may use this file in |
11 | ** accordance with the commercial license agreement provided with the |
12 | ** Software or, alternatively, in accordance with the terms contained in |
13 | ** a written agreement between you and Digia. For licensing terms and |
14 | ** conditions see http://qt.digia.com/licensing. For further information |
15 | ** use the contact form at http://qt.digia.com/contact-us. |
16 | ** |
17 | ** GNU Lesser General Public License Usage |
18 | ** Alternatively, this file may be used under the terms of the GNU Lesser |
19 | ** General Public License version 2.1 as published by the Free Software |
20 | ** Foundation and appearing in the file LICENSE.LGPL included in the |
21 | ** packaging of this file. Please review the following information to |
22 | ** ensure the GNU Lesser General Public License version 2.1 requirements |
23 | ** will be met: http://www.gnu.org/licenses/old-licenses/lgpl-2.1.html. |
24 | ** |
25 | ** In addition, as a special exception, Digia gives you certain additional |
26 | ** rights. These rights are described in the Digia Qt LGPL Exception |
27 | ** version 1.1, included in the file LGPL_EXCEPTION.txt in this package. |
28 | ** |
29 | ** GNU General Public License Usage |
30 | ** Alternatively, this file may be used under the terms of the GNU |
31 | ** General Public License version 3.0 as published by the Free Software |
32 | ** Foundation and appearing in the file LICENSE.GPL included in the |
33 | ** packaging of this file. Please review the following information to |
34 | ** ensure the GNU General Public License version 3.0 requirements will be |
35 | ** met: http://www.gnu.org/copyleft/gpl.html. |
36 | ** |
37 | ** |
38 | ** $QT_END_LICENSE$ |
39 | ** |
40 | ****************************************************************************/ |
41 | |
42 | #include "qlocale_tools_p.h" |
43 | #include "qlocale_p.h" |
44 | #include "qstring.h" |
45 | |
46 | #include <ctype.h> |
47 | #include <float.h> |
48 | #include <limits.h> |
49 | #include <math.h> |
50 | #include <stdlib.h> |
51 | #include <time.h> |
52 | |
53 | #ifdef Q_OS_WINCE |
54 | # include "qfunctions_wince.h" // for _control87 |
55 | #endif |
56 | |
57 | #if defined(Q_OS_LINUX) && !defined(__UCLIBC__) |
58 | # include <fenv.h> |
59 | #endif |
60 | |
61 | // Sizes as defined by the ISO C99 standard - fallback |
62 | #ifndef LLONG_MAX |
63 | # define LLONG_MAX Q_INT64_C(0x7fffffffffffffff) |
64 | #endif |
65 | #ifndef LLONG_MIN |
66 | # define LLONG_MIN (-LLONG_MAX - Q_INT64_C(1)) |
67 | #endif |
68 | #ifndef ULLONG_MAX |
69 | # define ULLONG_MAX Q_UINT64_C(0xffffffffffffffff) |
70 | #endif |
71 | |
72 | QT_BEGIN_NAMESPACE |
73 | |
74 | #ifndef QT_QLOCALE_USES_FCVT |
75 | static char *_qdtoa( NEEDS_VOLATILE double d, int mode, int ndigits, int *decpt, |
76 | int *sign, char **rve, char **digits_str); |
77 | #endif |
78 | |
79 | QString qulltoa(qulonglong l, int base, const QChar _zero) |
80 | { |
81 | ushort buff[65]; // length of MAX_ULLONG in base 2 |
82 | ushort *p = buff + 65; |
83 | |
84 | if (base != 10 || _zero.unicode() == '0') { |
85 | while (l != 0) { |
86 | int c = l % base; |
87 | |
88 | --p; |
89 | |
90 | if (c < 10) |
91 | *p = '0' + c; |
92 | else |
93 | *p = c - 10 + 'a'; |
94 | |
95 | l /= base; |
96 | } |
97 | } |
98 | else { |
99 | while (l != 0) { |
100 | int c = l % base; |
101 | |
102 | *(--p) = _zero.unicode() + c; |
103 | |
104 | l /= base; |
105 | } |
106 | } |
107 | |
108 | return QString(reinterpret_cast<QChar *>(p), 65 - (p - buff)); |
109 | } |
110 | |
111 | QString qlltoa(qlonglong l, int base, const QChar zero) |
112 | { |
113 | return qulltoa(l < 0 ? -l : l, base, zero); |
114 | } |
115 | |
116 | QString &decimalForm(QChar zero, QChar decimal, QChar group, |
117 | QString &digits, int decpt, uint precision, |
118 | PrecisionMode pm, |
119 | bool always_show_decpt, |
120 | bool thousands_group) |
121 | { |
122 | if (decpt < 0) { |
123 | for (int i = 0; i < -decpt; ++i) |
124 | digits.prepend(zero); |
125 | decpt = 0; |
126 | } |
127 | else if (decpt > digits.length()) { |
128 | for (int i = digits.length(); i < decpt; ++i) |
129 | digits.append(zero); |
130 | } |
131 | |
132 | if (pm == PMDecimalDigits) { |
133 | uint decimal_digits = digits.length() - decpt; |
134 | for (uint i = decimal_digits; i < precision; ++i) |
135 | digits.append(zero); |
136 | } |
137 | else if (pm == PMSignificantDigits) { |
138 | for (uint i = digits.length(); i < precision; ++i) |
139 | digits.append(zero); |
140 | } |
141 | else { // pm == PMChopTrailingZeros |
142 | } |
143 | |
144 | if (always_show_decpt || decpt < digits.length()) |
145 | digits.insert(decpt, decimal); |
146 | |
147 | if (thousands_group) { |
148 | for (int i = decpt - 3; i > 0; i -= 3) |
149 | digits.insert(i, group); |
150 | } |
151 | |
152 | if (decpt == 0) |
153 | digits.prepend(zero); |
154 | |
155 | return digits; |
156 | } |
157 | |
158 | QString &exponentForm(QChar zero, QChar decimal, QChar exponential, |
159 | QChar group, QChar plus, QChar minus, |
160 | QString &digits, int decpt, uint precision, |
161 | PrecisionMode pm, |
162 | bool always_show_decpt) |
163 | { |
164 | int exp = decpt - 1; |
165 | |
166 | if (pm == PMDecimalDigits) { |
167 | for (uint i = digits.length(); i < precision + 1; ++i) |
168 | digits.append(zero); |
169 | } |
170 | else if (pm == PMSignificantDigits) { |
171 | for (uint i = digits.length(); i < precision; ++i) |
172 | digits.append(zero); |
173 | } |
174 | else { // pm == PMChopTrailingZeros |
175 | } |
176 | |
177 | if (always_show_decpt || digits.length() > 1) |
178 | digits.insert(1, decimal); |
179 | |
180 | digits.append(exponential); |
181 | digits.append(QLocalePrivate::longLongToString(zero, group, plus, minus, |
182 | exp, 2, 10, -1, QLocalePrivate::AlwaysShowSign)); |
183 | |
184 | return digits; |
185 | } |
186 | |
187 | // Removes thousand-group separators in "C" locale. |
188 | bool removeGroupSeparators(QLocalePrivate::CharBuff *num) |
189 | { |
190 | int group_cnt = 0; // counts number of group chars |
191 | int decpt_idx = -1; |
192 | |
193 | char *data = num->data(); |
194 | int l = qstrlen(data); |
195 | |
196 | // Find the decimal point and check if there are any group chars |
197 | int i = 0; |
198 | for (; i < l; ++i) { |
199 | char c = data[i]; |
200 | |
201 | if (c == ',') { |
202 | if (i == 0 || data[i - 1] < '0' || data[i - 1] > '9') |
203 | return false; |
204 | if (i == l - 1 || data[i + 1] < '0' || data[i + 1] > '9') |
205 | return false; |
206 | ++group_cnt; |
207 | } |
208 | else if (c == '.') { |
209 | // Fail if more than one decimal points |
210 | if (decpt_idx != -1) |
211 | return false; |
212 | decpt_idx = i; |
213 | } else if (c == 'e' || c == 'E') { |
214 | // an 'e' or 'E' - if we have not encountered a decimal |
215 | // point, this is where it "is". |
216 | if (decpt_idx == -1) |
217 | decpt_idx = i; |
218 | } |
219 | } |
220 | |
221 | // If no group chars, we're done |
222 | if (group_cnt == 0) |
223 | return true; |
224 | |
225 | // No decimal point means that it "is" at the end of the string |
226 | if (decpt_idx == -1) |
227 | decpt_idx = l; |
228 | |
229 | i = 0; |
230 | while (i < l && group_cnt > 0) { |
231 | char c = data[i]; |
232 | |
233 | if (c == ',') { |
234 | // Don't allow group chars after the decimal point |
235 | if (i > decpt_idx) |
236 | return false; |
237 | |
238 | // Check that it is placed correctly relative to the decpt |
239 | if ((decpt_idx - i) % 4 != 0) |
240 | return false; |
241 | |
242 | // Remove it |
243 | memmove(data + i, data + i + 1, l - i - 1); |
244 | data[--l] = '\0'; |
245 | |
246 | --group_cnt; |
247 | --decpt_idx; |
248 | } else { |
249 | // Check that we are not missing a separator |
250 | if (i < decpt_idx |
251 | && (decpt_idx - i) % 4 == 0 |
252 | && !(i == 0 && c == '-')) // check for negative sign at start of string |
253 | return false; |
254 | ++i; |
255 | } |
256 | } |
257 | |
258 | return true; |
259 | } |
260 | |
261 | #if defined(Q_CC_MWERKS) && defined(Q_OS_WIN32) |
262 | inline bool isascii(int c) |
263 | { |
264 | return (c >= 0 && c <=127); |
265 | } |
266 | #endif |
267 | |
268 | /*- |
269 | * Copyright (c) 1992, 1993 |
270 | * The Regents of the University of California. All rights reserved. |
271 | * |
272 | * Redistribution and use in source and binary forms, with or without |
273 | * modification, are permitted provided that the following conditions |
274 | * are met: |
275 | * 1. Redistributions of source code must retain the above copyright |
276 | * notice, this list of conditions and the following disclaimer. |
277 | * 2. Redistributions in binary form must reproduce the above copyright |
278 | * notice, this list of conditions and the following disclaimer in the |
279 | * documentation and/or other materials provided with the distribution. |
280 | * 3. All advertising materials mentioning features or use of this software |
281 | * must display the following acknowledgment: |
282 | * This product includes software developed by the University of |
283 | * California, Berkeley and its contributors. |
284 | * 4. Neither the name of the University nor the names of its contributors |
285 | * may be used to endorse or promote products derived from this software |
286 | * without specific prior written permission. |
287 | * |
288 | * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND |
289 | * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE |
290 | * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE |
291 | * ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE |
292 | * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL |
293 | * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS |
294 | * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) |
295 | * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT |
296 | * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY |
297 | * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF |
298 | * SUCH DAMAGE. |
299 | */ |
300 | |
301 | // static char sccsid[] = "@(#)strtouq.c 8.1 (Berkeley) 6/4/93"; |
302 | // "$FreeBSD: src/lib/libc/stdlib/strtoull.c,v 1.5.2.1 2001/03/02 09:45:20 obrien Exp $"; |
303 | |
304 | /* |
305 | * Convert a string to an unsigned long long integer. |
306 | * |
307 | * Ignores `locale' stuff. Assumes that the upper and lower case |
308 | * alphabets and digits are each contiguous. |
309 | */ |
310 | qulonglong qstrtoull(const char *nptr, const char **endptr, register int base, bool *ok) |
311 | { |
312 | register const char *s = nptr; |
313 | register qulonglong acc; |
314 | register unsigned char c; |
315 | register qulonglong qbase, cutoff; |
316 | register int any, cutlim; |
317 | |
318 | if (ok != 0) |
319 | *ok = true; |
320 | |
321 | /* |
322 | * See strtoq for comments as to the logic used. |
323 | */ |
324 | s = nptr; |
325 | do { |
326 | c = *s++; |
327 | } while (isspace(c)); |
328 | if (c == '-') { |
329 | if (ok != 0) |
330 | *ok = false; |
331 | if (endptr != 0) |
332 | *endptr = s - 1; |
333 | return 0; |
334 | } else { |
335 | if (c == '+') |
336 | c = *s++; |
337 | } |
338 | if ((base == 0 || base == 16) && |
339 | c == '0' && (*s == 'x' || *s == 'X')) { |
340 | c = s[1]; |
341 | s += 2; |
342 | base = 16; |
343 | } |
344 | if (base == 0) |
345 | base = c == '0' ? 8 : 10; |
346 | qbase = unsigned(base); |
347 | cutoff = qulonglong(ULLONG_MAX) / qbase; |
348 | cutlim = qulonglong(ULLONG_MAX) % qbase; |
349 | for (acc = 0, any = 0;; c = *s++) { |
350 | if (!isascii(c)) |
351 | break; |
352 | if (isdigit(c)) |
353 | c -= '0'; |
354 | else if (isalpha(c)) |
355 | c -= isupper(c) ? 'A' - 10 : 'a' - 10; |
356 | else |
357 | break; |
358 | if (c >= base) |
359 | break; |
360 | if (any < 0 || acc > cutoff || (acc == cutoff && c > cutlim)) |
361 | any = -1; |
362 | else { |
363 | any = 1; |
364 | acc *= qbase; |
365 | acc += c; |
366 | } |
367 | } |
368 | if (any == 0) { |
369 | if (ok != 0) |
370 | *ok = false; |
371 | } else if (any < 0) { |
372 | acc = ULLONG_MAX; |
373 | if (ok != 0) |
374 | *ok = false; |
375 | } |
376 | if (endptr != 0) |
377 | *endptr = (any ? s - 1 : nptr); |
378 | return acc; |
379 | } |
380 | |
381 | |
382 | // "$FreeBSD: src/lib/libc/stdlib/strtoll.c,v 1.5.2.1 2001/03/02 09:45:20 obrien Exp $"; |
383 | |
384 | |
385 | /* |
386 | * Convert a string to a long long integer. |
387 | * |
388 | * Ignores `locale' stuff. Assumes that the upper and lower case |
389 | * alphabets and digits are each contiguous. |
390 | */ |
391 | qlonglong qstrtoll(const char *nptr, const char **endptr, register int base, bool *ok) |
392 | { |
393 | register const char *s; |
394 | register qulonglong acc; |
395 | register unsigned char c; |
396 | register qulonglong qbase, cutoff; |
397 | register int neg, any, cutlim; |
398 | |
399 | /* |
400 | * Skip white space and pick up leading +/- sign if any. |
401 | * If base is 0, allow 0x for hex and 0 for octal, else |
402 | * assume decimal; if base is already 16, allow 0x. |
403 | */ |
404 | s = nptr; |
405 | do { |
406 | c = *s++; |
407 | } while (isspace(c)); |
408 | if (c == '-') { |
409 | neg = 1; |
410 | c = *s++; |
411 | } else { |
412 | neg = 0; |
413 | if (c == '+') |
414 | c = *s++; |
415 | } |
416 | if ((base == 0 || base == 16) && |
417 | c == '0' && (*s == 'x' || *s == 'X')) { |
418 | c = s[1]; |
419 | s += 2; |
420 | base = 16; |
421 | } |
422 | if (base == 0) |
423 | base = c == '0' ? 8 : 10; |
424 | |
425 | /* |
426 | * Compute the cutoff value between legal numbers and illegal |
427 | * numbers. That is the largest legal value, divided by the |
428 | * base. An input number that is greater than this value, if |
429 | * followed by a legal input character, is too big. One that |
430 | * is equal to this value may be valid or not; the limit |
431 | * between valid and invalid numbers is then based on the last |
432 | * digit. For instance, if the range for quads is |
433 | * [-9223372036854775808..9223372036854775807] and the input base |
434 | * is 10, cutoff will be set to 922337203685477580 and cutlim to |
435 | * either 7 (neg==0) or 8 (neg==1), meaning that if we have |
436 | * accumulated a value > 922337203685477580, or equal but the |
437 | * next digit is > 7 (or 8), the number is too big, and we will |
438 | * return a range error. |
439 | * |
440 | * Set any if any `digits' consumed; make it negative to indicate |
441 | * overflow. |
442 | */ |
443 | qbase = unsigned(base); |
444 | cutoff = neg ? qulonglong(0-(LLONG_MIN + LLONG_MAX)) + LLONG_MAX : LLONG_MAX; |
445 | cutlim = cutoff % qbase; |
446 | cutoff /= qbase; |
447 | for (acc = 0, any = 0;; c = *s++) { |
448 | if (!isascii(c)) |
449 | break; |
450 | if (isdigit(c)) |
451 | c -= '0'; |
452 | else if (isalpha(c)) |
453 | c -= isupper(c) ? 'A' - 10 : 'a' - 10; |
454 | else |
455 | break; |
456 | if (c >= base) |
457 | break; |
458 | if (any < 0 || acc > cutoff || (acc == cutoff && c > cutlim)) |
459 | any = -1; |
460 | else { |
461 | any = 1; |
462 | acc *= qbase; |
463 | acc += c; |
464 | } |
465 | } |
466 | if (any < 0) { |
467 | acc = neg ? LLONG_MIN : LLONG_MAX; |
468 | if (ok != 0) |
469 | *ok = false; |
470 | } else if (neg) { |
471 | acc = (~acc) + 1; |
472 | } |
473 | if (endptr != 0) |
474 | *endptr = (any >= 0 ? s - 1 : nptr); |
475 | |
476 | if (ok != 0) |
477 | *ok = any > 0; |
478 | |
479 | return acc; |
480 | } |
481 | |
482 | #ifndef QT_QLOCALE_USES_FCVT |
483 | |
484 | /* From: NetBSD: strtod.c,v 1.26 1998/02/03 18:44:21 perry Exp */ |
485 | /* $FreeBSD: src/lib/libc/stdlib/netbsd_strtod.c,v 1.2.2.2 2001/03/02 17:14:15 tegge Exp $ */ |
486 | |
487 | /* Please send bug reports to |
488 | David M. Gay |
489 | AT&T Bell Laboratories, Room 2C-463 |
490 | 600 Mountain Avenue |
491 | Murray Hill, NJ 07974-2070 |
492 | U.S.A. |
493 | dmg@research.att.com or research!dmg |
494 | */ |
495 | |
496 | /* strtod for IEEE-, VAX-, and IBM-arithmetic machines. |
497 | * |
498 | * This strtod returns a nearest machine number to the input decimal |
499 | * string (or sets errno to ERANGE). With IEEE arithmetic, ties are |
500 | * broken by the IEEE round-even rule. Otherwise ties are broken by |
501 | * biased rounding (add half and chop). |
502 | * |
503 | * Inspired loosely by William D. Clinger's paper "How to Read Floating |
504 | * Point Numbers Accurately" [Proc. ACM SIGPLAN '90, pp. 92-101]. |
505 | * |
506 | * Modifications: |
507 | * |
508 | * 1. We only require IEEE, IBM, or VAX double-precision |
509 | * arithmetic (not IEEE double-extended). |
510 | * 2. We get by with floating-point arithmetic in a case that |
511 | * Clinger missed -- when we're computing d * 10^n |
512 | * for a small integer d and the integer n is not too |
513 | * much larger than 22 (the maximum integer k for which |
514 | * we can represent 10^k exactly), we may be able to |
515 | * compute (d*10^k) * 10^(e-k) with just one roundoff. |
516 | * 3. Rather than a bit-at-a-time adjustment of the binary |
517 | * result in the hard case, we use floating-point |
518 | * arithmetic to determine the adjustment to within |
519 | * one bit; only in really hard cases do we need to |
520 | * compute a second residual. |
521 | * 4. Because of 3., we don't need a large table of powers of 10 |
522 | * for ten-to-e (just some small tables, e.g. of 10^k |
523 | * for 0 <= k <= 22). |
524 | */ |
525 | |
526 | /* |
527 | * #define IEEE_LITTLE_ENDIAN for IEEE-arithmetic machines where the least |
528 | * significant byte has the lowest address. |
529 | * #define IEEE_BIG_ENDIAN for IEEE-arithmetic machines where the most |
530 | * significant byte has the lowest address. |
531 | * #define Long int on machines with 32-bit ints and 64-bit longs. |
532 | * #define Sudden_Underflow for IEEE-format machines without gradual |
533 | * underflow (i.e., that flush to zero on underflow). |
534 | * #define IBM for IBM mainframe-style floating-point arithmetic. |
535 | * #define VAX for VAX-style floating-point arithmetic. |
536 | * #define Unsigned_Shifts if >> does treats its left operand as unsigned. |
537 | * #define No_leftright to omit left-right logic in fast floating-point |
538 | * computation of dtoa. |
539 | * #define Check_FLT_ROUNDS if FLT_ROUNDS can assume the values 2 or 3. |
540 | * #define RND_PRODQUOT to use rnd_prod and rnd_quot (assembly routines |
541 | * that use extended-precision instructions to compute rounded |
542 | * products and quotients) with IBM. |
543 | * #define ROUND_BIASED for IEEE-format with biased rounding. |
544 | * #define Inaccurate_Divide for IEEE-format with correctly rounded |
545 | * products but inaccurate quotients, e.g., for Intel i860. |
546 | * #define Just_16 to store 16 bits per 32-bit Long when doing high-precision |
547 | * integer arithmetic. Whether this speeds things up or slows things |
548 | * down depends on the machine and the number being converted. |
549 | * #define KR_headers for old-style C function headers. |
550 | * #define Bad_float_h if your system lacks a float.h or if it does not |
551 | * define some or all of DBL_DIG, DBL_MAX_10_EXP, DBL_MAX_EXP, |
552 | * FLT_RADIX, FLT_ROUNDS, and DBL_MAX. |
553 | * #define MALLOC your_malloc, where your_malloc(n) acts like malloc(n) |
554 | * if memory is available and otherwise does something you deem |
555 | * appropriate. If MALLOC is undefined, malloc will be invoked |
556 | * directly -- and assumed always to succeed. |
557 | */ |
558 | |
559 | #if defined(LIBC_SCCS) && !defined(lint) |
560 | __RCSID("$NetBSD: strtod.c,v 1.26 1998/02/03 18:44:21 perry Exp $" ); |
561 | #endif /* LIBC_SCCS and not lint */ |
562 | |
563 | /* |
564 | #if defined(__m68k__) || defined(__sparc__) || defined(__i386__) || \ |
565 | defined(__mips__) || defined(__ns32k__) || defined(__alpha__) || \ |
566 | defined(__powerpc__) || defined(Q_OS_WIN) || defined(Q_OS_DARWIN) || defined(Q_OS_MAC) || \ |
567 | defined(mips) || defined(Q_OS_AIX) || defined(Q_OS_SOLARIS) |
568 | # define IEEE_BIG_OR_LITTLE_ENDIAN 1 |
569 | #endif |
570 | */ |
571 | |
572 | // *All* of our architectures have IEEE arithmetic, don't they? |
573 | #define IEEE_BIG_OR_LITTLE_ENDIAN 1 |
574 | |
575 | #ifdef __arm32__ |
576 | /* |
577 | * Although the CPU is little endian the FP has different |
578 | * byte and word endianness. The byte order is still little endian |
579 | * but the word order is big endian. |
580 | */ |
581 | #define IEEE_BIG_OR_LITTLE_ENDIAN |
582 | #endif |
583 | |
584 | #ifdef vax |
585 | #define VAX |
586 | #endif |
587 | |
588 | #define Long qint32 |
589 | #define ULong quint32 |
590 | |
591 | #define MALLOC malloc |
592 | |
593 | #ifdef BSD_QDTOA_DEBUG |
594 | QT_BEGIN_INCLUDE_NAMESPACE |
595 | #include <stdio.h> |
596 | QT_END_INCLUDE_NAMESPACE |
597 | |
598 | #define Bug(x) {fprintf(stderr, "%s\n", x); exit(1);} |
599 | #endif |
600 | |
601 | #ifdef Unsigned_Shifts |
602 | #define Sign_Extend(a,b) if (b < 0) a |= 0xffff0000; |
603 | #else |
604 | #define Sign_Extend(a,b) /*no-op*/ |
605 | #endif |
606 | |
607 | #if (defined(IEEE_BIG_OR_LITTLE_ENDIAN) + defined(VAX) + defined(IBM)) != 1 |
608 | #error Exactly one of IEEE_BIG_OR_LITTLE_ENDIAN, VAX, or IBM should be defined. |
609 | #endif |
610 | |
611 | static inline ULong _getWord0(const NEEDS_VOLATILE double x) |
612 | { |
613 | const NEEDS_VOLATILE uchar *ptr = reinterpret_cast<const NEEDS_VOLATILE uchar *>(&x); |
614 | if (QSysInfo::ByteOrder == QSysInfo::BigEndian) { |
615 | return (ptr[0]<<24) + (ptr[1]<<16) + (ptr[2]<<8) + ptr[3]; |
616 | } else { |
617 | return (ptr[7]<<24) + (ptr[6]<<16) + (ptr[5]<<8) + ptr[4]; |
618 | } |
619 | } |
620 | |
621 | static inline void _setWord0(NEEDS_VOLATILE double *x, ULong l) |
622 | { |
623 | NEEDS_VOLATILE uchar *ptr = reinterpret_cast<NEEDS_VOLATILE uchar *>(x); |
624 | if (QSysInfo::ByteOrder == QSysInfo::BigEndian) { |
625 | ptr[0] = uchar(l>>24); |
626 | ptr[1] = uchar(l>>16); |
627 | ptr[2] = uchar(l>>8); |
628 | ptr[3] = uchar(l); |
629 | } else { |
630 | ptr[7] = uchar(l>>24); |
631 | ptr[6] = uchar(l>>16); |
632 | ptr[5] = uchar(l>>8); |
633 | ptr[4] = uchar(l); |
634 | } |
635 | } |
636 | |
637 | static inline ULong _getWord1(const NEEDS_VOLATILE double x) |
638 | { |
639 | const NEEDS_VOLATILE uchar *ptr = reinterpret_cast<const NEEDS_VOLATILE uchar *>(&x); |
640 | if (QSysInfo::ByteOrder == QSysInfo::BigEndian) { |
641 | return (ptr[4]<<24) + (ptr[5]<<16) + (ptr[6]<<8) + ptr[7]; |
642 | } else { |
643 | return (ptr[3]<<24) + (ptr[2]<<16) + (ptr[1]<<8) + ptr[0]; |
644 | } |
645 | } |
646 | static inline void _setWord1(NEEDS_VOLATILE double *x, ULong l) |
647 | { |
648 | NEEDS_VOLATILE uchar *ptr = reinterpret_cast<uchar NEEDS_VOLATILE *>(x); |
649 | if (QSysInfo::ByteOrder == QSysInfo::BigEndian) { |
650 | ptr[4] = uchar(l>>24); |
651 | ptr[5] = uchar(l>>16); |
652 | ptr[6] = uchar(l>>8); |
653 | ptr[7] = uchar(l); |
654 | } else { |
655 | ptr[3] = uchar(l>>24); |
656 | ptr[2] = uchar(l>>16); |
657 | ptr[1] = uchar(l>>8); |
658 | ptr[0] = uchar(l); |
659 | } |
660 | } |
661 | |
662 | static inline ULong getWord0(const NEEDS_VOLATILE double x) |
663 | { |
664 | #ifdef QT_ARMFPA |
665 | return _getWord1(x); |
666 | #else |
667 | return _getWord0(x); |
668 | #endif |
669 | } |
670 | |
671 | static inline void setWord0(NEEDS_VOLATILE double *x, ULong l) |
672 | { |
673 | #ifdef QT_ARMFPA |
674 | _setWord1(x, l); |
675 | #else |
676 | _setWord0(x, l); |
677 | #endif |
678 | } |
679 | |
680 | static inline ULong getWord1(const NEEDS_VOLATILE double x) |
681 | { |
682 | #ifdef QT_ARMFPA |
683 | return _getWord0(x); |
684 | #else |
685 | return _getWord1(x); |
686 | #endif |
687 | } |
688 | |
689 | static inline void setWord1(NEEDS_VOLATILE double *x, ULong l) |
690 | { |
691 | #ifdef QT_ARMFPA |
692 | _setWord0(x, l); |
693 | #else |
694 | _setWord1(x, l); |
695 | #endif |
696 | } |
697 | |
698 | static inline void Storeinc(ULong *&a, const ULong &b, const ULong &c) |
699 | { |
700 | |
701 | *a = (ushort(b) << 16) | ushort(c); |
702 | ++a; |
703 | } |
704 | |
705 | /* #define P DBL_MANT_DIG */ |
706 | /* Ten_pmax = floor(P*log(2)/log(5)) */ |
707 | /* Bletch = (highest power of 2 < DBL_MAX_10_EXP) / 16 */ |
708 | /* Quick_max = floor((P-1)*log(FLT_RADIX)/log(10) - 1) */ |
709 | /* Int_max = floor(P*log(FLT_RADIX)/log(10) - 1) */ |
710 | |
711 | #if defined(IEEE_BIG_OR_LITTLE_ENDIAN) |
712 | #define Exp_shift 20 |
713 | #define Exp_shift1 20 |
714 | #define Exp_msk1 0x100000 |
715 | #define Exp_msk11 0x100000 |
716 | #define Exp_mask 0x7ff00000 |
717 | #define P 53 |
718 | #define Bias 1023 |
719 | #define IEEE_Arith |
720 | #define Emin (-1022) |
721 | #define Exp_1 0x3ff00000 |
722 | #define Exp_11 0x3ff00000 |
723 | #define Ebits 11 |
724 | #define Frac_mask 0xfffff |
725 | #define Frac_mask1 0xfffff |
726 | #define Ten_pmax 22 |
727 | #define Bletch 0x10 |
728 | #define Bndry_mask 0xfffff |
729 | #define Bndry_mask1 0xfffff |
730 | #if defined(LSB) && defined(Q_OS_VXWORKS) |
731 | #undef LSB |
732 | #endif |
733 | #define LSB 1 |
734 | #define Sign_bit 0x80000000 |
735 | #define Log2P 1 |
736 | #define Tiny0 0 |
737 | #define Tiny1 1 |
738 | #define Quick_max 14 |
739 | #define Int_max 14 |
740 | #define Infinite(x) (getWord0(x) == 0x7ff00000) /* sufficient test for here */ |
741 | #else |
742 | #undef Sudden_Underflow |
743 | #define Sudden_Underflow |
744 | #ifdef IBM |
745 | #define Exp_shift 24 |
746 | #define Exp_shift1 24 |
747 | #define Exp_msk1 0x1000000 |
748 | #define Exp_msk11 0x1000000 |
749 | #define Exp_mask 0x7f000000 |
750 | #define P 14 |
751 | #define Bias 65 |
752 | #define Exp_1 0x41000000 |
753 | #define Exp_11 0x41000000 |
754 | #define Ebits 8 /* exponent has 7 bits, but 8 is the right value in b2d */ |
755 | #define Frac_mask 0xffffff |
756 | #define Frac_mask1 0xffffff |
757 | #define Bletch 4 |
758 | #define Ten_pmax 22 |
759 | #define Bndry_mask 0xefffff |
760 | #define Bndry_mask1 0xffffff |
761 | #define LSB 1 |
762 | #define Sign_bit 0x80000000 |
763 | #define Log2P 4 |
764 | #define Tiny0 0x100000 |
765 | #define Tiny1 0 |
766 | #define Quick_max 14 |
767 | #define Int_max 15 |
768 | #else /* VAX */ |
769 | #define Exp_shift 23 |
770 | #define Exp_shift1 7 |
771 | #define Exp_msk1 0x80 |
772 | #define Exp_msk11 0x800000 |
773 | #define Exp_mask 0x7f80 |
774 | #define P 56 |
775 | #define Bias 129 |
776 | #define Exp_1 0x40800000 |
777 | #define Exp_11 0x4080 |
778 | #define Ebits 8 |
779 | #define Frac_mask 0x7fffff |
780 | #define Frac_mask1 0xffff007f |
781 | #define Ten_pmax 24 |
782 | #define Bletch 2 |
783 | #define Bndry_mask 0xffff007f |
784 | #define Bndry_mask1 0xffff007f |
785 | #define LSB 0x10000 |
786 | #define Sign_bit 0x8000 |
787 | #define Log2P 1 |
788 | #define Tiny0 0x80 |
789 | #define Tiny1 0 |
790 | #define Quick_max 15 |
791 | #define Int_max 15 |
792 | #endif |
793 | #endif |
794 | |
795 | #ifndef IEEE_Arith |
796 | #define ROUND_BIASED |
797 | #endif |
798 | |
799 | #ifdef RND_PRODQUOT |
800 | #define rounded_product(a,b) a = rnd_prod(a, b) |
801 | #define rounded_quotient(a,b) a = rnd_quot(a, b) |
802 | extern double rnd_prod(double, double), rnd_quot(double, double); |
803 | #else |
804 | #define rounded_product(a,b) a *= b |
805 | #define rounded_quotient(a,b) a /= b |
806 | #endif |
807 | |
808 | #define Big0 (Frac_mask1 | Exp_msk1*(DBL_MAX_EXP+Bias-1)) |
809 | #define Big1 0xffffffff |
810 | |
811 | #ifndef Just_16 |
812 | /* When Pack_32 is not defined, we store 16 bits per 32-bit Long. |
813 | * This makes some inner loops simpler and sometimes saves work |
814 | * during multiplications, but it often seems to make things slightly |
815 | * slower. Hence the default is now to store 32 bits per Long. |
816 | */ |
817 | #ifndef Pack_32 |
818 | #define Pack_32 |
819 | #endif |
820 | #endif |
821 | |
822 | #define Kmax 15 |
823 | |
824 | struct |
825 | Bigint { |
826 | struct Bigint *next; |
827 | int k, maxwds, sign, wds; |
828 | ULong x[1]; |
829 | }; |
830 | |
831 | typedef struct Bigint Bigint; |
832 | |
833 | static Bigint *Balloc(int k) |
834 | { |
835 | int x; |
836 | Bigint *rv; |
837 | |
838 | x = 1 << k; |
839 | rv = static_cast<Bigint *>(MALLOC(sizeof(Bigint) + (x-1)*sizeof(Long))); |
840 | Q_CHECK_PTR(rv); |
841 | rv->k = k; |
842 | rv->maxwds = x; |
843 | rv->sign = rv->wds = 0; |
844 | return rv; |
845 | } |
846 | |
847 | static void Bfree(Bigint *v) |
848 | { |
849 | free(v); |
850 | } |
851 | |
852 | #define Bcopy(x,y) memcpy(reinterpret_cast<char *>(&x->sign), reinterpret_cast<char *>(&y->sign), \ |
853 | y->wds*sizeof(Long) + 2*sizeof(int)) |
854 | |
855 | /* multiply by m and add a */ |
856 | static Bigint *multadd(Bigint *b, int m, int a) |
857 | { |
858 | int i, wds; |
859 | ULong *x, y; |
860 | #ifdef Pack_32 |
861 | ULong xi, z; |
862 | #endif |
863 | Bigint *b1; |
864 | |
865 | wds = b->wds; |
866 | x = b->x; |
867 | i = 0; |
868 | do { |
869 | #ifdef Pack_32 |
870 | xi = *x; |
871 | y = (xi & 0xffff) * m + a; |
872 | z = (xi >> 16) * m + (y >> 16); |
873 | a = (z >> 16); |
874 | *x++ = (z << 16) + (y & 0xffff); |
875 | #else |
876 | y = *x * m + a; |
877 | a = (y >> 16); |
878 | *x++ = y & 0xffff; |
879 | #endif |
880 | } |
881 | while(++i < wds); |
882 | if (a) { |
883 | if (wds >= b->maxwds) { |
884 | b1 = Balloc(b->k+1); |
885 | Bcopy(b1, b); |
886 | Bfree(b); |
887 | b = b1; |
888 | } |
889 | b->x[wds++] = a; |
890 | b->wds = wds; |
891 | } |
892 | return b; |
893 | } |
894 | |
895 | static Bigint *s2b(const char *s, int nd0, int nd, ULong y9) |
896 | { |
897 | Bigint *b; |
898 | int i, k; |
899 | Long x, y; |
900 | |
901 | x = (nd + 8) / 9; |
902 | for(k = 0, y = 1; x > y; y <<= 1, k++) ; |
903 | #ifdef Pack_32 |
904 | b = Balloc(k); |
905 | b->x[0] = y9; |
906 | b->wds = 1; |
907 | #else |
908 | b = Balloc(k+1); |
909 | b->x[0] = y9 & 0xffff; |
910 | b->wds = (b->x[1] = y9 >> 16) ? 2 : 1; |
911 | #endif |
912 | |
913 | i = 9; |
914 | if (9 < nd0) { |
915 | s += 9; |
916 | do b = multadd(b, 10, *s++ - '0'); |
917 | while(++i < nd0); |
918 | s++; |
919 | } |
920 | else |
921 | s += 10; |
922 | for(; i < nd; i++) |
923 | b = multadd(b, 10, *s++ - '0'); |
924 | return b; |
925 | } |
926 | |
927 | static int hi0bits(ULong x) |
928 | { |
929 | int k = 0; |
930 | |
931 | if (!(x & 0xffff0000)) { |
932 | k = 16; |
933 | x <<= 16; |
934 | } |
935 | if (!(x & 0xff000000)) { |
936 | k += 8; |
937 | x <<= 8; |
938 | } |
939 | if (!(x & 0xf0000000)) { |
940 | k += 4; |
941 | x <<= 4; |
942 | } |
943 | if (!(x & 0xc0000000)) { |
944 | k += 2; |
945 | x <<= 2; |
946 | } |
947 | if (!(x & 0x80000000)) { |
948 | k++; |
949 | if (!(x & 0x40000000)) |
950 | return 32; |
951 | } |
952 | return k; |
953 | } |
954 | |
955 | static int lo0bits(ULong *y) |
956 | { |
957 | int k; |
958 | ULong x = *y; |
959 | |
960 | if (x & 7) { |
961 | if (x & 1) |
962 | return 0; |
963 | if (x & 2) { |
964 | *y = x >> 1; |
965 | return 1; |
966 | } |
967 | *y = x >> 2; |
968 | return 2; |
969 | } |
970 | k = 0; |
971 | if (!(x & 0xffff)) { |
972 | k = 16; |
973 | x >>= 16; |
974 | } |
975 | if (!(x & 0xff)) { |
976 | k += 8; |
977 | x >>= 8; |
978 | } |
979 | if (!(x & 0xf)) { |
980 | k += 4; |
981 | x >>= 4; |
982 | } |
983 | if (!(x & 0x3)) { |
984 | k += 2; |
985 | x >>= 2; |
986 | } |
987 | if (!(x & 1)) { |
988 | k++; |
989 | x >>= 1; |
990 | if (!x & 1) |
991 | return 32; |
992 | } |
993 | *y = x; |
994 | return k; |
995 | } |
996 | |
997 | static Bigint *i2b(int i) |
998 | { |
999 | Bigint *b; |
1000 | |
1001 | b = Balloc(1); |
1002 | b->x[0] = i; |
1003 | b->wds = 1; |
1004 | return b; |
1005 | } |
1006 | |
1007 | static Bigint *mult(Bigint *a, Bigint *b) |
1008 | { |
1009 | Bigint *c; |
1010 | int k, wa, wb, wc; |
1011 | ULong carry, y, z; |
1012 | ULong *x, *xa, *xae, *xb, *xbe, *xc, *xc0; |
1013 | #ifdef Pack_32 |
1014 | ULong z2; |
1015 | #endif |
1016 | |
1017 | if (a->wds < b->wds) { |
1018 | c = a; |
1019 | a = b; |
1020 | b = c; |
1021 | } |
1022 | k = a->k; |
1023 | wa = a->wds; |
1024 | wb = b->wds; |
1025 | wc = wa + wb; |
1026 | if (wc > a->maxwds) |
1027 | k++; |
1028 | c = Balloc(k); |
1029 | for(x = c->x, xa = x + wc; x < xa; x++) |
1030 | *x = 0; |
1031 | xa = a->x; |
1032 | xae = xa + wa; |
1033 | xb = b->x; |
1034 | xbe = xb + wb; |
1035 | xc0 = c->x; |
1036 | #ifdef Pack_32 |
1037 | for(; xb < xbe; xb++, xc0++) { |
1038 | if ((y = *xb & 0xffff) != 0) { |
1039 | x = xa; |
1040 | xc = xc0; |
1041 | carry = 0; |
1042 | do { |
1043 | z = (*x & 0xffff) * y + (*xc & 0xffff) + carry; |
1044 | carry = z >> 16; |
1045 | z2 = (*x++ >> 16) * y + (*xc >> 16) + carry; |
1046 | carry = z2 >> 16; |
1047 | Storeinc(xc, z2, z); |
1048 | } |
1049 | while(x < xae); |
1050 | *xc = carry; |
1051 | } |
1052 | if ((y = *xb >> 16) != 0) { |
1053 | x = xa; |
1054 | xc = xc0; |
1055 | carry = 0; |
1056 | z2 = *xc; |
1057 | do { |
1058 | z = (*x & 0xffff) * y + (*xc >> 16) + carry; |
1059 | carry = z >> 16; |
1060 | Storeinc(xc, z, z2); |
1061 | z2 = (*x++ >> 16) * y + (*xc & 0xffff) + carry; |
1062 | carry = z2 >> 16; |
1063 | } |
1064 | while(x < xae); |
1065 | *xc = z2; |
1066 | } |
1067 | } |
1068 | #else |
1069 | for(; xb < xbe; xc0++) { |
1070 | if (y = *xb++) { |
1071 | x = xa; |
1072 | xc = xc0; |
1073 | carry = 0; |
1074 | do { |
1075 | z = *x++ * y + *xc + carry; |
1076 | carry = z >> 16; |
1077 | *xc++ = z & 0xffff; |
1078 | } |
1079 | while(x < xae); |
1080 | *xc = carry; |
1081 | } |
1082 | } |
1083 | #endif |
1084 | for(xc0 = c->x, xc = xc0 + wc; wc > 0 && !*--xc; --wc) ; |
1085 | c->wds = wc; |
1086 | return c; |
1087 | } |
1088 | |
1089 | static Bigint *p5s; |
1090 | |
1091 | struct p5s_deleter |
1092 | { |
1093 | ~p5s_deleter() |
1094 | { |
1095 | while (p5s) { |
1096 | Bigint *next = p5s->next; |
1097 | Bfree(p5s); |
1098 | p5s = next; |
1099 | } |
1100 | } |
1101 | }; |
1102 | |
1103 | static Bigint *pow5mult(Bigint *b, int k) |
1104 | { |
1105 | Bigint *b1, *p5, *p51; |
1106 | int i; |
1107 | static const int p05[3] = { 5, 25, 125 }; |
1108 | |
1109 | if ((i = k & 3) != 0) |
1110 | #if defined(Q_OS_IRIX) && defined(Q_CC_GNU) |
1111 | { |
1112 | // work around a bug on 64 bit IRIX gcc |
1113 | int *p = (int *) p05; |
1114 | b = multadd(b, p[i-1], 0); |
1115 | } |
1116 | #else |
1117 | b = multadd(b, p05[i-1], 0); |
1118 | #endif |
1119 | |
1120 | if (!(k >>= 2)) |
1121 | return b; |
1122 | if (!(p5 = p5s)) { |
1123 | /* first time */ |
1124 | static p5s_deleter deleter; |
1125 | p5 = p5s = i2b(625); |
1126 | p5->next = 0; |
1127 | } |
1128 | for(;;) { |
1129 | if (k & 1) { |
1130 | b1 = mult(b, p5); |
1131 | Bfree(b); |
1132 | b = b1; |
1133 | } |
1134 | if (!(k >>= 1)) |
1135 | break; |
1136 | if (!(p51 = p5->next)) { |
1137 | p51 = p5->next = mult(p5,p5); |
1138 | p51->next = 0; |
1139 | } |
1140 | p5 = p51; |
1141 | } |
1142 | return b; |
1143 | } |
1144 | |
1145 | static Bigint *lshift(Bigint *b, int k) |
1146 | { |
1147 | int i, k1, n, n1; |
1148 | Bigint *b1; |
1149 | ULong *x, *x1, *xe, z; |
1150 | |
1151 | #ifdef Pack_32 |
1152 | n = k >> 5; |
1153 | #else |
1154 | n = k >> 4; |
1155 | #endif |
1156 | k1 = b->k; |
1157 | n1 = n + b->wds + 1; |
1158 | for(i = b->maxwds; n1 > i; i <<= 1) |
1159 | k1++; |
1160 | b1 = Balloc(k1); |
1161 | x1 = b1->x; |
1162 | for(i = 0; i < n; i++) |
1163 | *x1++ = 0; |
1164 | x = b->x; |
1165 | xe = x + b->wds; |
1166 | #ifdef Pack_32 |
1167 | if (k &= 0x1f) { |
1168 | k1 = 32 - k; |
1169 | z = 0; |
1170 | do { |
1171 | *x1++ = *x << k | z; |
1172 | z = *x++ >> k1; |
1173 | } |
1174 | while(x < xe); |
1175 | if ((*x1 = z) != 0) |
1176 | ++n1; |
1177 | } |
1178 | #else |
1179 | if (k &= 0xf) { |
1180 | k1 = 16 - k; |
1181 | z = 0; |
1182 | do { |
1183 | *x1++ = *x << k & 0xffff | z; |
1184 | z = *x++ >> k1; |
1185 | } |
1186 | while(x < xe); |
1187 | if (*x1 = z) |
1188 | ++n1; |
1189 | } |
1190 | #endif |
1191 | else do |
1192 | *x1++ = *x++; |
1193 | while(x < xe); |
1194 | b1->wds = n1 - 1; |
1195 | Bfree(b); |
1196 | return b1; |
1197 | } |
1198 | |
1199 | static int cmp(Bigint *a, Bigint *b) |
1200 | { |
1201 | ULong *xa, *xa0, *xb, *xb0; |
1202 | int i, j; |
1203 | |
1204 | i = a->wds; |
1205 | j = b->wds; |
1206 | #ifdef BSD_QDTOA_DEBUG |
1207 | if (i > 1 && !a->x[i-1]) |
1208 | Bug("cmp called with a->x[a->wds-1] == 0" ); |
1209 | if (j > 1 && !b->x[j-1]) |
1210 | Bug("cmp called with b->x[b->wds-1] == 0" ); |
1211 | #endif |
1212 | if (i -= j) |
1213 | return i; |
1214 | xa0 = a->x; |
1215 | xa = xa0 + j; |
1216 | xb0 = b->x; |
1217 | xb = xb0 + j; |
1218 | for(;;) { |
1219 | if (*--xa != *--xb) |
1220 | return *xa < *xb ? -1 : 1; |
1221 | if (xa <= xa0) |
1222 | break; |
1223 | } |
1224 | return 0; |
1225 | } |
1226 | |
1227 | static Bigint *diff(Bigint *a, Bigint *b) |
1228 | { |
1229 | Bigint *c; |
1230 | int i, wa, wb; |
1231 | Long borrow, y; /* We need signed shifts here. */ |
1232 | ULong *xa, *xae, *xb, *xbe, *xc; |
1233 | #ifdef Pack_32 |
1234 | Long z; |
1235 | #endif |
1236 | |
1237 | i = cmp(a,b); |
1238 | if (!i) { |
1239 | c = Balloc(0); |
1240 | c->wds = 1; |
1241 | c->x[0] = 0; |
1242 | return c; |
1243 | } |
1244 | if (i < 0) { |
1245 | c = a; |
1246 | a = b; |
1247 | b = c; |
1248 | i = 1; |
1249 | } |
1250 | else |
1251 | i = 0; |
1252 | c = Balloc(a->k); |
1253 | c->sign = i; |
1254 | wa = a->wds; |
1255 | xa = a->x; |
1256 | xae = xa + wa; |
1257 | wb = b->wds; |
1258 | xb = b->x; |
1259 | xbe = xb + wb; |
1260 | xc = c->x; |
1261 | borrow = 0; |
1262 | #ifdef Pack_32 |
1263 | do { |
1264 | y = (*xa & 0xffff) - (*xb & 0xffff) + borrow; |
1265 | borrow = y >> 16; |
1266 | Sign_Extend(borrow, y); |
1267 | z = (*xa++ >> 16) - (*xb++ >> 16) + borrow; |
1268 | borrow = z >> 16; |
1269 | Sign_Extend(borrow, z); |
1270 | Storeinc(xc, z, y); |
1271 | } |
1272 | while(xb < xbe); |
1273 | while(xa < xae) { |
1274 | y = (*xa & 0xffff) + borrow; |
1275 | borrow = y >> 16; |
1276 | Sign_Extend(borrow, y); |
1277 | z = (*xa++ >> 16) + borrow; |
1278 | borrow = z >> 16; |
1279 | Sign_Extend(borrow, z); |
1280 | Storeinc(xc, z, y); |
1281 | } |
1282 | #else |
1283 | do { |
1284 | y = *xa++ - *xb++ + borrow; |
1285 | borrow = y >> 16; |
1286 | Sign_Extend(borrow, y); |
1287 | *xc++ = y & 0xffff; |
1288 | } |
1289 | while(xb < xbe); |
1290 | while(xa < xae) { |
1291 | y = *xa++ + borrow; |
1292 | borrow = y >> 16; |
1293 | Sign_Extend(borrow, y); |
1294 | *xc++ = y & 0xffff; |
1295 | } |
1296 | #endif |
1297 | while(!*--xc) |
1298 | wa--; |
1299 | c->wds = wa; |
1300 | return c; |
1301 | } |
1302 | |
1303 | static double ulp(double x) |
1304 | { |
1305 | Long L; |
1306 | double a; |
1307 | |
1308 | L = (getWord0(x) & Exp_mask) - (P-1)*Exp_msk1; |
1309 | #ifndef Sudden_Underflow |
1310 | if (L > 0) { |
1311 | #endif |
1312 | #ifdef IBM |
1313 | L |= Exp_msk1 >> 4; |
1314 | #endif |
1315 | setWord0(&a, L); |
1316 | setWord1(&a, 0); |
1317 | #ifndef Sudden_Underflow |
1318 | } |
1319 | else { |
1320 | L = -L >> Exp_shift; |
1321 | if (L < Exp_shift) { |
1322 | setWord0(&a, 0x80000 >> L); |
1323 | setWord1(&a, 0); |
1324 | } |
1325 | else { |
1326 | setWord0(&a, 0); |
1327 | L -= Exp_shift; |
1328 | setWord1(&a, L >= 31 ? 1U : 1U << (31 - L)); |
1329 | } |
1330 | } |
1331 | #endif |
1332 | return a; |
1333 | } |
1334 | |
1335 | static double b2d(Bigint *a, int *e) |
1336 | { |
1337 | ULong *xa, *xa0, w, y, z; |
1338 | int k; |
1339 | double d; |
1340 | |
1341 | xa0 = a->x; |
1342 | xa = xa0 + a->wds; |
1343 | y = *--xa; |
1344 | #ifdef BSD_QDTOA_DEBUG |
1345 | if (!y) Bug("zero y in b2d" ); |
1346 | #endif |
1347 | k = hi0bits(y); |
1348 | *e = 32 - k; |
1349 | #ifdef Pack_32 |
1350 | if (k < Ebits) { |
1351 | setWord0(&d, Exp_1 | y >> (Ebits - k)); |
1352 | w = xa > xa0 ? *--xa : 0; |
1353 | setWord1(&d, y << ((32-Ebits) + k) | w >> (Ebits - k)); |
1354 | goto ret_d; |
1355 | } |
1356 | z = xa > xa0 ? *--xa : 0; |
1357 | if (k -= Ebits) { |
1358 | setWord0(&d, Exp_1 | y << k | z >> (32 - k)); |
1359 | y = xa > xa0 ? *--xa : 0; |
1360 | setWord1(&d, z << k | y >> (32 - k)); |
1361 | } |
1362 | else { |
1363 | setWord0(&d, Exp_1 | y); |
1364 | setWord1(&d, z); |
1365 | } |
1366 | #else |
1367 | if (k < Ebits + 16) { |
1368 | z = xa > xa0 ? *--xa : 0; |
1369 | setWord0(&d, Exp_1 | y << k - Ebits | z >> Ebits + 16 - k); |
1370 | w = xa > xa0 ? *--xa : 0; |
1371 | y = xa > xa0 ? *--xa : 0; |
1372 | setWord1(&d, z << k + 16 - Ebits | w << k - Ebits | y >> 16 + Ebits - k); |
1373 | goto ret_d; |
1374 | } |
1375 | z = xa > xa0 ? *--xa : 0; |
1376 | w = xa > xa0 ? *--xa : 0; |
1377 | k -= Ebits + 16; |
1378 | setWord0(&d, Exp_1 | y << k + 16 | z << k | w >> 16 - k); |
1379 | y = xa > xa0 ? *--xa : 0; |
1380 | setWord1(&d, w << k + 16 | y << k); |
1381 | #endif |
1382 | ret_d: |
1383 | return d; |
1384 | } |
1385 | |
1386 | static Bigint *d2b(double d, int *e, int *bits) |
1387 | { |
1388 | Bigint *b; |
1389 | int de, i, k; |
1390 | ULong *x, y, z; |
1391 | |
1392 | #ifdef Pack_32 |
1393 | b = Balloc(1); |
1394 | #else |
1395 | b = Balloc(2); |
1396 | #endif |
1397 | x = b->x; |
1398 | |
1399 | z = getWord0(d) & Frac_mask; |
1400 | setWord0(&d, getWord0(d) & 0x7fffffff); /* clear sign bit, which we ignore */ |
1401 | #ifdef Sudden_Underflow |
1402 | de = (int)(getWord0(d) >> Exp_shift); |
1403 | #ifndef IBM |
1404 | z |= Exp_msk11; |
1405 | #endif |
1406 | #else |
1407 | if ((de = int(getWord0(d) >> Exp_shift)) != 0) |
1408 | z |= Exp_msk1; |
1409 | #endif |
1410 | #ifdef Pack_32 |
1411 | if ((y = getWord1(d)) != 0) { |
1412 | if ((k = lo0bits(&y)) != 0) { |
1413 | x[0] = y | z << (32 - k); |
1414 | z >>= k; |
1415 | } |
1416 | else |
1417 | x[0] = y; |
1418 | i = b->wds = (x[1] = z) ? 2 : 1; |
1419 | } |
1420 | else { |
1421 | #ifdef BSD_QDTOA_DEBUG |
1422 | if (!z) |
1423 | Bug("Zero passed to d2b" ); |
1424 | #endif |
1425 | k = lo0bits(&z); |
1426 | x[0] = z; |
1427 | i = b->wds = 1; |
1428 | k += 32; |
1429 | } |
1430 | #else |
1431 | if (y = getWord1(d)) { |
1432 | if (k = lo0bits(&y)) |
1433 | if (k >= 16) { |
1434 | x[0] = y | z << 32 - k & 0xffff; |
1435 | x[1] = z >> k - 16 & 0xffff; |
1436 | x[2] = z >> k; |
1437 | i = 2; |
1438 | } |
1439 | else { |
1440 | x[0] = y & 0xffff; |
1441 | x[1] = y >> 16 | z << 16 - k & 0xffff; |
1442 | x[2] = z >> k & 0xffff; |
1443 | x[3] = z >> k+16; |
1444 | i = 3; |
1445 | } |
1446 | else { |
1447 | x[0] = y & 0xffff; |
1448 | x[1] = y >> 16; |
1449 | x[2] = z & 0xffff; |
1450 | x[3] = z >> 16; |
1451 | i = 3; |
1452 | } |
1453 | } |
1454 | else { |
1455 | #ifdef BSD_QDTOA_DEBUG |
1456 | if (!z) |
1457 | Bug("Zero passed to d2b" ); |
1458 | #endif |
1459 | k = lo0bits(&z); |
1460 | if (k >= 16) { |
1461 | x[0] = z; |
1462 | i = 0; |
1463 | } |
1464 | else { |
1465 | x[0] = z & 0xffff; |
1466 | x[1] = z >> 16; |
1467 | i = 1; |
1468 | } |
1469 | k += 32; |
1470 | } |
1471 | while(!x[i]) |
1472 | --i; |
1473 | b->wds = i + 1; |
1474 | #endif |
1475 | #ifndef Sudden_Underflow |
1476 | if (de) { |
1477 | #endif |
1478 | #ifdef IBM |
1479 | *e = (de - Bias - (P-1) << 2) + k; |
1480 | *bits = 4*P + 8 - k - hi0bits(getWord0(d) & Frac_mask); |
1481 | #else |
1482 | *e = de - Bias - (P-1) + k; |
1483 | *bits = P - k; |
1484 | #endif |
1485 | #ifndef Sudden_Underflow |
1486 | } |
1487 | else { |
1488 | *e = de - Bias - (P-1) + 1 + k; |
1489 | #ifdef Pack_32 |
1490 | *bits = 32*i - hi0bits(x[i-1]); |
1491 | #else |
1492 | *bits = (i+2)*16 - hi0bits(x[i]); |
1493 | #endif |
1494 | } |
1495 | #endif |
1496 | return b; |
1497 | } |
1498 | |
1499 | static double ratio(Bigint *a, Bigint *b) |
1500 | { |
1501 | double da, db; |
1502 | int k, ka, kb; |
1503 | |
1504 | da = b2d(a, &ka); |
1505 | db = b2d(b, &kb); |
1506 | #ifdef Pack_32 |
1507 | k = ka - kb + 32*(a->wds - b->wds); |
1508 | #else |
1509 | k = ka - kb + 16*(a->wds - b->wds); |
1510 | #endif |
1511 | #ifdef IBM |
1512 | if (k > 0) { |
1513 | setWord0(&da, getWord0(da) + (k >> 2)*Exp_msk1); |
1514 | if (k &= 3) |
1515 | da *= 1 << k; |
1516 | } |
1517 | else { |
1518 | k = -k; |
1519 | setWord0(&db, getWord0(db) + (k >> 2)*Exp_msk1); |
1520 | if (k &= 3) |
1521 | db *= 1 << k; |
1522 | } |
1523 | #else |
1524 | if (k > 0) |
1525 | setWord0(&da, getWord0(da) + k*Exp_msk1); |
1526 | else { |
1527 | k = -k; |
1528 | setWord0(&db, getWord0(db) + k*Exp_msk1); |
1529 | } |
1530 | #endif |
1531 | return da / db; |
1532 | } |
1533 | |
1534 | static const double tens[] = { |
1535 | 1e0, 1e1, 1e2, 1e3, 1e4, 1e5, 1e6, 1e7, 1e8, 1e9, |
1536 | 1e10, 1e11, 1e12, 1e13, 1e14, 1e15, 1e16, 1e17, 1e18, 1e19, |
1537 | 1e20, 1e21, 1e22 |
1538 | #ifdef VAX |
1539 | , 1e23, 1e24 |
1540 | #endif |
1541 | }; |
1542 | |
1543 | #ifdef IEEE_Arith |
1544 | static const double bigtens[] = { 1e16, 1e32, 1e64, 1e128, 1e256 }; |
1545 | static const double tinytens[] = { 1e-16, 1e-32, 1e-64, 1e-128, 1e-256 }; |
1546 | #define n_bigtens 5 |
1547 | #else |
1548 | #ifdef IBM |
1549 | static const double bigtens[] = { 1e16, 1e32, 1e64 }; |
1550 | static const double tinytens[] = { 1e-16, 1e-32, 1e-64 }; |
1551 | #define n_bigtens 3 |
1552 | #else |
1553 | static const double bigtens[] = { 1e16, 1e32 }; |
1554 | static const double tinytens[] = { 1e-16, 1e-32 }; |
1555 | #define n_bigtens 2 |
1556 | #endif |
1557 | #endif |
1558 | |
1559 | /* |
1560 | The pre-release gcc3.3 shipped with SuSE 8.2 has a bug which causes |
1561 | the comparison 1e-100 == 0.0 to return true. As a workaround, we |
1562 | compare it to a global variable containing 0.0, which produces |
1563 | correct assembler output. |
1564 | |
1565 | ### consider detecting the broken compilers and using the static |
1566 | ### double for these, and use a #define for all working compilers |
1567 | */ |
1568 | static double g_double_zero = 0.0; |
1569 | |
1570 | Q_CORE_EXPORT double qstrtod(const char *s00, const char **se, bool *ok) |
1571 | { |
1572 | int bb2, bb5, bbe, bd2, bd5, bbbits, bs2, c, dsign, |
1573 | e, e1, esign, i, j, k, nd, nd0, nf, nz, nz0, sign; |
1574 | const char *s, *s0, *s1; |
1575 | double aadj, aadj1, adj, rv, rv0; |
1576 | Long L; |
1577 | ULong y, z; |
1578 | Bigint *bb1, *bd0; |
1579 | Bigint *bb = NULL, *bd = NULL, *bs = NULL, *delta = NULL;/* pacify gcc */ |
1580 | |
1581 | /* |
1582 | #ifndef KR_headers |
1583 | const char decimal_point = localeconv()->decimal_point[0]; |
1584 | #else |
1585 | const char decimal_point = '.'; |
1586 | #endif */ |
1587 | if (ok != 0) |
1588 | *ok = true; |
1589 | |
1590 | const char decimal_point = '.'; |
1591 | |
1592 | sign = nz0 = nz = 0; |
1593 | rv = 0.; |
1594 | |
1595 | |
1596 | for(s = s00; isspace(uchar(*s)); s++) |
1597 | ; |
1598 | |
1599 | if (*s == '-') { |
1600 | sign = 1; |
1601 | s++; |
1602 | } else if (*s == '+') { |
1603 | s++; |
1604 | } |
1605 | |
1606 | if (*s == '\0') { |
1607 | s = s00; |
1608 | goto ret; |
1609 | } |
1610 | |
1611 | if (*s == '0') { |
1612 | nz0 = 1; |
1613 | while(*++s == '0') ; |
1614 | if (!*s) |
1615 | goto ret; |
1616 | } |
1617 | s0 = s; |
1618 | y = z = 0; |
1619 | for(nd = nf = 0; (c = *s) >= '0' && c <= '9'; nd++, s++) |
1620 | if (nd < 9) |
1621 | y = 10*y + c - '0'; |
1622 | else if (nd < 16) |
1623 | z = 10*z + c - '0'; |
1624 | nd0 = nd; |
1625 | if (c == decimal_point) { |
1626 | c = *++s; |
1627 | if (!nd) { |
1628 | for(; c == '0'; c = *++s) |
1629 | nz++; |
1630 | if (c > '0' && c <= '9') { |
1631 | s0 = s; |
1632 | nf += nz; |
1633 | nz = 0; |
1634 | goto have_dig; |
1635 | } |
1636 | goto dig_done; |
1637 | } |
1638 | for(; c >= '0' && c <= '9'; c = *++s) { |
1639 | have_dig: |
1640 | nz++; |
1641 | if (c -= '0') { |
1642 | nf += nz; |
1643 | for(i = 1; i < nz; i++) |
1644 | if (nd++ < 9) |
1645 | y *= 10; |
1646 | else if (nd <= DBL_DIG + 1) |
1647 | z *= 10; |
1648 | if (nd++ < 9) |
1649 | y = 10*y + c; |
1650 | else if (nd <= DBL_DIG + 1) |
1651 | z = 10*z + c; |
1652 | nz = 0; |
1653 | } |
1654 | } |
1655 | } |
1656 | dig_done: |
1657 | e = 0; |
1658 | if (c == 'e' || c == 'E') { |
1659 | if (!nd && !nz && !nz0) { |
1660 | s = s00; |
1661 | goto ret; |
1662 | } |
1663 | s00 = s; |
1664 | esign = 0; |
1665 | switch(c = *++s) { |
1666 | case '-': |
1667 | esign = 1; |
1668 | case '+': |
1669 | c = *++s; |
1670 | } |
1671 | if (c >= '0' && c <= '9') { |
1672 | while(c == '0') |
1673 | c = *++s; |
1674 | if (c > '0' && c <= '9') { |
1675 | L = c - '0'; |
1676 | s1 = s; |
1677 | while((c = *++s) >= '0' && c <= '9') |
1678 | L = 10*L + c - '0'; |
1679 | if (s - s1 > 8 || L > 19999) |
1680 | /* Avoid confusion from exponents |
1681 | * so large that e might overflow. |
1682 | */ |
1683 | e = 19999; /* safe for 16 bit ints */ |
1684 | else |
1685 | e = int(L); |
1686 | if (esign) |
1687 | e = -e; |
1688 | } |
1689 | else |
1690 | e = 0; |
1691 | } |
1692 | else |
1693 | s = s00; |
1694 | } |
1695 | if (!nd) { |
1696 | if (!nz && !nz0) |
1697 | s = s00; |
1698 | goto ret; |
1699 | } |
1700 | e1 = e -= nf; |
1701 | |
1702 | /* Now we have nd0 digits, starting at s0, followed by a |
1703 | * decimal point, followed by nd-nd0 digits. The number we're |
1704 | * after is the integer represented by those digits times |
1705 | * 10**e */ |
1706 | |
1707 | if (!nd0) |
1708 | nd0 = nd; |
1709 | k = nd < DBL_DIG + 1 ? nd : DBL_DIG + 1; |
1710 | rv = y; |
1711 | if (k > 9) |
1712 | #if defined(Q_OS_IRIX) && defined(Q_CC_GNU) |
1713 | { |
1714 | // work around a bug on 64 bit IRIX gcc |
1715 | double *t = (double *) tens; |
1716 | rv = t[k - 9] * rv + z; |
1717 | } |
1718 | #else |
1719 | rv = tens[k - 9] * rv + z; |
1720 | #endif |
1721 | |
1722 | bd0 = 0; |
1723 | if (nd <= DBL_DIG |
1724 | #ifndef RND_PRODQUOT |
1725 | && FLT_ROUNDS == 1 |
1726 | #endif |
1727 | ) { |
1728 | if (!e) |
1729 | goto ret; |
1730 | if (e > 0) { |
1731 | if (e <= Ten_pmax) { |
1732 | #ifdef VAX |
1733 | goto vax_ovfl_check; |
1734 | #else |
1735 | /* rv = */ rounded_product(rv, tens[e]); |
1736 | goto ret; |
1737 | #endif |
1738 | } |
1739 | i = DBL_DIG - nd; |
1740 | if (e <= Ten_pmax + i) { |
1741 | /* A fancier test would sometimes let us do |
1742 | * this for larger i values. |
1743 | */ |
1744 | e -= i; |
1745 | rv *= tens[i]; |
1746 | #ifdef VAX |
1747 | /* VAX exponent range is so narrow we must |
1748 | * worry about overflow here... |
1749 | */ |
1750 | vax_ovfl_check: |
1751 | setWord0(&rv, getWord0(rv) - P*Exp_msk1); |
1752 | /* rv = */ rounded_product(rv, tens[e]); |
1753 | if ((getWord0(rv) & Exp_mask) |
1754 | > Exp_msk1*(DBL_MAX_EXP+Bias-1-P)) |
1755 | goto ovfl; |
1756 | setWord0(&rv, getWord0(rv) + P*Exp_msk1); |
1757 | #else |
1758 | /* rv = */ rounded_product(rv, tens[e]); |
1759 | #endif |
1760 | goto ret; |
1761 | } |
1762 | } |
1763 | #ifndef Inaccurate_Divide |
1764 | else if (e >= -Ten_pmax) { |
1765 | /* rv = */ rounded_quotient(rv, tens[-e]); |
1766 | goto ret; |
1767 | } |
1768 | #endif |
1769 | } |
1770 | e1 += nd - k; |
1771 | |
1772 | /* Get starting approximation = rv * 10**e1 */ |
1773 | |
1774 | if (e1 > 0) { |
1775 | if ((i = e1 & 15) != 0) |
1776 | rv *= tens[i]; |
1777 | if (e1 &= ~15) { |
1778 | if (e1 > DBL_MAX_10_EXP) { |
1779 | ovfl: |
1780 | // errno = ERANGE; |
1781 | if (ok != 0) |
1782 | *ok = false; |
1783 | #ifdef __STDC__ |
1784 | rv = HUGE_VAL; |
1785 | #else |
1786 | /* Can't trust HUGE_VAL */ |
1787 | #ifdef IEEE_Arith |
1788 | setWord0(&rv, Exp_mask); |
1789 | setWord1(&rv, 0); |
1790 | #else |
1791 | setWord0(&rv, Big0); |
1792 | setWord1(&rv, Big1); |
1793 | #endif |
1794 | #endif |
1795 | if (bd0) |
1796 | goto retfree; |
1797 | goto ret; |
1798 | } |
1799 | if (e1 >>= 4) { |
1800 | for(j = 0; e1 > 1; j++, e1 >>= 1) |
1801 | if (e1 & 1) |
1802 | rv *= bigtens[j]; |
1803 | /* The last multiplication could overflow. */ |
1804 | setWord0(&rv, getWord0(rv) - P*Exp_msk1); |
1805 | rv *= bigtens[j]; |
1806 | if ((z = getWord0(rv) & Exp_mask) |
1807 | > Exp_msk1*(DBL_MAX_EXP+Bias-P)) |
1808 | goto ovfl; |
1809 | if (z > Exp_msk1*(DBL_MAX_EXP+Bias-1-P)) { |
1810 | /* set to largest number */ |
1811 | /* (Can't trust DBL_MAX) */ |
1812 | setWord0(&rv, Big0); |
1813 | setWord1(&rv, Big1); |
1814 | } |
1815 | else |
1816 | setWord0(&rv, getWord0(rv) + P*Exp_msk1); |
1817 | } |
1818 | |
1819 | } |
1820 | } |
1821 | else if (e1 < 0) { |
1822 | e1 = -e1; |
1823 | if ((i = e1 & 15) != 0) |
1824 | rv /= tens[i]; |
1825 | if (e1 &= ~15) { |
1826 | e1 >>= 4; |
1827 | if (e1 >= 1 << n_bigtens) |
1828 | goto undfl; |
1829 | for(j = 0; e1 > 1; j++, e1 >>= 1) |
1830 | if (e1 & 1) |
1831 | rv *= tinytens[j]; |
1832 | /* The last multiplication could underflow. */ |
1833 | rv0 = rv; |
1834 | rv *= tinytens[j]; |
1835 | if (rv == g_double_zero) |
1836 | { |
1837 | rv = 2.*rv0; |
1838 | rv *= tinytens[j]; |
1839 | if (rv == g_double_zero) |
1840 | { |
1841 | undfl: |
1842 | rv = 0.; |
1843 | // errno = ERANGE; |
1844 | if (ok != 0) |
1845 | *ok = false; |
1846 | if (bd0) |
1847 | goto retfree; |
1848 | goto ret; |
1849 | } |
1850 | setWord0(&rv, Tiny0); |
1851 | setWord1(&rv, Tiny1); |
1852 | /* The refinement below will clean |
1853 | * this approximation up. |
1854 | */ |
1855 | } |
1856 | } |
1857 | } |
1858 | |
1859 | /* Now the hard part -- adjusting rv to the correct value.*/ |
1860 | |
1861 | /* Put digits into bd: true value = bd * 10^e */ |
1862 | |
1863 | bd0 = s2b(s0, nd0, nd, y); |
1864 | |
1865 | for(;;) { |
1866 | bd = Balloc(bd0->k); |
1867 | Bcopy(bd, bd0); |
1868 | bb = d2b(rv, &bbe, &bbbits); /* rv = bb * 2^bbe */ |
1869 | bs = i2b(1); |
1870 | |
1871 | if (e >= 0) { |
1872 | bb2 = bb5 = 0; |
1873 | bd2 = bd5 = e; |
1874 | } |
1875 | else { |
1876 | bb2 = bb5 = -e; |
1877 | bd2 = bd5 = 0; |
1878 | } |
1879 | if (bbe >= 0) |
1880 | bb2 += bbe; |
1881 | else |
1882 | bd2 -= bbe; |
1883 | bs2 = bb2; |
1884 | #ifdef Sudden_Underflow |
1885 | #ifdef IBM |
1886 | j = 1 + 4*P - 3 - bbbits + ((bbe + bbbits - 1) & 3); |
1887 | #else |
1888 | j = P + 1 - bbbits; |
1889 | #endif |
1890 | #else |
1891 | i = bbe + bbbits - 1; /* logb(rv) */ |
1892 | if (i < Emin) /* denormal */ |
1893 | j = bbe + (P-Emin); |
1894 | else |
1895 | j = P + 1 - bbbits; |
1896 | #endif |
1897 | bb2 += j; |
1898 | bd2 += j; |
1899 | i = bb2 < bd2 ? bb2 : bd2; |
1900 | if (i > bs2) |
1901 | i = bs2; |
1902 | if (i > 0) { |
1903 | bb2 -= i; |
1904 | bd2 -= i; |
1905 | bs2 -= i; |
1906 | } |
1907 | if (bb5 > 0) { |
1908 | bs = pow5mult(bs, bb5); |
1909 | bb1 = mult(bs, bb); |
1910 | Bfree(bb); |
1911 | bb = bb1; |
1912 | } |
1913 | if (bb2 > 0) |
1914 | bb = lshift(bb, bb2); |
1915 | if (bd5 > 0) |
1916 | bd = pow5mult(bd, bd5); |
1917 | if (bd2 > 0) |
1918 | bd = lshift(bd, bd2); |
1919 | if (bs2 > 0) |
1920 | bs = lshift(bs, bs2); |
1921 | delta = diff(bb, bd); |
1922 | dsign = delta->sign; |
1923 | delta->sign = 0; |
1924 | i = cmp(delta, bs); |
1925 | if (i < 0) { |
1926 | /* Error is less than half an ulp -- check for |
1927 | * special case of mantissa a power of two. |
1928 | */ |
1929 | if (dsign || getWord1(rv) || getWord0(rv) & Bndry_mask) |
1930 | break; |
1931 | delta = lshift(delta,Log2P); |
1932 | if (cmp(delta, bs) > 0) |
1933 | goto drop_down; |
1934 | break; |
1935 | } |
1936 | if (i == 0) { |
1937 | /* exactly half-way between */ |
1938 | if (dsign) { |
1939 | if ((getWord0(rv) & Bndry_mask1) == Bndry_mask1 |
1940 | && getWord1(rv) == 0xffffffff) { |
1941 | /*boundary case -- increment exponent*/ |
1942 | setWord0(&rv, (getWord0(rv) & Exp_mask) |
1943 | + Exp_msk1 |
1944 | #ifdef IBM |
1945 | | Exp_msk1 >> 4 |
1946 | #endif |
1947 | ); |
1948 | setWord1(&rv, 0); |
1949 | break; |
1950 | } |
1951 | } |
1952 | else if (!(getWord0(rv) & Bndry_mask) && !getWord1(rv)) { |
1953 | drop_down: |
1954 | /* boundary case -- decrement exponent */ |
1955 | #ifdef Sudden_Underflow |
1956 | L = getWord0(rv) & Exp_mask; |
1957 | #ifdef IBM |
1958 | if (L < Exp_msk1) |
1959 | #else |
1960 | if (L <= Exp_msk1) |
1961 | #endif |
1962 | goto undfl; |
1963 | L -= Exp_msk1; |
1964 | #else |
1965 | L = (getWord0(rv) & Exp_mask) - Exp_msk1; |
1966 | #endif |
1967 | setWord0(&rv, L | Bndry_mask1); |
1968 | setWord1(&rv, 0xffffffff); |
1969 | #ifdef IBM |
1970 | goto cont; |
1971 | #else |
1972 | break; |
1973 | #endif |
1974 | } |
1975 | #ifndef ROUND_BIASED |
1976 | if (!(getWord1(rv) & LSB)) |
1977 | break; |
1978 | #endif |
1979 | if (dsign) |
1980 | rv += ulp(rv); |
1981 | #ifndef ROUND_BIASED |
1982 | else { |
1983 | rv -= ulp(rv); |
1984 | #ifndef Sudden_Underflow |
1985 | if (rv == g_double_zero) |
1986 | goto undfl; |
1987 | #endif |
1988 | } |
1989 | #endif |
1990 | break; |
1991 | } |
1992 | if ((aadj = ratio(delta, bs)) <= 2.) { |
1993 | if (dsign) |
1994 | aadj = aadj1 = 1.; |
1995 | else if (getWord1(rv) || getWord0(rv) & Bndry_mask) { |
1996 | #ifndef Sudden_Underflow |
1997 | if (getWord1(rv) == Tiny1 && !getWord0(rv)) |
1998 | goto undfl; |
1999 | #endif |
2000 | aadj = 1.; |
2001 | aadj1 = -1.; |
2002 | } |
2003 | else { |
2004 | /* special case -- power of FLT_RADIX to be */ |
2005 | /* rounded down... */ |
2006 | |
2007 | if (aadj < 2./FLT_RADIX) |
2008 | aadj = 1./FLT_RADIX; |
2009 | else |
2010 | aadj *= 0.5; |
2011 | aadj1 = -aadj; |
2012 | } |
2013 | } |
2014 | else { |
2015 | aadj *= 0.5; |
2016 | aadj1 = dsign ? aadj : -aadj; |
2017 | #ifdef Check_FLT_ROUNDS |
2018 | switch(FLT_ROUNDS) { |
2019 | case 2: /* towards +infinity */ |
2020 | aadj1 -= 0.5; |
2021 | break; |
2022 | case 0: /* towards 0 */ |
2023 | case 3: /* towards -infinity */ |
2024 | aadj1 += 0.5; |
2025 | } |
2026 | #else |
2027 | if (FLT_ROUNDS == 0) |
2028 | aadj1 += 0.5; |
2029 | #endif |
2030 | } |
2031 | y = getWord0(rv) & Exp_mask; |
2032 | |
2033 | /* Check for overflow */ |
2034 | |
2035 | if (y == Exp_msk1*(DBL_MAX_EXP+Bias-1)) { |
2036 | rv0 = rv; |
2037 | setWord0(&rv, getWord0(rv) - P*Exp_msk1); |
2038 | adj = aadj1 * ulp(rv); |
2039 | rv += adj; |
2040 | if ((getWord0(rv) & Exp_mask) >= |
2041 | Exp_msk1*(DBL_MAX_EXP+Bias-P)) { |
2042 | if (getWord0(rv0) == Big0 && getWord1(rv0) == Big1) |
2043 | goto ovfl; |
2044 | setWord0(&rv, Big0); |
2045 | setWord1(&rv, Big1); |
2046 | goto cont; |
2047 | } |
2048 | else |
2049 | setWord0(&rv, getWord0(rv) + P*Exp_msk1); |
2050 | } |
2051 | else { |
2052 | #ifdef Sudden_Underflow |
2053 | if ((getWord0(rv) & Exp_mask) <= P*Exp_msk1) { |
2054 | rv0 = rv; |
2055 | setWord0(&rv, getWord0(rv) + P*Exp_msk1); |
2056 | adj = aadj1 * ulp(rv); |
2057 | rv += adj; |
2058 | #ifdef IBM |
2059 | if ((getWord0(rv) & Exp_mask) < P*Exp_msk1) |
2060 | #else |
2061 | if ((getWord0(rv) & Exp_mask) <= P*Exp_msk1) |
2062 | #endif |
2063 | { |
2064 | if (getWord0(rv0) == Tiny0 |
2065 | && getWord1(rv0) == Tiny1) |
2066 | goto undfl; |
2067 | setWord0(&rv, Tiny0); |
2068 | setWord1(&rv, Tiny1); |
2069 | goto cont; |
2070 | } |
2071 | else |
2072 | setWord0(&rv, getWord0(rv) - P*Exp_msk1); |
2073 | } |
2074 | else { |
2075 | adj = aadj1 * ulp(rv); |
2076 | rv += adj; |
2077 | } |
2078 | #else |
2079 | /* Compute adj so that the IEEE rounding rules will |
2080 | * correctly round rv + adj in some half-way cases. |
2081 | * If rv * ulp(rv) is denormalized (i.e., |
2082 | * y <= (P-1)*Exp_msk1), we must adjust aadj to avoid |
2083 | * trouble from bits lost to denormalization; |
2084 | * example: 1.2e-307 . |
2085 | */ |
2086 | if (y <= (P-1)*Exp_msk1 && aadj >= 1.) { |
2087 | aadj1 = int(aadj + 0.5); |
2088 | if (!dsign) |
2089 | aadj1 = -aadj1; |
2090 | } |
2091 | adj = aadj1 * ulp(rv); |
2092 | rv += adj; |
2093 | #endif |
2094 | } |
2095 | z = getWord0(rv) & Exp_mask; |
2096 | if (y == z) { |
2097 | /* Can we stop now? */ |
2098 | L = Long(aadj); |
2099 | aadj -= L; |
2100 | /* The tolerances below are conservative. */ |
2101 | if (dsign || getWord1(rv) || getWord0(rv) & Bndry_mask) { |
2102 | if (aadj < .4999999 || aadj > .5000001) |
2103 | break; |
2104 | } |
2105 | else if (aadj < .4999999/FLT_RADIX) |
2106 | break; |
2107 | } |
2108 | cont: |
2109 | Bfree(bb); |
2110 | Bfree(bd); |
2111 | Bfree(bs); |
2112 | Bfree(delta); |
2113 | } |
2114 | retfree: |
2115 | Bfree(bb); |
2116 | Bfree(bd); |
2117 | Bfree(bs); |
2118 | Bfree(bd0); |
2119 | Bfree(delta); |
2120 | ret: |
2121 | if (se) |
2122 | *se = s; |
2123 | return sign ? -rv : rv; |
2124 | } |
2125 | |
2126 | static int quorem(Bigint *b, Bigint *S) |
2127 | { |
2128 | int n; |
2129 | Long borrow, y; |
2130 | ULong carry, q, ys; |
2131 | ULong *bx, *bxe, *sx, *sxe; |
2132 | #ifdef Pack_32 |
2133 | Long z; |
2134 | ULong si, zs; |
2135 | #endif |
2136 | |
2137 | n = S->wds; |
2138 | #ifdef BSD_QDTOA_DEBUG |
2139 | /*debug*/ if (b->wds > n) |
2140 | /*debug*/ Bug("oversize b in quorem" ); |
2141 | #endif |
2142 | if (b->wds < n) |
2143 | return 0; |
2144 | sx = S->x; |
2145 | sxe = sx + --n; |
2146 | bx = b->x; |
2147 | bxe = bx + n; |
2148 | q = *bxe / (*sxe + 1); /* ensure q <= true quotient */ |
2149 | #ifdef BSD_QDTOA_DEBUG |
2150 | /*debug*/ if (q > 9) |
2151 | /*debug*/ Bug("oversized quotient in quorem" ); |
2152 | #endif |
2153 | if (q) { |
2154 | borrow = 0; |
2155 | carry = 0; |
2156 | do { |
2157 | #ifdef Pack_32 |
2158 | si = *sx++; |
2159 | ys = (si & 0xffff) * q + carry; |
2160 | zs = (si >> 16) * q + (ys >> 16); |
2161 | carry = zs >> 16; |
2162 | y = (*bx & 0xffff) - (ys & 0xffff) + borrow; |
2163 | borrow = y >> 16; |
2164 | Sign_Extend(borrow, y); |
2165 | z = (*bx >> 16) - (zs & 0xffff) + borrow; |
2166 | borrow = z >> 16; |
2167 | Sign_Extend(borrow, z); |
2168 | Storeinc(bx, z, y); |
2169 | #else |
2170 | ys = *sx++ * q + carry; |
2171 | carry = ys >> 16; |
2172 | y = *bx - (ys & 0xffff) + borrow; |
2173 | borrow = y >> 16; |
2174 | Sign_Extend(borrow, y); |
2175 | *bx++ = y & 0xffff; |
2176 | #endif |
2177 | } |
2178 | while(sx <= sxe); |
2179 | if (!*bxe) { |
2180 | bx = b->x; |
2181 | while(--bxe > bx && !*bxe) |
2182 | --n; |
2183 | b->wds = n; |
2184 | } |
2185 | } |
2186 | if (cmp(b, S) >= 0) { |
2187 | q++; |
2188 | borrow = 0; |
2189 | carry = 0; |
2190 | bx = b->x; |
2191 | sx = S->x; |
2192 | do { |
2193 | #ifdef Pack_32 |
2194 | si = *sx++; |
2195 | ys = (si & 0xffff) + carry; |
2196 | zs = (si >> 16) + (ys >> 16); |
2197 | carry = zs >> 16; |
2198 | y = (*bx & 0xffff) - (ys & 0xffff) + borrow; |
2199 | borrow = y >> 16; |
2200 | Sign_Extend(borrow, y); |
2201 | z = (*bx >> 16) - (zs & 0xffff) + borrow; |
2202 | borrow = z >> 16; |
2203 | Sign_Extend(borrow, z); |
2204 | Storeinc(bx, z, y); |
2205 | #else |
2206 | ys = *sx++ + carry; |
2207 | carry = ys >> 16; |
2208 | y = *bx - (ys & 0xffff) + borrow; |
2209 | borrow = y >> 16; |
2210 | Sign_Extend(borrow, y); |
2211 | *bx++ = y & 0xffff; |
2212 | #endif |
2213 | } |
2214 | while(sx <= sxe); |
2215 | bx = b->x; |
2216 | bxe = bx + n; |
2217 | if (!*bxe) { |
2218 | while(--bxe > bx && !*bxe) |
2219 | --n; |
2220 | b->wds = n; |
2221 | } |
2222 | } |
2223 | return q; |
2224 | } |
2225 | |
2226 | /* dtoa for IEEE arithmetic (dmg): convert double to ASCII string. |
2227 | * |
2228 | * Inspired by "How to Print Floating-Point Numbers Accurately" by |
2229 | * Guy L. Steele, Jr. and Jon L. White [Proc. ACM SIGPLAN '90, pp. 92-101]. |
2230 | * |
2231 | * Modifications: |
2232 | * 1. Rather than iterating, we use a simple numeric overestimate |
2233 | * to determine k = floor(log10(d)). We scale relevant |
2234 | * quantities using O(log2(k)) rather than O(k) multiplications. |
2235 | * 2. For some modes > 2 (corresponding to ecvt and fcvt), we don't |
2236 | * try to generate digits strictly left to right. Instead, we |
2237 | * compute with fewer bits and propagate the carry if necessary |
2238 | * when rounding the final digit up. This is often faster. |
2239 | * 3. Under the assumption that input will be rounded nearest, |
2240 | * mode 0 renders 1e23 as 1e23 rather than 9.999999999999999e22. |
2241 | * That is, we allow equality in stopping tests when the |
2242 | * round-nearest rule will give the same floating-point value |
2243 | * as would satisfaction of the stopping test with strict |
2244 | * inequality. |
2245 | * 4. We remove common factors of powers of 2 from relevant |
2246 | * quantities. |
2247 | * 5. When converting floating-point integers less than 1e16, |
2248 | * we use floating-point arithmetic rather than resorting |
2249 | * to multiple-precision integers. |
2250 | * 6. When asked to produce fewer than 15 digits, we first try |
2251 | * to get by with floating-point arithmetic; we resort to |
2252 | * multiple-precision integer arithmetic only if we cannot |
2253 | * guarantee that the floating-point calculation has given |
2254 | * the correctly rounded result. For k requested digits and |
2255 | * "uniformly" distributed input, the probability is |
2256 | * something like 10^(k-15) that we must resort to the Long |
2257 | * calculation. |
2258 | */ |
2259 | |
2260 | #if defined(Q_OS_WIN) && defined (Q_CC_GNU) && !defined(_clear87) // See QTBUG-7576 |
2261 | extern "C" { |
2262 | __attribute__ ((dllimport)) unsigned int __cdecl __MINGW_NOTHROW _control87 (unsigned int unNew, unsigned int unMask); |
2263 | __attribute__ ((dllimport)) unsigned int __cdecl __MINGW_NOTHROW _clearfp (void); /* Clear the FPU status word */ |
2264 | } |
2265 | # define _clear87 _clearfp |
2266 | #endif |
2267 | |
2268 | /* This actually sometimes returns a pointer to a string literal |
2269 | cast to a char*. Do NOT try to modify the return value. */ |
2270 | |
2271 | Q_CORE_EXPORT char *qdtoa ( double d, int mode, int ndigits, int *decpt, int *sign, char **rve, char **resultp) |
2272 | { |
2273 | // Some values of the floating-point control word can cause _qdtoa to crash with an underflow. |
2274 | // We set a safe value here. |
2275 | #ifdef Q_OS_WIN |
2276 | _clear87(); |
2277 | unsigned int oldbits = _control87(0, 0); |
2278 | #ifndef MCW_EM |
2279 | # ifdef _MCW_EM |
2280 | # define MCW_EM _MCW_EM |
2281 | # else |
2282 | # define MCW_EM 0x0008001F |
2283 | # endif |
2284 | #endif |
2285 | _control87(MCW_EM, MCW_EM); |
2286 | #endif |
2287 | |
2288 | #if defined(Q_OS_LINUX) && !defined(__UCLIBC__) |
2289 | fenv_t envp; |
2290 | feholdexcept(&envp); |
2291 | #endif |
2292 | |
2293 | char *s = _qdtoa(d, mode, ndigits, decpt, sign, rve, resultp); |
2294 | |
2295 | #ifdef Q_OS_WIN |
2296 | _clear87(); |
2297 | #ifndef _M_X64 |
2298 | _control87(oldbits, 0xFFFFF); |
2299 | #else |
2300 | _control87(oldbits, _MCW_EM|_MCW_DN|_MCW_RC); |
2301 | #endif //_M_X64 |
2302 | #endif //Q_OS_WIN |
2303 | |
2304 | #if defined(Q_OS_LINUX) && !defined(__UCLIBC__) |
2305 | fesetenv(&envp); |
2306 | #endif |
2307 | |
2308 | return s; |
2309 | } |
2310 | |
2311 | static char *_qdtoa( NEEDS_VOLATILE double d, int mode, int ndigits, int *decpt, int *sign, char **rve, char **resultp) |
2312 | { |
2313 | /* |
2314 | Arguments ndigits, decpt, sign are similar to those |
2315 | of ecvt and fcvt; trailing zeros are suppressed from |
2316 | the returned string. If not null, *rve is set to point |
2317 | to the end of the return value. If d is +-Infinity or NaN, |
2318 | then *decpt is set to 9999. |
2319 | |
2320 | mode: |
2321 | 0 ==> shortest string that yields d when read in |
2322 | and rounded to nearest. |
2323 | 1 ==> like 0, but with Steele & White stopping rule; |
2324 | e.g. with IEEE P754 arithmetic , mode 0 gives |
2325 | 1e23 whereas mode 1 gives 9.999999999999999e22. |
2326 | 2 ==> max(1,ndigits) significant digits. This gives a |
2327 | return value similar to that of ecvt, except |
2328 | that trailing zeros are suppressed. |
2329 | 3 ==> through ndigits past the decimal point. This |
2330 | gives a return value similar to that from fcvt, |
2331 | except that trailing zeros are suppressed, and |
2332 | ndigits can be negative. |
2333 | 4-9 should give the same return values as 2-3, i.e., |
2334 | 4 <= mode <= 9 ==> same return as mode |
2335 | 2 + (mode & 1). These modes are mainly for |
2336 | debugging; often they run slower but sometimes |
2337 | faster than modes 2-3. |
2338 | 4,5,8,9 ==> left-to-right digit generation. |
2339 | 6-9 ==> don't try fast floating-point estimate |
2340 | (if applicable). |
2341 | |
2342 | Values of mode other than 0-9 are treated as mode 0. |
2343 | |
2344 | Sufficient space is allocated to the return value |
2345 | to hold the suppressed trailing zeros. |
2346 | */ |
2347 | |
2348 | int bbits, b2, b5, be, dig, i, ieps, ilim0, |
2349 | j, j1, k, k0, k_check, leftright, m2, m5, s2, s5, |
2350 | try_quick; |
2351 | int ilim = 0, ilim1 = 0, spec_case = 0; /* pacify gcc */ |
2352 | Long L; |
2353 | #ifndef Sudden_Underflow |
2354 | int denorm; |
2355 | ULong x; |
2356 | #endif |
2357 | Bigint *b, *b1, *delta, *mhi, *S; |
2358 | Bigint *mlo = NULL; /* pacify gcc */ |
2359 | double d2; |
2360 | double ds, eps; |
2361 | char *s, *s0; |
2362 | |
2363 | if (getWord0(d) & Sign_bit) { |
2364 | /* set sign for everything, including 0's and NaNs */ |
2365 | *sign = 1; |
2366 | setWord0(&d, getWord0(d) & ~Sign_bit); /* clear sign bit */ |
2367 | } |
2368 | else |
2369 | *sign = 0; |
2370 | |
2371 | #if defined(IEEE_Arith) + defined(VAX) |
2372 | #ifdef IEEE_Arith |
2373 | if ((getWord0(d) & Exp_mask) == Exp_mask) |
2374 | #else |
2375 | if (getWord0(d) == 0x8000) |
2376 | #endif |
2377 | { |
2378 | /* Infinity or NaN */ |
2379 | *decpt = 9999; |
2380 | s = |
2381 | #ifdef IEEE_Arith |
2382 | !getWord1(d) && !(getWord0(d) & 0xfffff) ? const_cast<char*>("Infinity" ) : |
2383 | #endif |
2384 | const_cast<char*>("NaN" ); |
2385 | if (rve) |
2386 | *rve = |
2387 | #ifdef IEEE_Arith |
2388 | s[3] ? s + 8 : |
2389 | #endif |
2390 | s + 3; |
2391 | return s; |
2392 | } |
2393 | #endif |
2394 | #ifdef IBM |
2395 | d += 0; /* normalize */ |
2396 | #endif |
2397 | if (d == g_double_zero) |
2398 | { |
2399 | *decpt = 1; |
2400 | s = const_cast<char*>("0" ); |
2401 | if (rve) |
2402 | *rve = s + 1; |
2403 | return s; |
2404 | } |
2405 | |
2406 | b = d2b(d, &be, &bbits); |
2407 | #ifdef Sudden_Underflow |
2408 | i = (int)(getWord0(d) >> Exp_shift1 & (Exp_mask>>Exp_shift1)); |
2409 | #else |
2410 | if ((i = int(getWord0(d) >> Exp_shift1 & (Exp_mask>>Exp_shift1))) != 0) { |
2411 | #endif |
2412 | d2 = d; |
2413 | setWord0(&d2, getWord0(d2) & Frac_mask1); |
2414 | setWord0(&d2, getWord0(d2) | Exp_11); |
2415 | #ifdef IBM |
2416 | if (j = 11 - hi0bits(getWord0(d2) & Frac_mask)) |
2417 | d2 /= 1 << j; |
2418 | #endif |
2419 | |
2420 | /* log(x) ~=~ log(1.5) + (x-1.5)/1.5 |
2421 | * log10(x) = log(x) / log(10) |
2422 | * ~=~ log(1.5)/log(10) + (x-1.5)/(1.5*log(10)) |
2423 | * log10(d) = (i-Bias)*log(2)/log(10) + log10(d2) |
2424 | * |
2425 | * This suggests computing an approximation k to log10(d) by |
2426 | * |
2427 | * k = (i - Bias)*0.301029995663981 |
2428 | * + ( (d2-1.5)*0.289529654602168 + 0.176091259055681 ); |
2429 | * |
2430 | * We want k to be too large rather than too small. |
2431 | * The error in the first-order Taylor series approximation |
2432 | * is in our favor, so we just round up the constant enough |
2433 | * to compensate for any error in the multiplication of |
2434 | * (i - Bias) by 0.301029995663981; since |i - Bias| <= 1077, |
2435 | * and 1077 * 0.30103 * 2^-52 ~=~ 7.2e-14, |
2436 | * adding 1e-13 to the constant term more than suffices. |
2437 | * Hence we adjust the constant term to 0.1760912590558. |
2438 | * (We could get a more accurate k by invoking log10, |
2439 | * but this is probably not worthwhile.) |
2440 | */ |
2441 | |
2442 | i -= Bias; |
2443 | #ifdef IBM |
2444 | i <<= 2; |
2445 | i += j; |
2446 | #endif |
2447 | #ifndef Sudden_Underflow |
2448 | denorm = 0; |
2449 | } |
2450 | else { |
2451 | /* d is denormalized */ |
2452 | |
2453 | i = bbits + be + (Bias + (P-1) - 1); |
2454 | x = i > 32 ? getWord0(d) << (64 - i) | getWord1(d) >> (i - 32) |
2455 | : getWord1(d) << (32 - i); |
2456 | d2 = x; |
2457 | setWord0(&d2, getWord0(d2) - 31*Exp_msk1); /* adjust exponent */ |
2458 | i -= (Bias + (P-1) - 1) + 1; |
2459 | denorm = 1; |
2460 | } |
2461 | #endif |
2462 | ds = (d2-1.5)*0.289529654602168 + 0.1760912590558 + i*0.301029995663981; |
2463 | k = int(ds); |
2464 | if (ds < 0. && ds != k) |
2465 | k--; /* want k = floor(ds) */ |
2466 | k_check = 1; |
2467 | if (k >= 0 && k <= Ten_pmax) { |
2468 | if (d < tens[k]) |
2469 | k--; |
2470 | k_check = 0; |
2471 | } |
2472 | j = bbits - i - 1; |
2473 | if (j >= 0) { |
2474 | b2 = 0; |
2475 | s2 = j; |
2476 | } |
2477 | else { |
2478 | b2 = -j; |
2479 | s2 = 0; |
2480 | } |
2481 | if (k >= 0) { |
2482 | b5 = 0; |
2483 | s5 = k; |
2484 | s2 += k; |
2485 | } |
2486 | else { |
2487 | b2 -= k; |
2488 | b5 = -k; |
2489 | s5 = 0; |
2490 | } |
2491 | if (mode < 0 || mode > 9) |
2492 | mode = 0; |
2493 | try_quick = 1; |
2494 | if (mode > 5) { |
2495 | mode -= 4; |
2496 | try_quick = 0; |
2497 | } |
2498 | leftright = 1; |
2499 | switch(mode) { |
2500 | case 0: |
2501 | case 1: |
2502 | ilim = ilim1 = -1; |
2503 | i = 18; |
2504 | ndigits = 0; |
2505 | break; |
2506 | case 2: |
2507 | leftright = 0; |
2508 | /* no break */ |
2509 | case 4: |
2510 | if (ndigits <= 0) |
2511 | ndigits = 1; |
2512 | ilim = ilim1 = i = ndigits; |
2513 | break; |
2514 | case 3: |
2515 | leftright = 0; |
2516 | /* no break */ |
2517 | case 5: |
2518 | i = ndigits + k + 1; |
2519 | ilim = i; |
2520 | ilim1 = i - 1; |
2521 | if (i <= 0) |
2522 | i = 1; |
2523 | } |
2524 | QT_TRY { |
2525 | *resultp = static_cast<char *>(malloc(i + 1)); |
2526 | Q_CHECK_PTR(*resultp); |
2527 | } QT_CATCH(...) { |
2528 | Bfree(b); |
2529 | QT_RETHROW; |
2530 | } |
2531 | s = s0 = *resultp; |
2532 | |
2533 | if (ilim >= 0 && ilim <= Quick_max && try_quick) { |
2534 | |
2535 | /* Try to get by with floating-point arithmetic. */ |
2536 | |
2537 | i = 0; |
2538 | d2 = d; |
2539 | k0 = k; |
2540 | ilim0 = ilim; |
2541 | ieps = 2; /* conservative */ |
2542 | if (k > 0) { |
2543 | ds = tens[k&0xf]; |
2544 | j = k >> 4; |
2545 | if (j & Bletch) { |
2546 | /* prevent overflows */ |
2547 | j &= Bletch - 1; |
2548 | d /= bigtens[n_bigtens-1]; |
2549 | ieps++; |
2550 | } |
2551 | for(; j; j >>= 1, i++) |
2552 | if (j & 1) { |
2553 | ieps++; |
2554 | ds *= bigtens[i]; |
2555 | } |
2556 | d /= ds; |
2557 | } |
2558 | else if ((j1 = -k) != 0) { |
2559 | d *= tens[j1 & 0xf]; |
2560 | for(j = j1 >> 4; j; j >>= 1, i++) |
2561 | if (j & 1) { |
2562 | ieps++; |
2563 | d *= bigtens[i]; |
2564 | } |
2565 | } |
2566 | if (k_check && d < 1. && ilim > 0) { |
2567 | if (ilim1 <= 0) |
2568 | goto fast_failed; |
2569 | ilim = ilim1; |
2570 | k--; |
2571 | d *= 10.; |
2572 | ieps++; |
2573 | } |
2574 | eps = ieps*d + 7.; |
2575 | setWord0(&eps, getWord0(eps) - (P-1)*Exp_msk1); |
2576 | if (ilim == 0) { |
2577 | S = mhi = 0; |
2578 | d -= 5.; |
2579 | if (d > eps) |
2580 | goto one_digit; |
2581 | if (d < -eps) |
2582 | goto no_digits; |
2583 | goto fast_failed; |
2584 | } |
2585 | #ifndef No_leftright |
2586 | if (leftright) { |
2587 | /* Use Steele & White method of only |
2588 | * generating digits needed. |
2589 | */ |
2590 | eps = 0.5/tens[ilim-1] - eps; |
2591 | for(i = 0;;) { |
2592 | L = Long(d); |
2593 | d -= L; |
2594 | *s++ = '0' + int(L); |
2595 | if (d < eps) |
2596 | goto ret1; |
2597 | if (1. - d < eps) |
2598 | goto bump_up; |
2599 | if (++i >= ilim) |
2600 | break; |
2601 | eps *= 10.; |
2602 | d *= 10.; |
2603 | } |
2604 | } |
2605 | else { |
2606 | #endif |
2607 | /* Generate ilim digits, then fix them up. */ |
2608 | #if defined(Q_OS_IRIX) && defined(Q_CC_GNU) |
2609 | // work around a bug on 64 bit IRIX gcc |
2610 | double *t = (double *) tens; |
2611 | eps *= t[ilim-1]; |
2612 | #else |
2613 | eps *= tens[ilim-1]; |
2614 | #endif |
2615 | for(i = 1;; i++, d *= 10.) { |
2616 | L = Long(d); |
2617 | d -= L; |
2618 | *s++ = '0' + int(L); |
2619 | if (i == ilim) { |
2620 | if (d > 0.5 + eps) |
2621 | goto bump_up; |
2622 | else if (d < 0.5 - eps) { |
2623 | while(*--s == '0') {} |
2624 | s++; |
2625 | goto ret1; |
2626 | } |
2627 | break; |
2628 | } |
2629 | } |
2630 | #ifndef No_leftright |
2631 | } |
2632 | #endif |
2633 | fast_failed: |
2634 | s = s0; |
2635 | d = d2; |
2636 | k = k0; |
2637 | ilim = ilim0; |
2638 | } |
2639 | |
2640 | /* Do we have a "small" integer? */ |
2641 | |
2642 | if (be >= 0 && k <= Int_max) { |
2643 | /* Yes. */ |
2644 | ds = tens[k]; |
2645 | if (ndigits < 0 && ilim <= 0) { |
2646 | S = mhi = 0; |
2647 | if (ilim < 0 || d <= 5*ds) |
2648 | goto no_digits; |
2649 | goto one_digit; |
2650 | } |
2651 | for(i = 1;; i++) { |
2652 | L = Long(d / ds); |
2653 | d -= L*ds; |
2654 | #ifdef Check_FLT_ROUNDS |
2655 | /* If FLT_ROUNDS == 2, L will usually be high by 1 */ |
2656 | if (d < 0) { |
2657 | L--; |
2658 | d += ds; |
2659 | } |
2660 | #endif |
2661 | *s++ = '0' + int(L); |
2662 | if (i == ilim) { |
2663 | d += d; |
2664 | if (d > ds || (d == ds && L & 1)) { |
2665 | bump_up: |
2666 | while(*--s == '9') |
2667 | if (s == s0) { |
2668 | k++; |
2669 | *s = '0'; |
2670 | break; |
2671 | } |
2672 | ++*s++; |
2673 | } |
2674 | break; |
2675 | } |
2676 | if ((d *= 10.) == g_double_zero) |
2677 | break; |
2678 | } |
2679 | goto ret1; |
2680 | } |
2681 | |
2682 | m2 = b2; |
2683 | m5 = b5; |
2684 | mhi = mlo = 0; |
2685 | if (leftright) { |
2686 | if (mode < 2) { |
2687 | i = |
2688 | #ifndef Sudden_Underflow |
2689 | denorm ? be + (Bias + (P-1) - 1 + 1) : |
2690 | #endif |
2691 | #ifdef IBM |
2692 | 1 + 4*P - 3 - bbits + ((bbits + be - 1) & 3); |
2693 | #else |
2694 | 1 + P - bbits; |
2695 | #endif |
2696 | } |
2697 | else { |
2698 | j = ilim - 1; |
2699 | if (m5 >= j) |
2700 | m5 -= j; |
2701 | else { |
2702 | s5 += j -= m5; |
2703 | b5 += j; |
2704 | m5 = 0; |
2705 | } |
2706 | if ((i = ilim) < 0) { |
2707 | m2 -= i; |
2708 | i = 0; |
2709 | } |
2710 | } |
2711 | b2 += i; |
2712 | s2 += i; |
2713 | mhi = i2b(1); |
2714 | } |
2715 | if (m2 > 0 && s2 > 0) { |
2716 | i = m2 < s2 ? m2 : s2; |
2717 | b2 -= i; |
2718 | m2 -= i; |
2719 | s2 -= i; |
2720 | } |
2721 | if (b5 > 0) { |
2722 | if (leftright) { |
2723 | if (m5 > 0) { |
2724 | mhi = pow5mult(mhi, m5); |
2725 | b1 = mult(mhi, b); |
2726 | Bfree(b); |
2727 | b = b1; |
2728 | } |
2729 | if ((j = b5 - m5) != 0) |
2730 | b = pow5mult(b, j); |
2731 | } |
2732 | else |
2733 | b = pow5mult(b, b5); |
2734 | } |
2735 | S = i2b(1); |
2736 | if (s5 > 0) |
2737 | S = pow5mult(S, s5); |
2738 | |
2739 | /* Check for special case that d is a normalized power of 2. */ |
2740 | |
2741 | if (mode < 2) { |
2742 | if (!getWord1(d) && !(getWord0(d) & Bndry_mask) |
2743 | #ifndef Sudden_Underflow |
2744 | && getWord0(d) & Exp_mask |
2745 | #endif |
2746 | ) { |
2747 | /* The special case */ |
2748 | b2 += Log2P; |
2749 | s2 += Log2P; |
2750 | spec_case = 1; |
2751 | } |
2752 | else |
2753 | spec_case = 0; |
2754 | } |
2755 | |
2756 | /* Arrange for convenient computation of quotients: |
2757 | * shift left if necessary so divisor has 4 leading 0 bits. |
2758 | * |
2759 | * Perhaps we should just compute leading 28 bits of S once |
2760 | * and for all and pass them and a shift to quorem, so it |
2761 | * can do shifts and ors to compute the numerator for q. |
2762 | */ |
2763 | #ifdef Pack_32 |
2764 | if ((i = ((s5 ? 32 - hi0bits(S->x[S->wds-1]) : 1) + s2) & 0x1f) != 0) |
2765 | i = 32 - i; |
2766 | #else |
2767 | if (i = ((s5 ? 32 - hi0bits(S->x[S->wds-1]) : 1) + s2) & 0xf) |
2768 | i = 16 - i; |
2769 | #endif |
2770 | if (i > 4) { |
2771 | i -= 4; |
2772 | b2 += i; |
2773 | m2 += i; |
2774 | s2 += i; |
2775 | } |
2776 | else if (i < 4) { |
2777 | i += 28; |
2778 | b2 += i; |
2779 | m2 += i; |
2780 | s2 += i; |
2781 | } |
2782 | if (b2 > 0) |
2783 | b = lshift(b, b2); |
2784 | if (s2 > 0) |
2785 | S = lshift(S, s2); |
2786 | if (k_check) { |
2787 | if (cmp(b,S) < 0) { |
2788 | k--; |
2789 | b = multadd(b, 10, 0); /* we botched the k estimate */ |
2790 | if (leftright) |
2791 | mhi = multadd(mhi, 10, 0); |
2792 | ilim = ilim1; |
2793 | } |
2794 | } |
2795 | if (ilim <= 0 && mode > 2) { |
2796 | if (ilim < 0 || cmp(b,S = multadd(S,5,0)) <= 0) { |
2797 | /* no digits, fcvt style */ |
2798 | no_digits: |
2799 | k = -1 - ndigits; |
2800 | goto ret; |
2801 | } |
2802 | one_digit: |
2803 | *s++ = '1'; |
2804 | k++; |
2805 | goto ret; |
2806 | } |
2807 | if (leftright) { |
2808 | if (m2 > 0) |
2809 | mhi = lshift(mhi, m2); |
2810 | |
2811 | /* Compute mlo -- check for special case |
2812 | * that d is a normalized power of 2. |
2813 | */ |
2814 | |
2815 | mlo = mhi; |
2816 | if (spec_case) { |
2817 | mhi = Balloc(mhi->k); |
2818 | Bcopy(mhi, mlo); |
2819 | mhi = lshift(mhi, Log2P); |
2820 | } |
2821 | |
2822 | for(i = 1;;i++) { |
2823 | dig = quorem(b,S) + '0'; |
2824 | /* Do we yet have the shortest decimal string |
2825 | * that will round to d? |
2826 | */ |
2827 | j = cmp(b, mlo); |
2828 | delta = diff(S, mhi); |
2829 | j1 = delta->sign ? 1 : cmp(b, delta); |
2830 | Bfree(delta); |
2831 | #ifndef ROUND_BIASED |
2832 | if (j1 == 0 && !mode && !(getWord1(d) & 1)) { |
2833 | if (dig == '9') |
2834 | goto round_9_up; |
2835 | if (j > 0) |
2836 | dig++; |
2837 | *s++ = dig; |
2838 | goto ret; |
2839 | } |
2840 | #endif |
2841 | if (j < 0 || (j == 0 && !mode |
2842 | #ifndef ROUND_BIASED |
2843 | && !(getWord1(d) & 1) |
2844 | #endif |
2845 | )) { |
2846 | if (j1 > 0) { |
2847 | b = lshift(b, 1); |
2848 | j1 = cmp(b, S); |
2849 | if ((j1 > 0 || (j1 == 0 && dig & 1)) |
2850 | && dig++ == '9') |
2851 | goto round_9_up; |
2852 | } |
2853 | *s++ = dig; |
2854 | goto ret; |
2855 | } |
2856 | if (j1 > 0) { |
2857 | if (dig == '9') { /* possible if i == 1 */ |
2858 | round_9_up: |
2859 | *s++ = '9'; |
2860 | goto roundoff; |
2861 | } |
2862 | *s++ = dig + 1; |
2863 | goto ret; |
2864 | } |
2865 | *s++ = dig; |
2866 | if (i == ilim) |
2867 | break; |
2868 | b = multadd(b, 10, 0); |
2869 | if (mlo == mhi) |
2870 | mlo = mhi = multadd(mhi, 10, 0); |
2871 | else { |
2872 | mlo = multadd(mlo, 10, 0); |
2873 | mhi = multadd(mhi, 10, 0); |
2874 | } |
2875 | } |
2876 | } |
2877 | else |
2878 | for(i = 1;; i++) { |
2879 | *s++ = dig = quorem(b,S) + '0'; |
2880 | if (i >= ilim) |
2881 | break; |
2882 | b = multadd(b, 10, 0); |
2883 | } |
2884 | |
2885 | /* Round off last digit */ |
2886 | |
2887 | b = lshift(b, 1); |
2888 | j = cmp(b, S); |
2889 | if (j > 0 || (j == 0 && dig & 1)) { |
2890 | roundoff: |
2891 | while(*--s == '9') |
2892 | if (s == s0) { |
2893 | k++; |
2894 | *s++ = '1'; |
2895 | goto ret; |
2896 | } |
2897 | ++*s++; |
2898 | } |
2899 | else { |
2900 | while(*--s == '0') {} |
2901 | s++; |
2902 | } |
2903 | ret: |
2904 | Bfree(S); |
2905 | if (mhi) { |
2906 | if (mlo && mlo != mhi) |
2907 | Bfree(mlo); |
2908 | Bfree(mhi); |
2909 | } |
2910 | ret1: |
2911 | Bfree(b); |
2912 | if (s == s0) { /* don't return empty string */ |
2913 | *s++ = '0'; |
2914 | k = 0; |
2915 | } |
2916 | *s = 0; |
2917 | *decpt = k + 1; |
2918 | if (rve) |
2919 | *rve = s; |
2920 | return s0; |
2921 | } |
2922 | #else |
2923 | // NOT thread safe! |
2924 | |
2925 | #include <errno.h> |
2926 | |
2927 | Q_CORE_EXPORT char *qdtoa( double d, int mode, int ndigits, int *decpt, int *sign, char **rve, char **resultp) |
2928 | { |
2929 | if(rve) |
2930 | *rve = 0; |
2931 | |
2932 | char *res; |
2933 | if (mode == 0) |
2934 | ndigits = 80; |
2935 | |
2936 | if (mode == 3) |
2937 | res = fcvt(d, ndigits, decpt, sign); |
2938 | else |
2939 | res = ecvt(d, ndigits, decpt, sign); |
2940 | |
2941 | int n = qstrlen(res); |
2942 | if (mode == 0) { // remove trailing 0's |
2943 | const int stop = qMax(1, *decpt); |
2944 | int i; |
2945 | for (i = n-1; i >= stop; --i) { |
2946 | if (res[i] != '0') |
2947 | break; |
2948 | } |
2949 | n = i + 1; |
2950 | } |
2951 | *resultp = static_cast<char*>(malloc(n + 1)); |
2952 | Q_CHECK_PTR(resultp); |
2953 | qstrncpy(*resultp, res, n + 1); |
2954 | return *resultp; |
2955 | } |
2956 | |
2957 | Q_CORE_EXPORT double qstrtod(const char *s00, const char **se, bool *ok) |
2958 | { |
2959 | double ret = strtod((char*)s00, (char**)se); |
2960 | if (ok) { |
2961 | if((ret == 0.0l && errno == ERANGE) |
2962 | || ret == HUGE_VAL || ret == -HUGE_VAL) |
2963 | *ok = false; |
2964 | else |
2965 | *ok = true; // the result will be that we don't report underflow in this case |
2966 | } |
2967 | return ret; |
2968 | } |
2969 | |
2970 | #endif // QT_QLOCALE_USES_FCVT |
2971 | |
2972 | QT_END_NAMESPACE |
2973 | |