src/corelib/tools/qlocale_tools.cpp

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