1
/* -*- Mode: C; tab-width: 8; indent-tabs-mode: nil; c-basic-offset: 4 -*-
3
* ***** BEGIN LICENSE BLOCK *****
4
* Version: MPL 1.1/GPL 2.0/LGPL 2.1
6
* The contents of this file are subject to the Mozilla Public License Version
7
* 1.1 (the "License"); you may not use this file except in compliance with
8
* the License. You may obtain a copy of the License at
9
* http://www.mozilla.org/MPL/
11
* Software distributed under the License is distributed on an "AS IS" basis,
12
* WITHOUT WARRANTY OF ANY KIND, either express or implied. See the License
13
* for the specific language governing rights and limitations under the
16
* The Original Code is Mozilla Communicator client code, released
19
* The Initial Developer of the Original Code is
20
* Netscape Communications Corporation.
21
* Portions created by the Initial Developer are Copyright (C) 1998
22
* the Initial Developer. All Rights Reserved.
26
* Alternatively, the contents of this file may be used under the terms of
27
* either of the GNU General Public License Version 2 or later (the "GPL"),
28
* or the GNU Lesser General Public License Version 2.1 or later (the "LGPL"),
29
* in which case the provisions of the GPL or the LGPL are applicable instead
30
* of those above. If you wish to allow use of your version of this file only
31
* under the terms of either the GPL or the LGPL, and not to allow others to
32
* use your version of this file under the terms of the MPL, indicate your
33
* decision by deleting the provisions above and replace them with the notice
34
* and other provisions required by the GPL or the LGPL. If you do not delete
35
* the provisions above, a recipient may use your version of this file under
36
* the terms of any one of the MPL, the GPL or the LGPL.
38
* ***** END LICENSE BLOCK ***** */
41
* Portable double to alphanumeric string and back converters.
44
#include "jslibmath.h"
48
#include "jsutil.h" /* Added by JSIFY */
56
/****************************************************************
58
* The author of this software is David M. Gay.
60
* Copyright (c) 1991 by Lucent Technologies.
62
* Permission to use, copy, modify, and distribute this software for any
63
* purpose without fee is hereby granted, provided that this entire notice
64
* is included in all copies of any software which is or includes a copy
65
* or modification of this software and in all copies of the supporting
66
* documentation for such software.
68
* THIS SOFTWARE IS BEING PROVIDED "AS IS", WITHOUT ANY EXPRESS OR IMPLIED
69
* WARRANTY. IN PARTICULAR, NEITHER THE AUTHOR NOR LUCENT MAKES ANY
70
* REPRESENTATION OR WARRANTY OF ANY KIND CONCERNING THE MERCHANTABILITY
71
* OF THIS SOFTWARE OR ITS FITNESS FOR ANY PARTICULAR PURPOSE.
73
***************************************************************/
75
/* Please send bug reports to
77
Bell Laboratories, Room 2C-463
79
Murray Hill, NJ 07974-0636
84
/* On a machine with IEEE extended-precision registers, it is
85
* necessary to specify double-precision (53-bit) rounding precision
86
* before invoking strtod or dtoa. If the machine uses (the equivalent
87
* of) Intel 80x87 arithmetic, the call
88
* _control87(PC_53, MCW_PC);
89
* does this with many compilers. Whether this or another call is
90
* appropriate depends on the compiler; for this to work, it may be
91
* necessary to #include "float.h" or another system-dependent header
95
/* strtod for IEEE-arithmetic machines.
97
* This strtod returns a nearest machine number to the input decimal
98
* string (or sets err to JS_DTOA_ERANGE or JS_DTOA_ENOMEM). With IEEE
99
* arithmetic, ties are broken by the IEEE round-even rule. Otherwise
100
* ties are broken by biased rounding (add half and chop).
102
* Inspired loosely by William D. Clinger's paper "How to Read Floating
103
* Point Numbers Accurately" [Proc. ACM SIGPLAN '90, pp. 92-101].
107
* 1. We only require IEEE double-precision
108
* arithmetic (not IEEE double-extended).
109
* 2. We get by with floating-point arithmetic in a case that
110
* Clinger missed -- when we're computing d * 10^n
111
* for a small integer d and the integer n is not too
112
* much larger than 22 (the maximum integer k for which
113
* we can represent 10^k exactly), we may be able to
114
* compute (d*10^k) * 10^(e-k) with just one roundoff.
115
* 3. Rather than a bit-at-a-time adjustment of the binary
116
* result in the hard case, we use floating-point
117
* arithmetic to determine the adjustment to within
118
* one bit; only in really hard cases do we need to
119
* compute a second residual.
120
* 4. Because of 3., we don't need a large table of powers of 10
121
* for ten-to-e (just some small tables, e.g. of 10^k
126
* #define IEEE_8087 for IEEE-arithmetic machines where the least
127
* significant byte has the lowest address.
128
* #define IEEE_MC68k for IEEE-arithmetic machines where the most
129
* significant byte has the lowest address.
130
* #define Long int on machines with 32-bit ints and 64-bit longs.
131
* #define Sudden_Underflow for IEEE-format machines without gradual
132
* underflow (i.e., that flush to zero on underflow).
133
* #define No_leftright to omit left-right logic in fast floating-point
134
* computation of js_dtoa.
135
* #define Check_FLT_ROUNDS if FLT_ROUNDS can assume the values 2 or 3.
136
* #define RND_PRODQUOT to use rnd_prod and rnd_quot (assembly routines
137
* that use extended-precision instructions to compute rounded
138
* products and quotients) with IBM.
139
* #define ROUND_BIASED for IEEE-format with biased rounding.
140
* #define Inaccurate_Divide for IEEE-format with correctly rounded
141
* products but inaccurate quotients, e.g., for Intel i860.
142
* #define JS_HAVE_LONG_LONG on machines that have a "long long"
143
* integer type (of >= 64 bits). If long long is available and the name is
144
* something other than "long long", #define Llong to be the name,
145
* and if "unsigned Llong" does not work as an unsigned version of
146
* Llong, #define #ULLong to be the corresponding unsigned type.
147
* #define Bad_float_h if your system lacks a float.h or if it does not
148
* define some or all of DBL_DIG, DBL_MAX_10_EXP, DBL_MAX_EXP,
149
* FLT_RADIX, FLT_ROUNDS, and DBL_MAX.
150
* #define MALLOC your_malloc, where your_malloc(n) acts like malloc(n)
151
* if memory is available and otherwise does something you deem
152
* appropriate. If MALLOC is undefined, malloc will be invoked
153
* directly -- and assumed always to succeed.
154
* #define Omit_Private_Memory to omit logic (added Jan. 1998) for making
155
* memory allocations from a private pool of memory when possible.
156
* When used, the private pool is PRIVATE_MEM bytes long: 2000 bytes,
157
* unless #defined to be a different length. This default length
158
* suffices to get rid of MALLOC calls except for unusual cases,
159
* such as decimal-to-binary conversion of a very long string of
161
* #define INFNAN_CHECK on IEEE systems to cause strtod to check for
162
* Infinity and NaN (case insensitively). On some systems (e.g.,
163
* some HP systems), it may be necessary to #define NAN_WORD0
164
* appropriately -- to the most significant word of a quiet NaN.
165
* (On HP Series 700/800 machines, -DNAN_WORD0=0x7ff40000 works.)
166
* #define MULTIPLE_THREADS if the system offers preemptively scheduled
167
* multiple threads. In this case, you must provide (or suitably
168
* #define) two locks, acquired by ACQUIRE_DTOA_LOCK() and released
169
* by RELEASE_DTOA_LOCK(). (The second lock, accessed
170
* in pow5mult, ensures lazy evaluation of only one copy of high
171
* powers of 5; omitting this lock would introduce a small
172
* probability of wasting memory, but would otherwise be harmless.)
173
* You must also invoke freedtoa(s) to free the value s returned by
174
* dtoa. You may do so whether or not MULTIPLE_THREADS is #defined.
175
* #define NO_IEEE_Scale to disable new (Feb. 1997) logic in strtod that
176
* avoids underflows on inputs whose result does not underflow.
178
#ifdef IS_LITTLE_ENDIAN
192
#define Bug(errorMessageString) JS_ASSERT(!errorMessageString)
198
extern void *MALLOC(size_t);
200
#define MALLOC malloc
203
#define Omit_Private_Memory
204
/* Private memory currently doesn't work with JS_THREADSAFE */
205
#ifndef Omit_Private_Memory
207
#define PRIVATE_MEM 2000
209
#define PRIVATE_mem ((PRIVATE_MEM+sizeof(double)-1)/sizeof(double))
210
static double private_mem[PRIVATE_mem], *pmem_next = private_mem;
217
#define DBL_MAX_10_EXP 308
218
#define DBL_MAX_EXP 1024
221
#define DBL_MAX 1.7976931348623157e+308
226
#define LONG_MAX 2147483647
229
#else /* ifndef Bad_float_h */
232
* MacOS 10.2 defines the macro FLT_ROUNDS to an internal function
233
* which does not exist on 10.1. We can safely #define it to 1 here
234
* to allow 10.2 builds to run on 10.1, since we can't use fesetround()
235
* (which does not exist on 10.1 either).
237
#if defined(MACOS_DEPLOYMENT_TARGET) && (MACOS_DEPLOYMENT_TARGET < 100200)
241
#endif /* Bad_float_h */
251
#if defined(IEEE_8087) + defined(IEEE_MC68k) != 1
252
Exactly one of IEEE_8087 or IEEE_MC68k should be defined.
255
#define word0(x) JSDOUBLE_HI32(x)
256
#define set_word0(x, y) JSDOUBLE_SET_HI32(x, y)
257
#define word1(x) JSDOUBLE_LO32(x)
258
#define set_word1(x, y) JSDOUBLE_SET_LO32(x, y)
260
#define Storeinc(a,b,c) (*(a)++ = (b) << 16 | (c) & 0xffff)
262
/* #define P DBL_MANT_DIG */
263
/* Ten_pmax = floor(P*log(2)/log(5)) */
264
/* Bletch = (highest power of 2 < DBL_MAX_10_EXP) / 16 */
265
/* Quick_max = floor((P-1)*log(FLT_RADIX)/log(10) - 1) */
266
/* Int_max = floor(P*log(FLT_RADIX)/log(10) - 1) */
269
#define Exp_shift1 20
270
#define Exp_msk1 0x100000
271
#define Exp_msk11 0x100000
272
#define Exp_mask 0x7ff00000
276
#define Exp_1 0x3ff00000
277
#define Exp_11 0x3ff00000
279
#define Frac_mask 0xfffff
280
#define Frac_mask1 0xfffff
283
#define Bndry_mask 0xfffff
284
#define Bndry_mask1 0xfffff
286
#define Sign_bit 0x80000000
292
#define Infinite(x) (word0(x) == 0x7ff00000) /* sufficient test for here */
293
#ifndef NO_IEEE_Scale
294
#define Avoid_Underflow
300
#define rounded_product(a,b) a = rnd_prod(a, b)
301
#define rounded_quotient(a,b) a = rnd_quot(a, b)
302
extern double rnd_prod(double, double), rnd_quot(double, double);
304
#define rounded_product(a,b) a *= b
305
#define rounded_quotient(a,b) a /= b
308
#define Big0 (Frac_mask1 | Exp_msk1*(DBL_MAX_EXP+Bias-1))
309
#define Big1 0xffffffff
311
#ifndef JS_HAVE_LONG_LONG
313
#else /* long long available */
315
#define Llong JSInt64
318
#define ULLong JSUint64
320
#endif /* JS_HAVE_LONG_LONG */
323
#define MULTIPLE_THREADS
324
static PRLock *freelist_lock;
325
#define ACQUIRE_DTOA_LOCK() \
329
PR_Lock(freelist_lock); \
331
#define RELEASE_DTOA_LOCK() PR_Unlock(freelist_lock)
333
#undef MULTIPLE_THREADS
334
#define ACQUIRE_DTOA_LOCK() /*nothing*/
335
#define RELEASE_DTOA_LOCK() /*nothing*/
341
struct Bigint *next; /* Free list link */
342
int32 k; /* lg2(maxwds) */
343
int32 maxwds; /* Number of words allocated for x */
344
int32 sign; /* Zero if positive, 1 if negative. Ignored by most Bigint routines! */
345
int32 wds; /* Actual number of words. If value is nonzero, the most significant word must be nonzero. */
346
ULong x[1]; /* wds words of number in little endian order */
349
#ifdef ENABLE_OOM_TESTING
350
/* Out-of-memory testing. Use a good testcase (over and over) and then use
351
* these routines to cause a memory failure on every possible Balloc allocation,
352
* to make sure that all out-of-memory paths can be followed. See bug 14044.
355
static int allocationNum; /* which allocation is next? */
356
static int desiredFailure; /* which allocation should fail? */
359
* js_BigintTestingReset
361
* Call at the beginning of a test run to set the allocation failure position.
362
* (Set to 0 to just have the engine count allocations without failing.)
365
js_BigintTestingReset(int newFailure)
368
desiredFailure = newFailure;
372
* js_BigintTestingWhere
374
* Report the current allocation position. This is really only useful when you
375
* want to learn how many allocations a test run has.
378
js_BigintTestingWhere()
380
return allocationNum;
385
* So here's what you do: Set up a fantastic test case that exercises the
386
* elements of the code you wish. Set the failure point at 0 and run the test,
387
* then get the allocation position. This number is the number of allocations
388
* your test makes. Now loop from 1 to that number, setting the failure point
389
* at each loop count, and run the test over and over, causing failures at each
390
* step. Any memory failure *should* cause a Out-Of-Memory exception; if it
391
* doesn't, then there's still an error here.
395
typedef struct Bigint Bigint;
397
static Bigint *freelist[Kmax+1];
400
* Allocate a Bigint with 2^k words.
401
* This is not threadsafe. The caller must use thread locks
403
static Bigint *Balloc(int32 k)
407
#ifndef Omit_Private_Memory
411
#ifdef ENABLE_OOM_TESTING
412
if (++allocationNum == desiredFailure) {
413
printf("Forced Failing Allocation number %d\n", allocationNum);
418
if ((rv = freelist[k]) != NULL)
419
freelist[k] = rv->next;
422
#ifdef Omit_Private_Memory
423
rv = (Bigint *)MALLOC(sizeof(Bigint) + (x-1)*sizeof(ULong));
425
len = (sizeof(Bigint) + (x-1)*sizeof(ULong) + sizeof(double) - 1)
427
if (pmem_next - private_mem + len <= PRIVATE_mem) {
428
rv = (Bigint*)pmem_next;
432
rv = (Bigint*)MALLOC(len*sizeof(double));
439
rv->sign = rv->wds = 0;
443
static void Bfree(Bigint *v)
446
v->next = freelist[v->k];
451
#define Bcopy(x,y) memcpy((char *)&x->sign, (char *)&y->sign, \
452
y->wds*sizeof(Long) + 2*sizeof(int32))
454
/* Return b*m + a. Deallocate the old b. Both a and m must be between 0 and
455
* 65535 inclusive. NOTE: old b is deallocated on memory failure.
457
static Bigint *multadd(Bigint *b, int32 m, int32 a)
469
#ifdef ENABLE_OOM_TESTING
470
if (++allocationNum == desiredFailure) {
471
/* Faux allocation, because I'm not getting all of the failure paths
474
printf("Forced Failing Allocation number %d\n", allocationNum);
486
y = *x * (ULLong)m + carry;
488
*x++ = (ULong)(y & 0xffffffffUL);
491
y = (xi & 0xffff) * m + carry;
492
z = (xi >> 16) * m + (y >> 16);
494
*x++ = (z << 16) + (y & 0xffff);
499
if (wds >= b->maxwds) {
509
b->x[wds++] = (ULong)carry;
515
static Bigint *s2b(CONST char *s, int32 nd0, int32 nd, ULong y9)
522
for(k = 0, y = 1; x > y; y <<= 1, k++) ;
533
b = multadd(b, 10, *s++ - '0');
542
b = multadd(b, 10, *s++ - '0');
550
/* Return the number (0 through 32) of most significant zero bits in x. */
551
static int32 hi0bits(register ULong x)
553
register int32 k = 0;
555
if (!(x & 0xffff0000)) {
559
if (!(x & 0xff000000)) {
563
if (!(x & 0xf0000000)) {
567
if (!(x & 0xc0000000)) {
571
if (!(x & 0x80000000)) {
573
if (!(x & 0x40000000))
580
/* Return the number (0 through 32) of least significant zero bits in y.
581
* Also shift y to the right past these 0 through 32 zeros so that y's
582
* least significant bit will be set unless y was originally zero. */
583
static int32 lo0bits(ULong *y)
586
register ULong x = *y;
625
/* Return a new Bigint with the given integer value, which must be nonnegative. */
626
static Bigint *i2b(int32 i)
638
/* Return a newly allocated product of a and b. */
639
static Bigint *mult(CONST Bigint *a, CONST Bigint *b)
645
ULong *xc, *xc0, *xce;
646
CONST ULong *x, *xa, *xae, *xb, *xbe;
654
if (a->wds < b->wds) {
668
for(xc = c->x, xce = xc + wc; xc < xce; xc++)
676
for(; xb < xbe; xc0++) {
677
if ((y = *xb++) != 0) {
682
z = *x++ * (ULLong)y + *xc + carry;
684
*xc++ = (ULong)(z & 0xffffffffUL);
691
for(; xb < xbe; xb++, xc0++) {
692
if ((y = *xb & 0xffff) != 0) {
697
z = (*x & 0xffff) * y + (*xc & 0xffff) + carry;
699
z2 = (*x++ >> 16) * y + (*xc >> 16) + carry;
706
if ((y = *xb >> 16) != 0) {
712
z = (*x & 0xffff) * y + (*xc >> 16) + carry;
715
z2 = (*x++ >> 16) * y + (*xc & 0xffff) + carry;
723
for(xc0 = c->x, xc = xc0 + wc; wc > 0 && !*--xc; --wc) ;
729
* 'p5s' points to a linked list of Bigints that are powers of 5.
730
* This list grows on demand, and it can only grow: it won't change
731
* in any other way. So if we read 'p5s' or the 'next' field of
732
* some Bigint on the list, and it is not NULL, we know it won't
733
* change to NULL or some other value. Only when the value of
734
* 'p5s' or 'next' is NULL do we need to acquire the lock and add
735
* a new Bigint to the list.
741
static PRLock *p5s_lock;
744
/* Return b * 5^k. Deallocate the old b. k must be nonnegative. */
745
/* NOTE: old b is deallocated on memory failure. */
746
static Bigint *pow5mult(Bigint *b, int32 k)
748
Bigint *b1, *p5, *p51;
750
static CONST int32 p05[3] = { 5, 25, 125 };
752
if ((i = k & 3) != 0) {
753
b = multadd(b, p05[i-1], 0);
763
* We take great care to not call i2b() and Bfree()
764
* while holding the lock.
766
Bigint *wasted_effort = NULL;
772
/* lock and check again */
779
/* some other thread just beat us */
785
Bfree(wasted_effort);
807
if (!(p51 = p5->next)) {
809
Bigint *wasted_effort = NULL;
825
Bfree(wasted_effort);
842
/* Return b * 2^k. Deallocate the old b. k must be nonnegative.
843
* NOTE: on memory failure, old b is deallocated. */
844
static Bigint *lshift(Bigint *b, int32 k)
848
ULong *x, *x1, *xe, z;
853
for(i = b->maxwds; n1 > i; i <<= 1)
859
for(i = 0; i < n; i++)
883
/* Return -1, 0, or 1 depending on whether a<b, a==b, or a>b, respectively. */
884
static int32 cmp(Bigint *a, Bigint *b)
886
ULong *xa, *xa0, *xb, *xb0;
892
if (i > 1 && !a->x[i-1])
893
Bug("cmp called with a->x[a->wds-1] == 0");
894
if (j > 1 && !b->x[j-1])
895
Bug("cmp called with b->x[b->wds-1] == 0");
905
return *xa < *xb ? -1 : 1;
912
static Bigint *diff(Bigint *a, Bigint *b)
916
ULong *xa, *xae, *xb, *xbe, *xc;
955
y = (ULLong)*xa++ - *xb++ - borrow;
956
borrow = y >> 32 & 1UL;
957
*xc++ = (ULong)(y & 0xffffffffUL);
962
borrow = y >> 32 & 1UL;
963
*xc++ = (ULong)(y & 0xffffffffUL);
967
y = (*xa & 0xffff) - (*xb & 0xffff) - borrow;
968
borrow = (y & 0x10000) >> 16;
969
z = (*xa++ >> 16) - (*xb++ >> 16) - borrow;
970
borrow = (z & 0x10000) >> 16;
975
y = (*xa & 0xffff) - borrow;
976
borrow = (y & 0x10000) >> 16;
977
z = (*xa++ >> 16) - borrow;
978
borrow = (z & 0x10000) >> 16;
988
/* Return the absolute difference between x and the adjacent greater-magnitude double number (ignoring exponent overflows). */
989
static double ulp(double x)
994
L = (word0(x) & Exp_mask) - (P-1)*Exp_msk1;
995
#ifndef Sudden_Underflow
1000
#ifndef Sudden_Underflow
1003
L = -L >> Exp_shift;
1004
if (L < Exp_shift) {
1005
set_word0(a, 0x80000 >> L);
1011
set_word1(a, L >= 31 ? 1 : 1 << (31 - L));
1019
static double b2d(Bigint *a, int32 *e)
1021
ULong *xa, *xa0, w, y, z;
1026
#define set_d0(x) set_word0(d, x)
1027
#define set_d1(x) set_word1(d, x)
1033
if (!y) Bug("zero y in b2d");
1038
set_d0(Exp_1 | y >> (Ebits - k));
1039
w = xa > xa0 ? *--xa : 0;
1040
set_d1(y << (32-Ebits + k) | w >> (Ebits - k));
1043
z = xa > xa0 ? *--xa : 0;
1045
set_d0(Exp_1 | y << k | z >> (32 - k));
1046
y = xa > xa0 ? *--xa : 0;
1047
set_d1(z << k | y >> (32 - k));
1062
/* Convert d into the form b*2^e, where b is an odd integer. b is the returned
1063
* Bigint and e is the returned binary exponent. Return the number of significant
1064
* bits in b in bits. d must be finite and nonzero. */
1065
static Bigint *d2b(double d, int32 *e, int32 *bits)
1072
#define set_d0(x) set_word0(d, x)
1073
#define set_d1(x) set_word1(d, x)
1081
set_d0(d0 & 0x7fffffff); /* clear sign bit, which we ignore */
1082
#ifdef Sudden_Underflow
1083
de = (int32)(d0 >> Exp_shift);
1086
if ((de = (int32)(d0 >> Exp_shift)) != 0)
1089
if ((y = d1) != 0) {
1090
if ((k = lo0bits(&y)) != 0) {
1091
x[0] = y | z << (32 - k);
1096
i = b->wds = (x[1] = z) ? 2 : 1;
1105
#ifndef Sudden_Underflow
1108
*e = de - Bias - (P-1) + k;
1110
#ifndef Sudden_Underflow
1113
*e = de - Bias - (P-1) + 1 + k;
1114
*bits = 32*i - hi0bits(x[i-1]);
1125
static double ratio(Bigint *a, Bigint *b)
1132
k = ka - kb + 32*(a->wds - b->wds);
1134
set_word0(da, word0(da) + k*Exp_msk1);
1137
set_word0(db, word0(db) + k*Exp_msk1);
1144
1e0, 1e1, 1e2, 1e3, 1e4, 1e5, 1e6, 1e7, 1e8, 1e9,
1145
1e10, 1e11, 1e12, 1e13, 1e14, 1e15, 1e16, 1e17, 1e18, 1e19,
1149
static CONST double bigtens[] = { 1e16, 1e32, 1e64, 1e128, 1e256 };
1150
static CONST double tinytens[] = { 1e-16, 1e-32, 1e-64, 1e-128,
1151
#ifdef Avoid_Underflow
1152
9007199254740992.e-256
1157
/* The factor of 2^53 in tinytens[4] helps us avoid setting the underflow */
1158
/* flag unnecessarily. It leads to a song and dance at the end of strtod. */
1159
#define Scale_Bit 0x10
1166
#define NAN_WORD0 0x7ff80000
1173
static int match(CONST char **sp, char *t)
1176
CONST char *s = *sp;
1179
if ((c = *++s) >= 'A' && c <= 'Z')
1187
#endif /* INFNAN_CHECK */
1190
#ifdef JS_THREADSAFE
1191
static JSBool initialized = JS_FALSE;
1193
/* hacked replica of nspr _PR_InitDtoa */
1194
static void InitDtoa(void)
1196
freelist_lock = PR_NewLock();
1197
p5s_lock = PR_NewLock();
1198
initialized = JS_TRUE;
1202
void js_FinishDtoa(void)
1207
#ifdef JS_THREADSAFE
1208
if (initialized == JS_TRUE) {
1209
PR_DestroyLock(freelist_lock);
1210
PR_DestroyLock(p5s_lock);
1211
initialized = JS_FALSE;
1215
/* clear down the freelist array and p5s */
1217
/* static Bigint *freelist[Kmax+1]; */
1218
for (count = 0; count <= Kmax; count++) {
1219
Bigint **listp = &freelist[count];
1220
while ((temp = *listp) != NULL) {
1221
*listp = temp->next;
1224
freelist[count] = NULL;
1227
/* static Bigint *p5s; */
1235
/* nspr2 watcom bug ifdef omitted */
1237
JS_FRIEND_API(double)
1238
JS_strtod(CONST char *s00, char **se, int *err)
1241
int32 bb2, bb5, bbe, bd2, bd5, bbbits, bs2, c, dsign,
1242
e, e1, esign, i, j, k, nd, nd0, nf, nz, nz0, sign;
1243
CONST char *s, *s0, *s1;
1244
double aadj, aadj1, adj, rv, rv0;
1247
Bigint *bb, *bb1, *bd, *bd0, *bs, *delta;
1253
bb = bd = bs = delta = NULL;
1254
sign = nz0 = nz = 0;
1257
/* Locking for Balloc's shared buffers that will be used in this block */
1258
ACQUIRE_DTOA_LOCK();
1260
for(s = s00;;s++) switch(*s) {
1285
while(*++s == '0') ;
1291
for(nd = nf = 0; (c = *s) >= '0' && c <= '9'; nd++, s++)
1300
for(; c == '0'; c = *++s)
1302
if (c > '0' && c <= '9') {
1310
for(; c >= '0' && c <= '9'; c = *++s) {
1315
for(i = 1; i < nz; i++)
1318
else if (nd <= DBL_DIG + 1)
1322
else if (nd <= DBL_DIG + 1)
1330
if (c == 'e' || c == 'E') {
1331
if (!nd && !nz && !nz0) {
1343
if (c >= '0' && c <= '9') {
1346
if (c > '0' && c <= '9') {
1349
while((c = *++s) >= '0' && c <= '9')
1351
if (s - s1 > 8 || L > 19999)
1352
/* Avoid confusion from exponents
1353
* so large that e might overflow.
1355
e = 19999; /* safe for 16 bit ints */
1370
/* Check for Nan and Infinity */
1374
if (match(&s,"nfinity")) {
1375
word0(rv) = 0x7ff00000;
1382
if (match(&s, "an")) {
1383
word0(rv) = NAN_WORD0;
1384
word1(rv) = NAN_WORD1;
1388
#endif /* INFNAN_CHECK */
1395
/* Now we have nd0 digits, starting at s0, followed by a
1396
* decimal point, followed by nd-nd0 digits. The number we're
1397
* after is the integer represented by those digits times
1402
k = nd < DBL_DIG + 1 ? nd : DBL_DIG + 1;
1405
rv = tens[k - 9] * rv + z;
1408
#ifndef RND_PRODQUOT
1415
if (e <= Ten_pmax) {
1416
/* rv = */ rounded_product(rv, tens[e]);
1420
if (e <= Ten_pmax + i) {
1421
/* A fancier test would sometimes let us do
1422
* this for larger i values.
1426
/* rv = */ rounded_product(rv, tens[e]);
1430
#ifndef Inaccurate_Divide
1431
else if (e >= -Ten_pmax) {
1432
/* rv = */ rounded_quotient(rv, tens[-e]);
1441
/* Get starting approximation = rv * 10**e1 */
1444
if ((i = e1 & 15) != 0)
1447
if (e1 > DBL_MAX_10_EXP) {
1449
*err = JS_DTOA_ERANGE;
1453
/* Can't trust HUGE_VAL */
1454
word0(rv) = Exp_mask;
1462
for(j = 0; e1 > 1; j++, e1 >>= 1)
1465
/* The last multiplication could overflow. */
1466
set_word0(rv, word0(rv) - P*Exp_msk1);
1468
if ((z = word0(rv) & Exp_mask) > Exp_msk1*(DBL_MAX_EXP+Bias-P))
1470
if (z > Exp_msk1*(DBL_MAX_EXP+Bias-1-P)) {
1471
/* set to largest number */
1472
/* (Can't trust DBL_MAX) */
1473
set_word0(rv, Big0);
1474
set_word1(rv, Big1);
1477
set_word0(rv, word0(rv) + P*Exp_msk1);
1482
if ((i = e1 & 15) != 0)
1486
if (e1 >= 1 << n_bigtens)
1488
#ifdef Avoid_Underflow
1491
for(j = 0; e1 > 0; j++, e1 >>= 1)
1494
if (scale && (j = P + 1 - ((word0(rv) & Exp_mask)
1495
>> Exp_shift)) > 0) {
1496
/* scaled rv is denormal; zap j low bits */
1499
set_word0(rv, word0(rv) & (0xffffffff << (j-32)));
1504
set_word1(rv, word1(rv) & (0xffffffff << j));
1507
for(j = 0; e1 > 1; j++, e1 >>= 1)
1510
/* The last multiplication could underflow. */
1520
*err = JS_DTOA_ERANGE;
1525
#ifndef Avoid_Underflow
1526
set_word0(rv, Tiny0);
1527
set_word1(rv, Tiny1);
1528
/* The refinement below will clean
1529
* this approximation up.
1536
/* Now the hard part -- adjusting rv to the correct value.*/
1538
/* Put digits into bd: true value = bd * 10^e */
1540
bd0 = s2b(s0, nd0, nd, y);
1545
bd = Balloc(bd0->k);
1549
bb = d2b(rv, &bbe, &bbbits); /* rv = bb * 2^bbe */
1569
#ifdef Sudden_Underflow
1572
#ifdef Avoid_Underflow
1577
i = j + bbbits - 1; /* logb(rv) */
1578
if (i < Emin) /* denormal */
1585
#ifdef Avoid_Underflow
1588
i = bb2 < bd2 ? bb2 : bd2;
1597
bs = pow5mult(bs, bb5);
1607
bb = lshift(bb, bb2);
1612
bd = pow5mult(bd, bd5);
1617
bd = lshift(bd, bd2);
1622
bs = lshift(bs, bs2);
1626
delta = diff(bb, bd);
1629
dsign = delta->sign;
1633
/* Error is less than half an ulp -- check for
1634
* special case of mantissa a power of two.
1636
if (dsign || word1(rv) || word0(rv) & Bndry_mask
1637
#ifdef Avoid_Underflow
1638
|| (word0(rv) & Exp_mask) <= Exp_msk1 + P*Exp_msk1
1640
|| (word0(rv) & Exp_mask) <= Exp_msk1
1643
#ifdef Avoid_Underflow
1644
if (!delta->x[0] && delta->wds == 1)
1649
delta = lshift(delta,Log2P);
1652
if (cmp(delta, bs) > 0)
1657
/* exactly half-way between */
1659
if ((word0(rv) & Bndry_mask1) == Bndry_mask1
1660
&& word1(rv) == 0xffffffff) {
1661
/*boundary case -- increment exponent*/
1662
set_word0(rv, (word0(rv) & Exp_mask) + Exp_msk1);
1664
#ifdef Avoid_Underflow
1670
else if (!(word0(rv) & Bndry_mask) && !word1(rv)) {
1671
#ifdef Avoid_Underflow
1675
/* boundary case -- decrement exponent */
1676
#ifdef Sudden_Underflow
1677
L = word0(rv) & Exp_mask;
1682
L = (word0(rv) & Exp_mask) - Exp_msk1;
1684
set_word0(rv, L | Bndry_mask1);
1685
set_word1(rv, 0xffffffff);
1688
#ifndef ROUND_BIASED
1689
if (!(word1(rv) & LSB))
1694
#ifndef ROUND_BIASED
1697
#ifndef Sudden_Underflow
1702
#ifdef Avoid_Underflow
1708
if ((aadj = ratio(delta, bs)) <= 2.) {
1711
else if (word1(rv) || word0(rv) & Bndry_mask) {
1712
#ifndef Sudden_Underflow
1713
if (word1(rv) == Tiny1 && !word0(rv))
1720
/* special case -- power of FLT_RADIX to be */
1721
/* rounded down... */
1723
if (aadj < 2./FLT_RADIX)
1724
aadj = 1./FLT_RADIX;
1732
aadj1 = dsign ? aadj : -aadj;
1733
#ifdef Check_FLT_ROUNDS
1734
switch(FLT_ROUNDS) {
1735
case 2: /* towards +infinity */
1738
case 0: /* towards 0 */
1739
case 3: /* towards -infinity */
1743
if (FLT_ROUNDS == 0)
1747
y = word0(rv) & Exp_mask;
1749
/* Check for overflow */
1751
if (y == Exp_msk1*(DBL_MAX_EXP+Bias-1)) {
1753
set_word0(rv, word0(rv) - P*Exp_msk1);
1754
adj = aadj1 * ulp(rv);
1756
if ((word0(rv) & Exp_mask) >=
1757
Exp_msk1*(DBL_MAX_EXP+Bias-P)) {
1758
if (word0(rv0) == Big0 && word1(rv0) == Big1)
1760
set_word0(rv, Big0);
1761
set_word1(rv, Big1);
1765
set_word0(rv, word0(rv) + P*Exp_msk1);
1768
#ifdef Sudden_Underflow
1769
if ((word0(rv) & Exp_mask) <= P*Exp_msk1) {
1771
set_word0(rv, word0(rv) + P*Exp_msk1);
1772
adj = aadj1 * ulp(rv);
1774
if ((word0(rv) & Exp_mask) <= P*Exp_msk1)
1776
if (word0(rv0) == Tiny0
1777
&& word1(rv0) == Tiny1)
1779
set_word0(rv, Tiny0);
1780
set_word1(rv, Tiny1);
1784
set_word0(rv, word0(rv) - P*Exp_msk1);
1787
adj = aadj1 * ulp(rv);
1791
/* Compute adj so that the IEEE rounding rules will
1792
* correctly round rv + adj in some half-way cases.
1793
* If rv * ulp(rv) is denormalized (i.e.,
1794
* y <= (P-1)*Exp_msk1), we must adjust aadj to avoid
1795
* trouble from bits lost to denormalization;
1796
* example: 1.2e-307 .
1798
#ifdef Avoid_Underflow
1799
if (y <= P*Exp_msk1 && aadj > 1.)
1801
if (y <= (P-1)*Exp_msk1 && aadj > 1.)
1804
aadj1 = (double)(int32)(aadj + 0.5);
1808
#ifdef Avoid_Underflow
1809
if (scale && y <= P*Exp_msk1)
1810
set_word0(aadj1, word0(aadj1) + (P+1)*Exp_msk1 - y);
1812
adj = aadj1 * ulp(rv);
1816
z = word0(rv) & Exp_mask;
1817
#ifdef Avoid_Underflow
1821
/* Can we stop now? */
1824
/* The tolerances below are conservative. */
1825
if (dsign || word1(rv) || word0(rv) & Bndry_mask) {
1826
if (aadj < .4999999 || aadj > .5000001)
1829
else if (aadj < .4999999/FLT_RADIX)
1837
bb = bd = bs = delta = NULL;
1839
#ifdef Avoid_Underflow
1841
set_word0(rv0, Exp_1 - P*Exp_msk1);
1843
if ((word0(rv) & Exp_mask) <= P*Exp_msk1
1847
#ifdef Sudden_Underflow
1848
/* rv will be 0, but this would give the */
1849
/* right result if only rv *= rv0 worked. */
1850
set_word0(rv, word0(rv) + P*Exp_msk1);
1851
set_word0(rv0, Exp_1 - 2*P*Exp_msk1);
1856
set_word1(rv, word1(rv) & ~1);
1860
#endif /* Avoid_Underflow */
1868
RELEASE_DTOA_LOCK();
1871
rv0 = sign ? -rv : rv;
1880
*err = JS_DTOA_ENOMEM;
1889
/* Return floor(b/2^k) and set b to be the remainder. The returned quotient must be less than 2^32. */
1890
static uint32 quorem2(Bigint *b, int32 k)
1909
JS_ASSERT(!(bxe[1] & ~mask));
1911
result |= bxe[1] << (32 - k);
1914
while (!*bxe && bxe != bx) {
1922
/* Return floor(b/S) and set b to be the remainder. As added restrictions, b must not have
1923
* more words than S, the most significant word of S must not start with a 1 bit, and the
1924
* returned quotient must be less than 36. */
1925
static int32 quorem(Bigint *b, Bigint *S)
1928
ULong *bx, *bxe, q, *sx, *sxe;
1930
ULLong borrow, carry, y, ys;
1932
ULong borrow, carry, y, ys;
1937
JS_ASSERT(b->wds <= n);
1944
JS_ASSERT(*sxe <= 0x7FFFFFFF);
1945
q = *bxe / (*sxe + 1); /* ensure q <= true quotient */
1952
ys = *sx++ * (ULLong)q + carry;
1954
y = *bx - (ys & 0xffffffffUL) - borrow;
1955
borrow = y >> 32 & 1UL;
1956
*bx++ = (ULong)(y & 0xffffffffUL);
1959
ys = (si & 0xffff) * q + carry;
1960
zs = (si >> 16) * q + (ys >> 16);
1962
y = (*bx & 0xffff) - (ys & 0xffff) - borrow;
1963
borrow = (y & 0x10000) >> 16;
1964
z = (*bx >> 16) - (zs & 0xffff) - borrow;
1965
borrow = (z & 0x10000) >> 16;
1972
while(--bxe > bx && !*bxe)
1977
if (cmp(b, S) >= 0) {
1987
y = *bx - (ys & 0xffffffffUL) - borrow;
1988
borrow = y >> 32 & 1UL;
1989
*bx++ = (ULong)(y & 0xffffffffUL);
1992
ys = (si & 0xffff) + carry;
1993
zs = (si >> 16) + (ys >> 16);
1995
y = (*bx & 0xffff) - (ys & 0xffff) - borrow;
1996
borrow = (y & 0x10000) >> 16;
1997
z = (*bx >> 16) - (zs & 0xffff) - borrow;
1998
borrow = (z & 0x10000) >> 16;
2005
while(--bxe > bx && !*bxe)
2013
/* dtoa for IEEE arithmetic (dmg): convert double to ASCII string.
2015
* Inspired by "How to Print Floating-Point Numbers Accurately" by
2016
* Guy L. Steele, Jr. and Jon L. White [Proc. ACM SIGPLAN '90, pp. 92-101].
2019
* 1. Rather than iterating, we use a simple numeric overestimate
2020
* to determine k = floor(log10(d)). We scale relevant
2021
* quantities using O(log2(k)) rather than O(k) multiplications.
2022
* 2. For some modes > 2 (corresponding to ecvt and fcvt), we don't
2023
* try to generate digits strictly left to right. Instead, we
2024
* compute with fewer bits and propagate the carry if necessary
2025
* when rounding the final digit up. This is often faster.
2026
* 3. Under the assumption that input will be rounded nearest,
2027
* mode 0 renders 1e23 as 1e23 rather than 9.999999999999999e22.
2028
* That is, we allow equality in stopping tests when the
2029
* round-nearest rule will give the same floating-point value
2030
* as would satisfaction of the stopping test with strict
2032
* 4. We remove common factors of powers of 2 from relevant
2034
* 5. When converting floating-point integers less than 1e16,
2035
* we use floating-point arithmetic rather than resorting
2036
* to multiple-precision integers.
2037
* 6. When asked to produce fewer than 15 digits, we first try
2038
* to get by with floating-point arithmetic; we resort to
2039
* multiple-precision integer arithmetic only if we cannot
2040
* guarantee that the floating-point calculation has given
2041
* the correctly rounded result. For k requested digits and
2042
* "uniformly" distributed input, the probability is
2043
* something like 10^(k-15) that we must resort to the Long
2047
/* Always emits at least one digit. */
2048
/* If biasUp is set, then rounding in modes 2 and 3 will round away from zero
2049
* when the number is exactly halfway between two representable values. For example,
2050
* rounding 2.5 to zero digits after the decimal point will return 3 and not 2.
2051
* 2.49 will still round to 2, and 2.51 will still round to 3. */
2052
/* bufsize should be at least 20 for modes 0 and 1. For the other modes,
2053
* bufsize should be two greater than the maximum number of output characters expected. */
2055
js_dtoa(double d, int mode, JSBool biasUp, int ndigits,
2056
int *decpt, int *sign, char **rve, char *buf, size_t bufsize)
2058
/* Arguments ndigits, decpt, sign are similar to those
2059
of ecvt and fcvt; trailing zeros are suppressed from
2060
the returned string. If not null, *rve is set to point
2061
to the end of the return value. If d is +-Infinity or NaN,
2062
then *decpt is set to 9999.
2065
0 ==> shortest string that yields d when read in
2066
and rounded to nearest.
2067
1 ==> like 0, but with Steele & White stopping rule;
2068
e.g. with IEEE P754 arithmetic , mode 0 gives
2069
1e23 whereas mode 1 gives 9.999999999999999e22.
2070
2 ==> max(1,ndigits) significant digits. This gives a
2071
return value similar to that of ecvt, except
2072
that trailing zeros are suppressed.
2073
3 ==> through ndigits past the decimal point. This
2074
gives a return value similar to that from fcvt,
2075
except that trailing zeros are suppressed, and
2076
ndigits can be negative.
2077
4-9 should give the same return values as 2-3, i.e.,
2078
4 <= mode <= 9 ==> same return as mode
2079
2 + (mode & 1). These modes are mainly for
2080
debugging; often they run slower but sometimes
2081
faster than modes 2-3.
2082
4,5,8,9 ==> left-to-right digit generation.
2083
6-9 ==> don't try fast floating-point estimate
2086
Values of mode other than 0-9 are treated as mode 0.
2088
Sufficient space is allocated to the return value
2089
to hold the suppressed trailing zeros.
2092
int32 bbits, b2, b5, be, dig, i, ieps, ilim, ilim0, ilim1,
2093
j, j1, k, k0, k_check, leftright, m2, m5, s2, s5,
2094
spec_case, try_quick;
2096
#ifndef Sudden_Underflow
2100
Bigint *b, *b1, *delta, *mlo, *mhi, *S;
2107
if (word0(d) & Sign_bit) {
2108
/* set sign for everything, including 0's and NaNs */
2110
set_word0(d, word0(d) & ~Sign_bit); /* clear sign bit */
2115
if ((word0(d) & Exp_mask) == Exp_mask) {
2116
/* Infinity or NaN */
2118
s = !word1(d) && !(word0(d) & Frac_mask) ? "Infinity" : "NaN";
2119
if ((s[0] == 'I' && bufsize < 9) || (s[0] == 'N' && bufsize < 4)) {
2120
JS_ASSERT(JS_FALSE);
2121
/* JS_SetError(JS_BUFFER_OVERFLOW_ERROR, 0); */
2127
*rve = buf[3] ? buf + 8 : buf + 3;
2128
JS_ASSERT(**rve == '\0');
2134
b = NULL; /* initialize for abort protection */
2142
JS_ASSERT(JS_FALSE);
2143
/* JS_SetError(JS_BUFFER_OVERFLOW_ERROR, 0); */
2147
buf[0] = '0'; buf[1] = '\0'; /* copy "0" to buffer */
2150
/* We might have jumped to "no_digits" from below, so we need
2151
* to be sure to free the potentially allocated Bigints to avoid
2162
b = d2b(d, &be, &bbits);
2165
#ifdef Sudden_Underflow
2166
i = (int32)(word0(d) >> Exp_shift1 & (Exp_mask>>Exp_shift1));
2168
if ((i = (int32)(word0(d) >> Exp_shift1 & (Exp_mask>>Exp_shift1))) != 0) {
2171
set_word0(d2, word0(d2) & Frac_mask1);
2172
set_word0(d2, word0(d2) | Exp_11);
2174
/* log(x) ~=~ log(1.5) + (x-1.5)/1.5
2175
* log10(x) = log(x) / log(10)
2176
* ~=~ log(1.5)/log(10) + (x-1.5)/(1.5*log(10))
2177
* log10(d) = (i-Bias)*log(2)/log(10) + log10(d2)
2179
* This suggests computing an approximation k to log10(d) by
2181
* k = (i - Bias)*0.301029995663981
2182
* + ( (d2-1.5)*0.289529654602168 + 0.176091259055681 );
2184
* We want k to be too large rather than too small.
2185
* The error in the first-order Taylor series approximation
2186
* is in our favor, so we just round up the constant enough
2187
* to compensate for any error in the multiplication of
2188
* (i - Bias) by 0.301029995663981; since |i - Bias| <= 1077,
2189
* and 1077 * 0.30103 * 2^-52 ~=~ 7.2e-14,
2190
* adding 1e-13 to the constant term more than suffices.
2191
* Hence we adjust the constant term to 0.1760912590558.
2192
* (We could get a more accurate k by invoking log10,
2193
* but this is probably not worthwhile.)
2197
#ifndef Sudden_Underflow
2201
/* d is denormalized */
2203
i = bbits + be + (Bias + (P-1) - 1);
2204
x = i > 32 ? word0(d) << (64 - i) | word1(d) >> (i - 32) : word1(d) << (32 - i);
2206
set_word0(d2, word0(d2) - 31*Exp_msk1); /* adjust exponent */
2207
i -= (Bias + (P-1) - 1) + 1;
2211
/* At this point d = f*2^i, where 1 <= f < 2. d2 is an approximation of f. */
2212
ds = (d2-1.5)*0.289529654602168 + 0.1760912590558 + i*0.301029995663981;
2214
if (ds < 0. && ds != k)
2215
k--; /* want k = floor(ds) */
2217
if (k >= 0 && k <= Ten_pmax) {
2222
/* At this point floor(log10(d)) <= k <= floor(log10(d))+1.
2223
If k_check is zero, we're guaranteed that k = floor(log10(d)). */
2225
/* At this point d = b/2^j, where b is an odd integer. */
2244
/* At this point d/10^k = (b * 2^b2 * 5^b5) / (2^s2 * 5^s5), where b is an odd integer,
2245
b2 >= 0, b5 >= 0, s2 >= 0, and s5 >= 0. */
2246
if (mode < 0 || mode > 9)
2268
ilim = ilim1 = i = ndigits;
2274
i = ndigits + k + 1;
2280
/* ilim is the maximum number of significant digits we want, based on k and ndigits. */
2281
/* ilim1 is the maximum number of significant digits we want, based on k and ndigits,
2282
when it turns out that k was computed too high by one. */
2284
/* Ensure space for at least i+1 characters, including trailing null. */
2285
if (bufsize <= (size_t)i) {
2287
JS_ASSERT(JS_FALSE);
2293
if (ilim >= 0 && ilim <= Quick_max && try_quick) {
2295
/* Try to get by with floating-point arithmetic. */
2301
ieps = 2; /* conservative */
2302
/* Divide d by 10^k, keeping track of the roundoff error and avoiding overflows. */
2307
/* prevent overflows */
2309
d /= bigtens[n_bigtens-1];
2312
for(; j; j >>= 1, i++)
2319
else if ((j1 = -k) != 0) {
2320
d *= tens[j1 & 0xf];
2321
for(j = j1 >> 4; j; j >>= 1, i++)
2327
/* Check that k was computed correctly. */
2328
if (k_check && d < 1. && ilim > 0) {
2336
/* eps bounds the cumulative error. */
2338
set_word0(eps, word0(eps) - (P-1)*Exp_msk1);
2348
#ifndef No_leftright
2350
/* Use Steele & White method of only
2351
* generating digits needed.
2353
eps = 0.5/tens[ilim-1] - eps;
2357
*s++ = '0' + (char)L;
2370
/* Generate ilim digits, then fix them up. */
2371
eps *= tens[ilim-1];
2372
for(i = 1;; i++, d *= 10.) {
2375
*s++ = '0' + (char)L;
2379
else if (d < 0.5 - eps) {
2380
while(*--s == '0') ;
2387
#ifndef No_leftright
2397
/* Do we have a "small" integer? */
2399
if (be >= 0 && k <= Int_max) {
2402
if (ndigits < 0 && ilim <= 0) {
2404
if (ilim < 0 || d < 5*ds || (!biasUp && d == 5*ds))
2409
L = (Long) (d / ds);
2411
#ifdef Check_FLT_ROUNDS
2412
/* If FLT_ROUNDS == 2, L will usually be high by 1 */
2418
*s++ = '0' + (char)L;
2421
if ((d > ds) || (d == ds && (L & 1 || biasUp))) {
2444
#ifndef Sudden_Underflow
2445
denorm ? be + (Bias + (P-1) - 1 + 1) :
2448
/* i is 1 plus the number of trailing zero bits in d's significand. Thus,
2449
(2^m2 * 5^m5) / (2^(s2+i) * 5^s5) = (1/2 lsb of d)/10^k. */
2460
if ((i = ilim) < 0) {
2464
/* (2^m2 * 5^m5) / (2^(s2+i) * 5^s5) = (1/2 * 10^(1-ilim))/10^k. */
2471
/* (mhi * 2^m2 * 5^m5) / (2^s2 * 5^s5) = one-half of last printed (when mode >= 2) or
2472
input (when mode < 2) significant digit, divided by 10^k. */
2474
/* We still have d/10^k = (b * 2^b2 * 5^b5) / (2^s2 * 5^s5). Reduce common factors in
2475
b2, m2, and s2 without changing the equalities. */
2476
if (m2 > 0 && s2 > 0) {
2477
i = m2 < s2 ? m2 : s2;
2483
/* Fold b5 into b and m5 into mhi. */
2487
mhi = pow5mult(mhi, m5);
2496
if ((j = b5 - m5) != 0) {
2503
b = pow5mult(b, b5);
2508
/* Now we have d/10^k = (b * 2^b2) / (2^s2 * 5^s5) and
2509
(mhi * 2^m2) / (2^s2 * 5^s5) = one-half of last printed or input significant digit, divided by 10^k. */
2515
S = pow5mult(S, s5);
2519
/* Now we have d/10^k = (b * 2^b2) / (S * 2^s2) and
2520
(mhi * 2^m2) / (S * 2^s2) = one-half of last printed or input significant digit, divided by 10^k. */
2522
/* Check for special case that d is a normalized power of 2. */
2525
if (!word1(d) && !(word0(d) & Bndry_mask)
2526
#ifndef Sudden_Underflow
2527
&& word0(d) & (Exp_mask & Exp_mask << 1)
2530
/* The special case. Here we want to be within a quarter of the last input
2531
significant digit instead of one half of it when the decimal output string's value is less than d. */
2538
/* Arrange for convenient computation of quotients:
2539
* shift left if necessary so divisor has 4 leading 0 bits.
2541
* Perhaps we should just compute leading 28 bits of S once
2542
* and for all and pass them and a shift to quorem, so it
2543
* can do shifts and ors to compute the numerator for q.
2545
if ((i = ((s5 ? 32 - hi0bits(S->x[S->wds-1]) : 1) + s2) & 0x1f) != 0)
2547
/* i is the number of leading zero bits in the most significant word of S*2^s2. */
2560
/* Now S*2^s2 has exactly four leading zero bits in its most significant word. */
2571
/* Now we have d/10^k = b/S and
2572
(mhi * 2^m2) / S = maximum acceptable error, divided by 10^k. */
2576
b = multadd(b, 10, 0); /* we botched the k estimate */
2580
mhi = multadd(mhi, 10, 0);
2587
/* At this point 1 <= d/10^k = b/S < 10. */
2589
if (ilim <= 0 && mode > 2) {
2590
/* We're doing fixed-mode output and d is less than the minimum nonzero output in this mode.
2591
Output either zero or the minimum nonzero output depending on which is closer to d. */
2598
if (i < 0 || (i == 0 && !biasUp)) {
2599
/* Always emit at least one digit. If the number appears to be zero
2600
using the current mode, then emit one '0' digit and set decpt to 1. */
2613
mhi = lshift(mhi, m2);
2618
/* Compute mlo -- check for special case
2619
* that d is a normalized power of 2.
2624
mhi = Balloc(mhi->k);
2628
mhi = lshift(mhi, Log2P);
2632
/* mlo/S = maximum acceptable error, divided by 10^k, if the output is less than d. */
2633
/* mhi/S = maximum acceptable error, divided by 10^k, if the output is greater than d. */
2636
dig = quorem(b,S) + '0';
2637
/* Do we yet have the shortest decimal string
2638
* that will round to d?
2641
/* j is b/S compared with mlo/S. */
2642
delta = diff(S, mhi);
2645
j1 = delta->sign ? 1 : cmp(b, delta);
2647
/* j1 is b/S compared with 1 - mhi/S. */
2648
#ifndef ROUND_BIASED
2649
if (j1 == 0 && !mode && !(word1(d) & 1)) {
2658
if ((j < 0) || (j == 0 && !mode
2659
#ifndef ROUND_BIASED
2664
/* Either dig or dig+1 would work here as the least significant decimal digit.
2665
Use whichever would produce a decimal value closer to d. */
2670
if (((j1 > 0) || (j1 == 0 && (dig & 1 || biasUp)))
2678
if (dig == '9') { /* possible if i == 1 */
2689
b = multadd(b, 10, 0);
2693
mlo = mhi = multadd(mhi, 10, 0);
2698
mlo = multadd(mlo, 10, 0);
2701
mhi = multadd(mhi, 10, 0);
2709
*s++ = (char)(dig = quorem(b,S) + '0');
2712
b = multadd(b, 10, 0);
2717
/* Round off last digit */
2723
if ((j > 0) || (j == 0 && (dig & 1 || biasUp))) {
2734
/* Strip trailing zeros */
2735
while(*--s == '0') ;
2741
if (mlo && mlo != mhi)
2747
JS_ASSERT(s < buf + bufsize);
2758
if (mlo && mlo != mhi)
2771
/* Mapping of JSDToStrMode -> js_dtoa mode */
2772
static const int dtoaModes[] = {
2773
0, /* DTOSTR_STANDARD */
2774
0, /* DTOSTR_STANDARD_EXPONENTIAL, */
2775
3, /* DTOSTR_FIXED, */
2776
2, /* DTOSTR_EXPONENTIAL, */
2777
2}; /* DTOSTR_PRECISION */
2779
JS_FRIEND_API(char *)
2780
JS_dtostr(char *buffer, size_t bufferSize, JSDToStrMode mode, int precision, double d)
2782
int decPt; /* Position of decimal point relative to first digit returned by js_dtoa */
2783
int sign; /* Nonzero if the sign bit was set in d */
2784
int nDigits; /* Number of significand digits returned by js_dtoa */
2785
char *numBegin = buffer+2; /* Pointer to the digits returned by js_dtoa; the +2 leaves space for */
2786
/* the sign and/or decimal point */
2787
char *numEnd; /* Pointer past the digits returned by js_dtoa */
2790
JS_ASSERT(bufferSize >= (size_t)(mode <= DTOSTR_STANDARD_EXPONENTIAL ? DTOSTR_STANDARD_BUFFER_SIZE :
2791
DTOSTR_VARIABLE_BUFFER_SIZE(precision)));
2793
if (mode == DTOSTR_FIXED && (d >= 1e21 || d <= -1e21))
2794
mode = DTOSTR_STANDARD; /* Change mode here rather than below because the buffer may not be large enough to hold a large integer. */
2796
/* Locking for Balloc's shared buffers */
2797
ACQUIRE_DTOA_LOCK();
2798
dtoaRet = js_dtoa(d, dtoaModes[mode], mode >= DTOSTR_FIXED, precision, &decPt, &sign, &numEnd, numBegin, bufferSize-2);
2799
RELEASE_DTOA_LOCK();
2803
nDigits = numEnd - numBegin;
2805
/* If Infinity, -Infinity, or NaN, return the string regardless of the mode. */
2806
if (decPt != 9999) {
2807
JSBool exponentialNotation = JS_FALSE;
2808
int minNDigits = 0; /* Minimum number of significand digits required by mode and precision */
2813
case DTOSTR_STANDARD:
2814
if (decPt < -5 || decPt > 21)
2815
exponentialNotation = JS_TRUE;
2822
minNDigits = decPt + precision;
2827
case DTOSTR_EXPONENTIAL:
2828
JS_ASSERT(precision > 0);
2829
minNDigits = precision;
2831
case DTOSTR_STANDARD_EXPONENTIAL:
2832
exponentialNotation = JS_TRUE;
2835
case DTOSTR_PRECISION:
2836
JS_ASSERT(precision > 0);
2837
minNDigits = precision;
2838
if (decPt < -5 || decPt > precision)
2839
exponentialNotation = JS_TRUE;
2843
/* If the number has fewer than minNDigits, pad it with zeros at the end */
2844
if (nDigits < minNDigits) {
2845
p = numBegin + minNDigits;
2846
nDigits = minNDigits;
2849
} while (numEnd != p);
2853
if (exponentialNotation) {
2854
/* Insert a decimal point if more than one significand digit */
2857
numBegin[0] = numBegin[1];
2860
JS_snprintf(numEnd, bufferSize - (numEnd - buffer), "e%+d", decPt-1);
2861
} else if (decPt != nDigits) {
2862
/* Some kind of a fraction in fixed notation */
2863
JS_ASSERT(decPt <= nDigits);
2865
/* dd...dd . dd...dd */
2873
/* 0 . 00...00dd...dd */
2875
numEnd += 1 - decPt;
2877
JS_ASSERT(numEnd < buffer + bufferSize);
2879
while (p != numBegin)
2881
for (p = numBegin + 1; p != q; p++)
2889
/* If negative and neither -0.0 nor NaN, output a leading '-'. */
2891
!(word0(d) == Sign_bit && word1(d) == 0) &&
2892
!((word0(d) & Exp_mask) == Exp_mask &&
2893
(word1(d) || (word0(d) & Frac_mask)))) {
2900
/* Let b = floor(b / divisor), and return the remainder. b must be nonnegative.
2901
* divisor must be between 1 and 65536.
2902
* This function cannot run out of memory. */
2904
divrem(Bigint *b, uint32 divisor)
2907
uint32 remainder = 0;
2911
JS_ASSERT(divisor > 0 && divisor <= 65536);
2914
return 0; /* b is zero */
2919
ULong dividend = remainder << 16 | a >> 16;
2920
ULong quotientHi = dividend / divisor;
2923
remainder = dividend - quotientHi*divisor;
2924
JS_ASSERT(quotientHi <= 0xFFFF && remainder < divisor);
2925
dividend = remainder << 16 | (a & 0xFFFF);
2926
quotientLo = dividend / divisor;
2927
remainder = dividend - quotientLo*divisor;
2928
JS_ASSERT(quotientLo <= 0xFFFF && remainder < divisor);
2929
*bp = quotientHi << 16 | quotientLo;
2931
/* Decrease the size of the number if its most significant word is now zero. */
2938
/* "-0.0000...(1073 zeros after decimal point)...0001\0" is the longest string that we could produce,
2939
* which occurs when printing -5e-324 in binary. We could compute a better estimate of the size of
2940
* the output string and malloc fewer bytes depending on d and base, but why bother? */
2941
#define DTOBASESTR_BUFFER_SIZE 1078
2942
#define BASEDIGIT(digit) ((char)(((digit) >= 10) ? 'a' - 10 + (digit) : '0' + (digit)))
2944
JS_FRIEND_API(char *)
2945
JS_dtobasestr(int base, double d)
2947
char *buffer; /* The output string */
2948
char *p; /* Pointer to current position in the buffer */
2949
char *pInt; /* Pointer to the beginning of the integer part of the string */
2952
double di; /* d truncated to an integer */
2953
double df; /* The fractional part of d */
2955
JS_ASSERT(base >= 2 && base <= 36);
2957
buffer = (char*) malloc(DTOBASESTR_BUFFER_SIZE);
2961
#if defined(XP_WIN) || defined(XP_OS2)
2962
&& !((word0(d) & Exp_mask) == Exp_mask && ((word0(d) & Frac_mask) || word1(d))) /* Visual C++ doesn't know how to compare against NaN */
2969
/* Check for Infinity and NaN */
2970
if ((word0(d) & Exp_mask) == Exp_mask) {
2971
strcpy(p, !word1(d) && !(word0(d) & Frac_mask) ? "Infinity" : "NaN");
2975
/* Locking for Balloc's shared buffers */
2976
ACQUIRE_DTOA_LOCK();
2978
/* Output the integer part of d with the digits in reverse order. */
2981
if (di <= 4294967295.0) {
2982
uint32 n = (uint32)di;
2985
uint32 m = n / base;
2988
JS_ASSERT(digit < (uint32)base);
2989
*p++ = BASEDIGIT(digit);
2994
int32 bits; /* Number of significant bits in di; not used. */
2995
Bigint *b = d2b(di, &e, &bits);
3005
digit = divrem(b, base);
3006
JS_ASSERT(digit < (uint32)base);
3007
*p++ = BASEDIGIT(digit);
3011
/* Reverse the digits of the integer part of d. */
3021
/* We have a fraction. */
3022
int32 e, bbits, s2, done;
3023
Bigint *b, *s, *mlo, *mhi;
3025
b = s = mlo = mhi = NULL;
3028
b = d2b(df, &e, &bbits);
3039
/* At this point df = b * 2^e. e must be less than zero because 0 < df < 1. */
3041
s2 = -(int32)(word0(d) >> Exp_shift1 & Exp_mask>>Exp_shift1);
3042
#ifndef Sudden_Underflow
3047
/* 1/2^s2 = (nextDouble(d) - d)/2 */
3053
if (!word1(d) && !(word0(d) & Bndry_mask)
3054
#ifndef Sudden_Underflow
3055
&& word0(d) & (Exp_mask & Exp_mask << 1)
3058
/* The special case. Here we want to be within a quarter of the last input
3059
significant digit instead of one half of it when the output string's value is less than d. */
3061
mhi = i2b(1<<Log2P);
3065
b = lshift(b, e + s2);
3074
/* At this point we have the following:
3076
* 1 > df = b/2^s2 > 0;
3077
* (d - prevDouble(d))/2 = mlo/2^s2;
3078
* (nextDouble(d) - d)/2 = mhi/2^s2. */
3085
b = multadd(b, base, 0);
3088
digit = quorem2(b, s2);
3090
mlo = mhi = multadd(mlo, base, 0);
3095
mlo = multadd(mlo, base, 0);
3098
mhi = multadd(mhi, base, 0);
3103
/* Do we yet have the shortest string that will round to d? */
3105
/* j is b/2^s2 compared with mlo/2^s2. */
3106
delta = diff(s, mhi);
3109
j1 = delta->sign ? 1 : cmp(b, delta);
3111
/* j1 is b/2^s2 compared with 1 - mhi/2^s2. */
3113
#ifndef ROUND_BIASED
3114
if (j1 == 0 && !(word1(d) & 1)) {
3120
if (j < 0 || (j == 0
3121
#ifndef ROUND_BIASED
3126
/* Either dig or dig+1 would work here as the least significant digit.
3127
Use whichever would produce an output value closer to d. */
3132
if (j1 > 0) /* The even test (|| (j1 == 0 && (digit & 1))) is not here because it messes up odd base output
3133
* such as 3.5 in base 3. */
3137
} else if (j1 > 0) {
3141
JS_ASSERT(digit < (uint32)base);
3142
*p++ = BASEDIGIT(digit);
3150
JS_ASSERT(p < buffer + DTOBASESTR_BUFFER_SIZE);
3152
RELEASE_DTOA_LOCK();
added
removed
src/extension/script/js/jsdtoa.c
