Bug Summary

File:platform/mac/avmshell/../../../core/MathUtils.cpp
Location:line 1145, column 47
Description:The left operand of '!=' is a garbage value

Annotated Source Code

1/* -*- Mode: C++; c-basic-offset: 4; indent-tabs-mode: nil; tab-width: 4 -*- */
2/* vi: set ts=4 sw=4 expandtab: (add to ~/.vimrc: set modeline modelines=5) */
3/* ***** BEGIN LICENSE BLOCK *****
4 * Version: MPL 1.1/GPL 2.0/LGPL 2.1
5 *
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/
10 *
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
14 * License.
15 *
16 * The Original Code is [Open Source Virtual Machine.].
17 *
18 * The Initial Developer of the Original Code is
19 * Adobe System Incorporated.
20 * Portions created by the Initial Developer are Copyright (C) 2004-2006
21 * the Initial Developer. All Rights Reserved.
22 *
23 * Contributor(s):
24 * Adobe AS3 Team
25 *
26 * Alternatively, the contents of this file may be used under the terms of
27 * either the GNU General Public License Version 2 or later (the "GPL"), or
28 * 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.
37 *
38 * ***** END LICENSE BLOCK ***** */
39
40#include "avmplus.h"
41#include "BigInteger.h"
42
43//GCC only allows intrinsics if sse2 is enabled
44#if (defined(_MSC_VER) || (defined(__GNUC__4) && defined(__SSE2__1))) && (defined(AVMPLUS_IA32) || defined(AVMPLUS_AMD64))
45 #include <emmintrin.h>
46#endif
47
48namespace avmplus
49{
50 const double kLog2_10 = 0.30102999566398119521373889472449;
51
52 // optimize pow(10,x) for commonly sized numbers;
53 static const double kPowersOfTen[23] = { 1, 1e1, 1e2, 1e3, 1e4, 1e5, 1e6, 1e7, 1e8, 1e9, 1e10,
54 1e11, 1e12, 1e13, 1e14, 1e15, 1e16, 1e17, 1e18, 1e19, 1e20,
55 1e21, 1e22 };
56
57 static inline double quickPowTen(int32_t exp)
58 {
59 if (exp < 23 && exp > 0) {
60 return kPowersOfTen[exp];
61 } else {
62 return MathUtils::pow(10,exp);
63 }
64 }
65
66 // skipSpaces: skip whitespace characters
67
68 static int32_t skipSpaces(const StringIndexer& s, int32_t index)
69 {
70 while (index < s->length())
71 {
72 uint32_t ch = s[index];
73 if (!(ch == ' ' || ch == '\t' || ch == '\r' || ch == '\n' || ch == '\f' || ch == '\v' ||
74 (ch >= 0x2000 && ch <=0x200b) || ch == 0x2028 || ch == 0x2029 || ch == 0x205f || ch == 0x3000 ))
75 break;
76 index++;
77 }
78 return index;
79 }
80
81 // handleSign: helper function which checks for sign indicator
82 //
83
84 int32_t handleSign(const StringIndexer& s, int32_t index, boolbool& negative)
85 {
86 negative = falsefalse;
87 if (index >= s->length())
88 return index;
89 uint32_t ch = s[index];
90 if (ch == '+') {
91 index++;
92 } else if (ch == '-') {
93 negative = truetrue;
94 index++;
95 }
96 return index;
97 }
98
99 boolbool MathUtils::equals(double x, double y)
100 {
101 // Infiniti-ness must match up
102 if (isInfinite(x) != isInfinite(y)) {
103 return falsefalse;
104 }
105 // Nothing is equal to NaN, not even itself.
106 if (isNaN(x) || isNaN(y)) {
107 return falsefalse;
108 }
109 return x == y;
110 }
111
112 static double _nan()
113 {
114 union { float f; uint32_t d; }; d = 0x7FFFFFFF;
115 return f;
116 }
117
118 static double _infinity()
119 {
120 union { float f; uint32_t d; }; d = 0x7F800000;
121 return f;
122 }
123
124 static double _neg_infinity()
125 {
126 union { float f; uint32_t d; }; d = 0xFF800000;
127 return f;
128 }
129
130
131 /*static*/ const double MathUtils::kNaN = _nan();
132 /*static*/ const double MathUtils::kInfinity = _infinity();
133 /*static*/ const double MathUtils::kNegInfinity = _neg_infinity();
134
135#ifdef UNIX
136 /*
137 * MathUtils::isNaN and MathUtils::isInfinite are modified versions of
138 * s_isNaN.c and s_isInfinite.c, freely usable code by Sun. Copyright follows:
139 *======================================================
140 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
141 *
142 * Developed at SunPro, a Sun Microsystems, Inc. business.
143 * Permission to use, copy, modify and distribute this
144 * software is freely granted, provided that this notice is preserved.
145 *======================================================
146 */
147
148 int32_t MathUtils::isInfinite(double x)
149 {
150 double_overlay u(x);
151
152 int32_t hx = u.words.msw;
153 int32_t lx = u.words.lsw;
154
155 lx |= (hx & 0x7FFFFFFF) ^ 0x7FF00000;
156 lx |= -lx;
157 return ~(lx >> 31) & (hx >> 30);
158 }
159
160 boolbool MathUtils::isNaN(double x)
161 {
162 // Infinity is NOT considered NaN in
163 // JavaScript
164 if (MathUtils::isInfinite(x) != 0) {
165 return falsefalse;
166 }
167
168 double_overlay u(x);
169
170 int32_t hx = u.words.msw;
171 int32_t lx = u.words.lsw;
172
173 hx &= 0x7FFFFFFF;
174 hx |= (uint32_t) (lx | (-lx)) >> 31;
175 hx = 0x7FF00000 - hx;
176 return (boolbool) (((uint32_t)(hx) >> 31) != 0);
177 }
178
179 boolbool MathUtils::isNegZero(double x)
180 {
181 double_overlay u(x);
182
183 int32_t hx = u.words.msw;
184 int32_t lx = u.words.lsw;
185
186 return (hx == (int32_t)0x80000000 && lx == (int32_t)0x0);
187 }
188
189#else
190
191 /*
192 * double precision special value tests:
193
194 inf = 7FF0 0000 0000 0000
195 FFF0 0000 0000 0000
196
197 nan = 7FFx xxxx xxxx xxxx where x != 0
198 FFFx xxxx xxxx xxxx
199 */
200
201 int32_t MathUtils::isInfinite(double x)
202 {
203 double_overlay u(x);
204
205 if ((u.bits64 & 0x7fffffffffffffffLL) != 0x7FF0000000000000LL)
206 return 0;
207 if (u.bits64 & 0x8000000000000000LL)
208 return -1;
209 else
210 return 1;
211 }
212
213 boolbool MathUtils::isNaN(double value)
214 {
215 double_overlay u(value);
216
217 return (u.bits64 & ~0x8000000000000000ULL) > 0x7ff0000000000000ULL;
218 }
219
220 boolbool MathUtils::isNegZero(double x)
221 {
222 double_overlay d(x);
223
224 return d.bits64 == 0x8000000000000000ULL;
225 }
226
227#endif // UNIX
228
229 int32_t MathUtils::nextPowerOfTwo(int32_t n)
230 {
231 int32_t i = 2;
232 while (i < n)
233 {
234 i <<= 1;
235 }
236
237 return (i);
238 }
239
240 static boolbool isHexNumber(const StringIndexer& s, int32_t index)
241 {
242 if (s->length() - index < 2 || s[index] != '0')
243 return falsefalse;
244 uint32_t ch = s[index+1];
245 return (ch == 'x' || ch == 'X');
246 }
247
248 // parseIntDigit: Helper for ParseInt and ConvertStringToInteger
249 //
250
251 static int32_t parseIntDigit(wchar ch)
252 {
253 if (ch >= '0' && ch <= '9') {
254 return (ch - '0');
255 } else if (ch >= 'a' && ch <= 'z') {
256 return (ch - 'a' + 10);
257 } else if (ch >= 'A' && ch <= 'Z') {
258 return (ch - 'A' + 10);
259 } else {
260 return -1;
261 }
262 }
263
264 // parseInt: Implementation of ECMA-262 parseInt function
265 double MathUtils::parseInt(Stringp inStr, int32_t radix /*=10*/, boolbool strict /*=false*/ )
266 {
267 boolbool negate;
268 boolbool gotDigits = falsefalse;
269 double result = 0;
270
271 StringIndexer s(inStr);
272 int32_t index = skipSpaces(s, 0); // leading and trailing whitespace is valid.
273 index = handleSign(s, index, negate);
274 if (isHexNumber(s, index) && (radix == 16 || radix == 0) ) {
275 index += 2;
276 if (radix == 0)
277 radix = 16;
278 } else if (radix == 0) {
279 // default radix is 10
280 radix = 10;
281 }
282
283 // Make sure radix is valid, and we have digits
284 if (radix >= 2 && radix <= 36 && index < s->length()) {
285 result = 0;
286 int32_t start = index;
287
288 // Read the digits, generate result
289 while (index < s->length()) {
290 int32_t v = parseIntDigit(s[index]);
291 if (v == -1 || v >= radix) {
292 break;
293 }
294 result = result * radix + v;
295 gotDigits = truetrue;
296 index++;
297 }
298
299 index = skipSpaces(s, index); // leading and trailing whitespace is valid.
300 if (strict && index < s->length()) {
301 return MathUtils::kNaN;
302 }
303
304 if ( result >= 0x20000000000000LL && // i.e. if the result may need at least 54 bits of mantissa
305 (radix == 2 || radix == 4 || radix == 8 || radix == 16 || radix == 32) ) {
306
307 // CN: we're here because we may have incurred roundoff error with the above.
308 // Error will creep in once we need more than the available 53 bits
309 // of precision in the mantissa portion of a double. No way to deduce
310 // this from the result, so we have to recalculate it more slowly.
311 result = 0;
312
313 int32_t powOf2 = 1;
314 for(int32_t x = radix; (x != 1); x >>= 1)
315 powOf2++;
316 powOf2--; // each word contains one less than this # of bits.
317 index = start;
318 int32_t v=0;
319
320 int32_t end,next;
321 // skip leading zeros
322 for(end=index; end < s->length() && s[end] == '0'; end++)
323 ;
324 if (end >= s->length())
325 return 0;
326
327 for (next=0; next*powOf2 <= 52; next++) { // read first 52 non-zero digits. Once charPosition*log2(radix) is > 53, we can have rounding issues
328 v = parseIntDigit(s[end++]);
329 if (v == -1 || v >= radix) {
330 v = 0;
331 break;
332 }
333 result = result * radix + v;
334 if (end >= s->length())
335 break;
336 }
337 if (next*powOf2 > 52) { // If number contains more than 53 bits of precision, may need to roundUp last digit processed.
338 boolbool roundUp = falsefalse;
339 int32_t bit53 = 0;
340 int32_t bit54 = 0;
341
342 double factor = 1;
343
344 switch(radix) {
345 case 32: // last word read contained digits 51,52,53,54,55
346 bit53 = v & (1 << 2);
347 bit54 = v & (1 << 1);
348 roundUp = (v & 1);
349 break;
350 case 16: // last word read contained digits 50,51,52,53
351 bit53 = v & (1 << 0);
352 v = parseIntDigit(s[end]);
353 if (v != -1 && v < radix) {
354 factor *= radix;
355 bit54 = v & (1 << 3);
356 roundUp = (v & 0x3) != 0; // check if any bit after bit54 is set
357 } else {
358 roundUp = bit53 != 0;
359 }
360 break;
361 case 8: // last work read contained digits 49,50,51, next word contains 52,53,54
362 v = parseIntDigit(s[end]);
363 if (v == -1 || v >= radix) {
364 v = 0;
365 }
366 factor *= radix;
367 bit53 = v & (1 << 1);
368 bit54 = v & (1 << 0);
369 break;
370 case 4: // 51,52 - 53,54
371 v = parseIntDigit(s[end]);
372 if (v == -1 || v >= radix) {
373 v = 0;
374 }
375 factor *= radix;
376 bit53 = v & (1 << 1);
377 bit54 = v & (1 << 0);
378 break;
379 case 2: // 52 - 53 - 54
380 /*
381 v = parseIntDigit(s[end++]);
382 result = result * radix; // add factor before round-off adjustment for 52 bit
383 */
384 bit53 = v & (1 << 0);
385 v = parseIntDigit(s[end]);
386 if (v != -1 && v < radix) {
387 factor *= radix;
388 bit54 = (v != -1 ? (v & (1 << 0)) : 0); // Might be there are only 53 digits.
389 }
390
391 break;
392 }
393
394 bit53 = !!bit53;
395 bit54 = !!bit54;
396
397
398 while(++end < s->length()) {
399 v = parseIntDigit(s[end]);
400 if (v == -1 || v >= radix) {
401 break;
402 }
403 roundUp |= (v != 0); // any trailing positive bit causes us to round up
404 factor *= radix;
405 }
406 roundUp = bit54 && (bit53 || roundUp);
407 result += (roundUp ? 1.0 : 0.0);
408 result *= factor;
409 }
410
411 }
412 /*
413 else if (radix == 10 && result >= 0x20000000000000)
414 // if there are more than 15 digits, roundoff error may affect us. Need to use exact integer rep instead of float
415 //int32_t numDigits = len - (s - sStart);
416 */
417 if (negate) {
418 result = -result;
419 }
420 }
421 return gotDigits ? result : MathUtils::kNaN;
422 }
423
424 // used by AvmCore::number() and numberAtom() for converting arbitrary string to a number.
425 double MathUtils::convertStringToNumber(Stringp inStr)
426 {
427 double value;
428
429 if (inStr->length() == 0) { // toNumber("") should be 0, not NaN
430 value = 0;
431 } else if (MathUtils::convertStringToDouble(inStr, &value, truetrue)) {
432 // nothing to do, value set by side effect.
433 } else {
434 value = MathUtils::parseInt(inStr, 0,truetrue); // 0 radix means figure it out from the string.
435 }
436 return value;
437 }
438
439 // apparently SunPro compiler doesn't like combining REALLY_INLINE with static functions.
440 /*static*/
441 REALLY_INLINEinline __attribute__((always_inline)) int32_t real_to_int(double value)
442 {
443#if defined(WIN32) && defined(AVMPLUS_AMD64)
444 int32_t intValue = _mm_cvttsd_si32(_mm_set_sd(value));
445#elif defined(WIN32) && defined(AVMPLUS_IA32)
446 int32_t intValue;
447 _asm fld [value];
448 _asm fistp [intValue];
449#elif defined(_MAC1) && (defined(AVMPLUS_IA32) || defined(AVMPLUS_AMD64))
450 int32_t intValue = _mm_cvttsd_si32(_mm_set_sd(value));
451#else
452 int32_t intValue = int32_t(value);
453#endif
454 return intValue;
455 }
456
457 // MathUtils::toInt is the ToInteger algorithm from
458 // ECMA-262 section 9.4
459 double MathUtils::toInt(double value)
460 {
461 int32_t intValue = real_to_int(value);
462
463 if ((value == (double)(intValue)) && ((uint32_t)intValue != 0x80000000))
464 return value;
465
466 if (MathUtils::isNaN(value)) {
467 return 0;
468 }
469 if (MathUtils::isInfinite(value) || value == 0) {
470 return value;
471 }
472 if (value < 0) {
473 return -MathUtils::floor(-value);
474 } else {
475 return MathUtils::floor(value);
476 }
477 }
478
479 //
480 // powerOfTen(exponent, value) returns the value
481 // value * 10 ^ exponent
482 //
483
484 double MathUtils::powerOfTen(int32_t exponent, double value)
485 {
486 double base = 10.0;
487
488 if (exponent < 0) {
489 exponent = -exponent;
490 while (exponent) {
491 if (exponent & 1) {
492 value /= base;
493 }
494 exponent >>= 1;
495 base *= base;
496 }
497 } else {
498 while (exponent) {
499 if (exponent & 1) {
500 value *= base;
501 }
502 exponent >>= 1;
503 base *= base;
504 }
505 }
506 return value;
507 }
508
509 //
510 // powInt(base, exponent) returns the value
511 // base ^ exponent
512 // where exponent is an integer.
513 //
514 static double powInt(double base, int32_t exponent)
515 {
516 double value = 1.0;
517 double original_base = base;
518 int32_t original_exponent = exponent;
519
520 if (MathUtils::isInfinite(base))
521 {
522 if (exponent < 0) // return -0 or 0 depending on sign of base
523 return (base < 0) ? 1.0/base : 0.0;
524
525 // if base is -Infinity, return Infinity or -Infinity depending whether the exponent is even or not.
526 if (base < 0)
527 return (MathUtils::mod(exponent, 2) != 0) ? base : 0-base;
528
529 return base;
530 }
531
532 if (exponent < 0) {
533 exponent = -exponent;
534 while (exponent) {
535 if (exponent & 1) {
536 value /= base;
537 // cn: max double value is magnitude 1e308 (base 10), but min double value is magnitude 1e-324
538 // That means we can't invert the exponent when the base ten value exponent is < 307.
539 // We could first check if powerOfTwo(base)*exp > 1074, but that would
540 // slow down all powInts calls. This limits the xtra work to just very small numbers.
541 if (value == 0 && base != 0) // double check by calling the real thing
542 return MathUtils::powInternal(original_base,(double)original_exponent);
543 }
544 exponent >>= 1;
545 base *= base;
546 }
547 } else {
548 while (exponent) {
549 if (exponent & 1) {
550 value *= base;
551 }
552 exponent >>= 1;
553 base *= base;
554 }
555 }
556 return value;
557 }
558
559 // x = base, y = exponent
560 double MathUtils::pow(double x, double y)
561 {
562 if (MathUtils::isNaN(y)) {
563 // x^NaN = NaN
564 return MathUtils::kNaN;
565 }
566
567 if (y == 0) {
568 // x^0 is 1 for all x
569 return 1;
570 }
571
572 int32_t infy = MathUtils::isInfinite(y);
573
574 if (infy == 0 && y == (int32_t)y) {
575 // If y is an integer, use multiplication
576 // algorithm which yields more accurate
577 // results
578 return powInt(x, (int32_t)y);
579 }
580
581 double absx = MathUtils::abs(x);
582
583 if (absx < 1) {
584 infy = -infy;
585 }
586 if (absx == 1 && infy != 0) {
587 // (+/-)1^(+/-)Infinity = NaN
588 return MathUtils::kNaN;
589 }
590 if (infy == 1) {
591 // x^Infinity = Infinity
592 return MathUtils::kInfinity;
593 } else if (infy == -1) {
594 // x^-Infinity = +0
595 return +0;
596 }
597
598 if (MathUtils::isInfinite(x))
599 {
600 if (y < 0) {
601 // y<0
602 // Infinity^y = 0
603 // -Infinity^y = 0
604 return 0;
605 } else {
606 if (y < 1.0) {
607 return MathUtils::kInfinity;
608 } else {
609 // y>0
610 // Infinity^y = Infinity
611 // -Infinity^y = Infinity
612 return x;
613 }
614 }
615 }
616
617 // R41
618 // Handle negative x.
619 if (x < 0.0) {
620 // If y is non-integer, return NaN.
621 if (y != MathUtils::floor(y)) {
622 return MathUtils::kNaN;
623 }
624 // Switch sign of base
625 x = -x;
626 // If exponent is odd, negate result
627 if (MathUtils::mod(y, 2) != 0) {
628 return -powInternal(x, y);
629 }
630 // Else do the usual thing.
631 }
632 // end R41
633
634 // Bug 33811 - Windows math returns NaN incorrectly when base = 0.0
635 if (x == 0.0)
636 {
637 if (y < 0.0)
638 return MathUtils::kInfinity;
639 else
640 return 0.0;
641 }
642
643 return powInternal(x, y);
644 }
645
646 int32_t MathUtils::nextDigit(double *value)
647 {
648 int32_t digit;
649
650 #if defined(WIN32) && defined(AVMPLUS_AMD64)
651 digit = _mm_cvttsd_si32(_mm_set_sd(*value));
652 #elif defined(WIN32) && defined(AVMPLUS_IA32)
653 // WARNING! nextDigit assumes that rounding mode is
654 // set to truncate.
655 _asm mov eax,[value];
656 _asm fld qword ptr [eax];
657 _asm fistp [digit];
658 #elif defined(_MAC1) && (defined(AVMPLUS_IA32) || defined(AVMPLUS_AMD64))
659 digit = _mm_cvttsd_si32(_mm_set_sd(*value));
660 #else
661 digit = (int32_t) *value;
662 #endif
663
664 *value -= digit;
665 *value *= 10;
666 return digit;
667 }
668
669 int32_t MathUtils::roundInt(double x)
670 {
671 if (x < 0) {
672 return (int32_t) (x - 0.5);
673 } else {
674 return (int32_t) (x + 0.5);
675 }
676 }
677
678 Stringp MathUtils::convertIntegerToStringRadix(AvmCore* core,
679 intptr_t value,
680 int32_t radix,
681 UnsignedTreatment treatAs)
682 {
683 MMgc::GC::AllocaAutoPtr _buffer;
684 char* buffer = (char*)VMPI_alloca(core, _buffer, kMinSizeForInt64_t_toString)(kMinSizeForInt64_t_toString > 4000 ? core->gc->allocaPush
(kMinSizeForInt64_t_toString, _buffer) : __builtin_alloca(kMinSizeForInt64_t_toString
));
685 int32_t len = kMinSizeForInt64_t_toString;
686 char* p = convertIntegerToStringBuffer(value, buffer, len, radix, treatAs);
687 return core->newStringLatin1(p, len);
688 }
689
690 Stringp MathUtils::convertIntegerToStringBase10(AvmCore* core,
691 int32_t value,
692 UnsignedTreatment treatAs)
693 {
694 char buffer[kMinSizeForInt32_t_base10_toString];
695 int32_t len = kMinSizeForInt32_t_base10_toString;
696 #ifdef AVMPLUS_64BIT
697 intptr_t wideVal = treatAs == kTreatAsUnsigned ? (intptr_t)(uint32_t)value : (intptr_t)value;
698 #else
699 intptr_t wideVal = (intptr_t) value;
700 #endif
701 char* p = convertIntegerToStringBuffer(wideVal, buffer, len, 10, treatAs);
702 return core->newStringLatin1(p, len);
703 }
704
705 char* MathUtils::convertIntegerToStringBuffer(intptr_t value,
706 char *buffer,
707 int32_t& len,
708 int32_t radix,
709 UnsignedTreatment treatAs)
710 {
711 #ifndef AVMPLUS_64BIT
712 // This routines does not work with this integer number since negating
713 // this value returns the same negative value and screws up our code
714 // in the while loop below.
715 if ( (uint32_t)value == 0x80000000 && treatAs == kTreatAsSigned)
716 // MathUtils::convertIntegerToString doesn't deal with this number because you can't negate it.
717 {
718 AvmAssert(len >= 12)do { } while (0);
719 if (len < 12)
720 return NULL__null;
721 VMPI_memcpy::memcpy (buffer, "-2147483648", 12);
722 len = 11;
723 return buffer;
724 }
725 #endif
726
727 if (radix < 2 || radix > 36)
728 {
729 AvmAssert( 0 )do { } while (0);
730 return NULL__null;
731 }
732
733 char *src = buffer + len - 1;
734 char *srcEnd = src;
735
736 *src-- = '\0';
737
738 if (value == 0)
739 {
740 *src-- = '0';
741 }
742 else
743 {
744 uintptr_t uVal;
745 boolbool negative=falsefalse;
746
747 if (treatAs == kTreatAsUnsigned)
748 {
749 uVal = (uintptr_t)value;
750 }
751 else
752 {
753 negative = (value < 0);
754 if (negative)
755 value = -value;
756 uVal = (uintptr_t)value;
757 }
758
759 while (uVal != 0)
760 {
761 // buffer too small?
762 AvmAssert(src >= buffer)do { } while (0);
763 uintptr_t j = uVal;
764 uVal = uVal / radix;
765 j -= (uVal * radix);
766
767 *src-- = (char)((j < 10) ? (j + '0') : (j + ('a' - 10)));
768 }
769
770 if (negative)
771 {
772 AvmAssert(src >= buffer)do { } while (0);
773 if (src < buffer)
774 // buffer too small
775 return NULL__null;
776 *src-- = '-';
777 }
778 }
779
780 src++;
781 len = (int32_t)(srcEnd-src);
782
783 return src;
784 }
785
786 Stringp MathUtils::convertDoubleToStringRadix(AvmCore *core,
787 double value,
788 int32_t radix)
789 {
790 if (radix < 2 || radix > 36)
791 {
792 AvmAssert( 0 )do { } while (0);
793 return NULL__null;
794 }
795
796 MMgc::GC::AllocaAutoPtr _tmp;
797 char* tmp = (char*)VMPI_alloca(core, _tmp, kMinSizeForDouble_base2_toString)(kMinSizeForDouble_base2_toString > 4000 ? core->gc->
allocaPush(kMinSizeForDouble_base2_toString, _tmp) : __builtin_alloca
(kMinSizeForDouble_base2_toString));
798 char *src = tmp + kMinSizeForDouble_base2_toString - 1;
799 char *srcEnd = src;
800
801 boolbool negative=falsefalse;
802
803 negative = (value < 0);
804 if (negative)
805 value = -value;
806
807 if (value < 1)
808 {
809 *src-- = '0';
810 }
811 else
812 {
813 double uVal = MathUtils::floor(value);
814
815 while (uVal != 0)
816 {
817 double j = uVal;
818 uVal = MathUtils::floor(uVal / radix);
819 j -= (uVal * radix);
820
821 *src-- = (char)((j < 10) ? ((int32_t)j + '0') : ((int32_t)j + ('a' - 10)));
822
823 AvmAssert(src >= tmp || 0 == uVal)do { } while (0);
824 }
825
826 if (negative)
827 *src-- = '-';
828 }
829
830 int32_t len = (int32_t)(srcEnd-src);
831 return core->newStringLatin1(src+1, len);
832 }
833
834 Stringp MathUtils::convertDoubleToString(AvmCore* core,
835 double value,
836 int32_t mode,
837 int32_t precision)
838 {
839 switch (isInfinite(value)) {
1
'Default' branch taken. Execution continues on line 845
840 case -1:
841 return core->newConstantStringLatin1("-Infinity");
842 case 1:
843 return core->newConstantStringLatin1("Infinity");
844 }
845 if (isNaN(value)) {
2
Taking false branch
846 return core->newConstantStringLatin1("NaN");
847 }
848
849 // If our double is really an integer, call our integer version which
850 // is much, much faster then our double version. (Skips a ton of _ftol
851 // which are slow).
852 if (mode == DTOSTR_NORMAL) {
3
Taking false branch
853#if defined(WIN32) && defined(AVMPLUS_AMD64)
854 int32_t intValue = _mm_cvttsd_si32(_mm_set_sd(value));
855#elif defined(WIN32) && defined(AVMPLUS_IA32)
856 int32_t intValue;
857 _asm fld [value];
858 _asm fistp [intValue];
859#elif defined(_MAC1) && (defined(AVMPLUS_IA32) || defined(AVMPLUS_AMD64))
860 int32_t intValue = _mm_cvttsd_si32(_mm_set_sd(value));
861#else
862 int32_t intValue = int32_t(value);
863#endif
864 if ((value == (double)(intValue)) && ((uint32_t)intValue != 0x80000000))
865 return convertIntegerToStringBase10(core, intValue, kTreatAsSigned);
866 }
867
868 int32_t i, len = 0;
869
870 MMgc::GC::AllocaAutoPtr _buffer;
871 char* buffer = (char*)VMPI_alloca(core, _buffer, kMinSizeForDouble_base10_toString)(kMinSizeForDouble_base10_toString > 4000 ? core->gc->
allocaPush(kMinSizeForDouble_base10_toString, _buffer) : __builtin_alloca
(kMinSizeForDouble_base10_toString));
872 char *s = buffer;
873 boolbool negative = falsefalse;
874
875 // Deal with negative numbers
876 if (value < 0) {
4
Taking false branch
877 value = -value;
878 negative = truetrue;
879 // make room for minus sign later (after applying sentinel)
880 s++;
881 }
882
883 // initialize d2a engine
884 D2A *d2a = mmfx_new(D2A(value, mode, precision))new (MMgc::kUseFixedMalloc) D2A(value, mode, precision);
885 int32_t exp10 = d2a->expBase10()-1;
886
887 // Sentinel is used for rounding
888 char *sentinel = s;
889
890 enum {
891 kNormal,
892 kExponential,
893 kFraction,
894 kFixedFraction
895 };
896 int32_t format;
897 int32_t digit;
898
899 switch (mode) {
5
Control jumps to 'case DTOSTR_EXPONENTIAL:' at line 921
900 case DTOSTR_FIXED:
901 {
902 if (exp10 < 0) {
903 format = kFixedFraction;
904 } else {
905 format = kNormal;
906 precision++;
907 }
908 }
909 break;
910 case DTOSTR_PRECISION:
911 {
912 if (exp10 < 0) {
913 format = kFraction;
914 } else if (exp10 >= precision) {
915 format = kExponential;
916 } else {
917 format = kNormal;
918 }
919 }
920 break;
921 case DTOSTR_EXPONENTIAL:
922 format = kExponential;
923 precision++;
924 break;
6
 Execution continues on line 951
925 default:
926 if (exp10 < 0 && exp10 > -7) {
927 // Number is of form 0.######
928 if (exp10 < -precision) {
929 exp10 = -precision-1;
930 }
931 format = kFraction;
932 } else if (exp10 > 20) { // ECMA spec 9.8.1
933 format = kExponential;
934 } else {
935 format = kNormal;
936 }
937 }
938
939 #if 0 // ifdef WIN32
940 // On Windows, set the FPU control word to round
941 // down for the rest of this operation.
942 volatile uint16_t oldcw;
943 volatile uint16_t newcw;
944 _asm fnstcw [oldcw];
945 _asm mov ax,[oldcw];
946 _asm or ax,0xc3f;
947 _asm mov [newcw],ax;
948 _asm fldcw [newcw];
949 #endif
950
951 boolbool wroteDecimal = falsefalse;
952 switch (format) {
7
Control jumps to 'case kExponential:' at line 1054
953 case kNormal:
954 {
955 int32_t digits = 0;
956 sentinel = s;
957 *s++ = '0';
958 digit = d2a->nextDigit();
959 if (digit > 0) {
960 *s++ = (char)(digit + '0');
961 }
962 while (exp10 > 0) {
963 digit = (d2a->finished) ? 0 : d2a->nextDigit();
964 *s++ = (char)(digit + '0');
965 exp10--;
966 digits++;
967 }
968 if (mode == DTOSTR_FIXED) {
969 digits = 0;
970 }
971 if (mode == DTOSTR_NORMAL) {
972 if (!d2a->finished) {
973 *s++ = '.';
974 wroteDecimal = truetrue;
975 while( !d2a->finished ) {
976 *s++ = (char)(d2a->nextDigit() + '0');
977 }
978 }
979 }
980 else if (digits < precision-1)
981 {
982 *s++ = '.';
983 wroteDecimal = truetrue;
984 for (; digits < precision-1; digits++) {
985 digit = d2a->finished ? 0 : d2a->nextDigit();
986 *s++ = (char)(digit + '0');
987 }
988 }
989 }
990 break;
991 case kFixedFraction:
992 {
993 *s++ = '0'; // Sentinel
994 *s++ = '0';
995 *s++ = '.';
996 wroteDecimal = truetrue;
997 // Write out leading zeros
998 int32_t digits = 0;
999 if (exp10 > 0) {
1000 while (++exp10 < 10 && digits < precision) {
1001 *s++ = '0';
1002 digits++;
1003 }
1004 } else if (exp10 < 0) {
1005 while ((++exp10 < 0) && (precision-- > 0))
1006 *s++ = '0';
1007 }
1008
1009 // Write out significand
1010 for ( ; digits<precision; digits++) {
1011 if (d2a->finished)
1012 {
1013 if (mode == DTOSTR_NORMAL)
1014 break;
1015 *s++ = '0';
1016 } else {
1017 *s++ = (char)(d2a->nextDigit() + '0');
1018 }
1019 }
1020 exp10 = 0;
1021 }
1022 break;
1023 case kFraction:
1024 {
1025 *s++ = '0'; // Sentinel
1026 *s++ = '0';
1027 *s++ = '.';
1028 wroteDecimal = truetrue;
1029
1030 // Write out leading zeros
1031 if (value)
1032 {
1033 for (i=exp10; i<-1; i++) {
1034 *s++ = '0';
1035 }
1036 }
1037
1038 // Write out significand
1039 i=0;
1040 while (!d2a->finished)
1041 {
1042 *s++ = (char)(d2a->nextDigit() + '0');
1043 if (mode != DTOSTR_NORMAL && ++i >= precision)
1044 break;
1045 }
1046 if (mode == DTOSTR_PRECISION)
1047 {
1048 while(i++ < precision)
1049 *s++ = (char)(d2a->finished ? '0' : d2a->nextDigit() + '0');
1050 }
1051 exp10 = 0;
1052 }
1053 break;
1054 case kExponential:
1055 {
1056 digit = d2a->finished ? 0 : d2a->nextDigit();
8
'?' condition is false
1057 *s++ = (char)(digit + '0');
1058 if ( ((mode == DTOSTR_NORMAL) && !d2a->finished) ||
9
Taking false branch
1059 ((mode != DTOSTR_NORMAL) && precision > 1) ) {
1060 *s++ = '.';
1061 wroteDecimal = truetrue;
1062 for (i=0; i<precision-1; i++) {
1063 if (d2a->finished)
1064 {
1065 if (mode == DTOSTR_NORMAL)
1066 break;
1067 *s++ = '0';
1068 } else {
1069 *s++ = (char)(d2a->nextDigit() + '0');
1070 }
1071 }
1072 }
1073 }
1074 break;
10
 Execution continues on line 1079
1075 }
1076
1077 // Rounding (todo: argh, was hoping to get rid of this, but we still need it for fastEstimate mode)
1078 // fix bug 121952: must also do rounding for fixed mode
1079 if (d2a->bFastEstimateOk || mode == DTOSTR_FIXED || mode == DTOSTR_PRECISION)
11
Taking false branch
1080 {
1081 i = d2a->nextDigit();
1082 if (i > 4) {
1083 char *ptr = s-1;
1084 while (ptr >= buffer) {
1085 if (*ptr < '0') {
1086 ptr--;
1087 continue;
1088 }
1089 (*ptr)++;
1090 if (*ptr != 0x3A) {
1091 break;
1092 }
1093 *ptr-- = '0';
1094 }
1095 }
1096 if (mode == MathUtils::DTOSTR_NORMAL && wroteDecimal) {
1097 // Remove trailing zeros
1098 while (*(s-1) == '0') {
1099 s--;
1100 }
1101 if (*(s-1) == '.') {
1102 s--;
1103 }
1104 }
1105 }
1106
1107 if (exp10) {
12
Taking false branch
1108 char *firstNonZero = buffer + int(negative);
1109 while (firstNonZero < s && *firstNonZero == '0') {
1110 firstNonZero++;
1111 }
1112 if (s == firstNonZero) {
1113 // Deal with case where all digits were
1114 // rounded away
1115 *s++ = '1';
1116 exp10++;
1117 } else {
1118 char *lastNonZero = s;
1119 while (lastNonZero > firstNonZero) {
1120 if (*--lastNonZero != '0') {
1121 break;
1122 }
1123 }
1124 if (value && (firstNonZero == lastNonZero) ) {
1125 // Watch out for extra zeros like
1126 // 10e+95
1127 exp10 += (int32_t) (s - firstNonZero - 1);
1128 s = lastNonZero+1;
1129 }
1130 }
1131 *s++ = 'e';
1132 if (exp10 > 0) {
1133 *s++ = '+';
1134 }
1135 char expstr[kMinSizeForInt32_t_base10_toString];
1136 len = kMinSizeForInt32_t_base10_toString;
1137 char* t = convertIntegerToStringBuffer(exp10, expstr, len, 10, kTreatAsSigned);
1138 while (*t) { *s++ = *t++; }
1139 }
1140
1141 len = (int32_t)(s-buffer);
1142 s = buffer;
1143 if (negative)
13
Taking false branch
1144 s++;
1145 if (sentinel && sentinel[0] == '0' && sentinel[1] != '.') {
14
The left operand of '!=' is a garbage value
1146 s = sentinel + 1;
1147 len--;
1148 }
1149 if (negative)
1150 *--s = '-';
1151
1152 #if 0 // def WIN32
1153 // Restore control word
1154 _asm fldcw [oldcw];
1155 #endif
1156
1157 mmfx_delete(d2a)::MMgcDestructTaggedScalarChecked(d2a);
1158
1159 return core->newStringLatin1(s, len);
1160 }
1161
1162 // convertStringToDouble: Converts an ASCII string in the form
1163 // [+-]######.######e[+-]###
1164 // to a double-precision floating-point number
1165 //
1166 boolbool MathUtils::convertStringToDouble(Stringp inStr,
1167 double *value,
1168 boolbool strict /*=false*/ )
1169 {
1170 boolbool negate;
1171 int32_t numDigits = 0;
1172 int32_t exp10 = 0;
1173 double result = 0;
1174 uint32_t ch = 0;
1175
1176 StringIndexer s(inStr);
1177 int32_t index = skipSpaces(s, 0);
1178
1179 if (index >= s->length()) // empty string or string of spaces == 0
1180 {
1181 *value = 0;
1182 return (strict ? truetrue : falsefalse); // odd use of strict, but Number(' ') is 0 and that uses strict == true. parseFloat(' ') is NaN and that uses strict == false
1183 }
1184
1185 index = handleSign(s, index, negate);
1186
1187 // Pass 1: Validate input and determine exponents.
1188 int32_t start = index;
1189 int32_t slength = s->length();
1190 while (index < slength) {
1191 ch = s[index];
1192 if (ch >= '0' && ch <= '9') {
1193 numDigits++;
1194 index++;
1195 }
1196 else
1197 {
1198 // https://bugzilla.mozilla.org/show_bug.cgi?id=564839 -- stop parsing at null char
1199 if (ch == 0)
1200 {
1201 slength = index;
1202 }
1203 break;
1204 }
1205 }
1206
1207 // If decimal point encountered, read past digits
1208 // to right of decimal point.
1209 if (ch == '.') {
1210 while (++index < slength) {
1211 ch = s[index];
1212 if (ch >= '0' && ch <= '9') {
1213 numDigits++;
1214 }
1215 else
1216 {
1217 // https://bugzilla.mozilla.org/show_bug.cgi?id=564839 -- stop parsing at null char
1218 if (ch == 0 && index < slength)
1219 {
1220 slength = index;
1221 }
1222 break;
1223 }
1224 }
1225 }
1226
1227 // Handle exponent
1228 if (index < slength && (s[index] == 'e' || s[index] == 'E')) {
1229 int32_t num = 0;
1230 boolbool negexp;
1231 index = handleSign(s, ++index, negexp);
1232 if (negexp && index >= slength) // fail if string ends in "e-" with no exponent value specified.
1233 return falsefalse;
1234
1235 while (index < slength) {
1236 ch = s[index];
1237 if (ch >= '0' && ch <= '9') {
1238 num = num * 10 + (ch - '0');
1239 index++;
1240 }
1241 else
1242 {
1243 // https://bugzilla.mozilla.org/show_bug.cgi?id=564839 -- stop parsing at null char
1244 if (ch == 0 && index < slength)
1245 {
1246 slength = index;
1247 }
1248 break;
1249 }
1250 }
1251 if (negexp) {
1252 num = -num;
1253 }
1254 exp10 += num;
1255 }
1256 // ecma3 compliant to allow leading and trailing whitespace
1257 index = skipSpaces(s, index);
1258
1259 // If we got no digits, check for Infinity/-Infinity, else fail
1260 if (numDigits == 0) {
1261 // check for the strings "Infinity" and "-Infinity"
1262 if (s->matchesLatin1("Infinity", 8, index)) {
1263 index += 8;
1264 // there may be trailing whitespace
1265 if (index < slength && skipSpaces(s, index) == index)
1266 return falsefalse;
1267 *value = (negate ? MathUtils::kNegInfinity : MathUtils::kInfinity);
1268 return truetrue;
1269 }
1270 return falsefalse;
1271 }
1272
1273 // If we got digits, but we're in strict mode and not at end of string, fail.
1274 if (index < slength && strict) {
1275 return falsefalse;
1276 }
1277
1278 // Pass 2: We've validated the input, now calculate the result.
1279
1280 // Read all the digits
1281 // cn: multiplyinging by negative powers of ten injects roundoff errors.
1282 /*
1283 while ((*s >= '0' && *s <= '9') || *s == '.') {
1284 if (*s != '.') {
1285 result += powerOfTen(exp10--, *s - '0');
1286 }
1287 s++;
1288 }
1289 */
1290 const int32_t limit = inStr->core()->currentBugCompatibility()->bugzilla513018 ? index : slength;
1291 index = start;
1292 if (numDigits > 15) // after 15 digits, we can run into roundoff error
1293 {
1294 BigInteger exactInt;
1295 exactInt.setFromInteger(0);
1296 int32_t decimalDigits= -1;
1297 while (index < limit) {
1298 ch = s[index];
1299 if ((ch >= '0' && ch <= '9') || ch == '.') {
1300 if (decimalDigits != -1) {
1301 decimalDigits++;
1302 }
1303 if (ch != '.') {
1304 exactInt.multAndIncrementBy(10, ch - '0');
1305 } else {
1306 decimalDigits=0;
1307 }
1308 index++;
1309 }
1310 else
1311 break;
1312 }
1313 if (decimalDigits > 0) {
1314 exp10 -= decimalDigits;
1315 }
1316 if (exp10 > 0) {
1317 exactInt.multBy(quickPowTen(exp10));
1318 exp10 = 0;
1319 }
1320 // This is an approximation which is at most off by 1. To gain complete accuracy, we need to implement
1321 // the algorithm commented out below
1322 result = exactInt.doubleValueOf();
1323 if (exp10 < 0) {
1324 if (exp10 < -307) { // max positive value is e308. min neg value is e-324 Avoid overflow
1325 int32_t diff = exp10 + 307;
1326 result /= quickPowTen(-diff);
1327 exp10 -= diff;
1328 }
1329 result /= quickPowTen(-exp10); // actually more accurate than multiplying by negative power of 10
1330 }
1331 }
1332 else // we can use double
1333 {
1334 int32_t decimalDigits= -1;
1335 while (index < limit)
1336 {
1337 ch = s[index];
1338 if ((ch >= '0' && ch <= '9') || ch == '.') {
1339 if (decimalDigits != -1) {
1340 decimalDigits++;
1341 }
1342 if (ch != '.') {
1343 result = result*10 + ch - '0';
1344 } else {
1345 decimalDigits=0;
1346 }
1347 index++;
1348 }
1349 else
1350 break;
1351 }
1352 if (decimalDigits > 0) {
1353 exp10 -= decimalDigits;
1354 }
1355 if (exp10 >= 0) {
1356 result *= quickPowTen(exp10);
1357
1358 } else {
1359 if (exp10 < -307) { // max positive value is e308. min neg value is e-324. Avoid overflow
1360 int32_t diff = exp10 + 307;
1361 result /= quickPowTen(-diff);
1362 exp10 -= diff;
1363 }
1364 result /= quickPowTen(-exp10); // actually more accurate than multiplying by negative power of 10
1365 }
1366 }
1367
1368 *value = negate ? -result : result;
1369 return truetrue;
1370 }
1371
1372 // cn: The following Alg is based on William Clinger's paper "How to Read Floating Point Numbers Accurately"
1373 //
1374 // Using BigIntegers gives us a result which is either the best approximation or the next best approximation.
1375 // The rest of this algorithm should be completed to adjust the result in order to make sure its the best.
1376 // Held off on finishing it as there are bigger fish to fry, the ecmaScript number tests pass with the above
1377 // and perhaps that's good enough. To finish
1378 // 1) BigInteger class needs to be extended to support negative numbers
1379 // 2) The section under findClosestFloat after the comment "Compare" needs to be fleshed out with the
1380 // details from Clinger's paper to figure out if the result, nextFloat(), or previousFloat() should
1381 // be returned. Note that no loop needs to be completed because the BigInteger
1382 // calclulation gaurantees we are at worst one estimate off the true value and the loop will
1383 // only require a single iteration.
1384 // 3) Implement nextFloat() / previousFloat()
1385 /*
1386 // Check if both the mantissa and exponent are guaranteed to accurately fit into a double
1387 // (mantissa is <= 2^53, exponent is <= 10^23) If so, we can get away with double math.
1388 if ( (exactInt->numWords == 1) ||
1389 (exactInt->numWords == 2 && wordBuffer[1] < 0x00200000) ) { // i.e. its < 0x0020 0000 0000 0000 == 2^53
1390
1391 if (exp10 >= 0 && exp10 < 23) { // log5-of-two^53 = ceil( log(pow(2,53)) / log(5) ) = ceil(22.82585758) = 23
1392 double f = exactInt->doubleValueOf();
1393 result = f * quickPowTen(exp10);
1394 }
1395 else if (exp10 < 0 && -exp10 < 23) {
1396 double f = exactInt->doubleValueOf();
1397 result = f / quickPowTen(exp10);
1398 }
1399 else { // exponent doesn't fit exactly into a double
1400 result = Bellerophon_MultiplyAndTest(exactInt, exp10, 0);
1401 }
1402 } else if (exactInt->numWords == 2) { // i.e. its < 0x0000 0001 0000 0000 0000 0000 == 2^64
1403 if ( (exp10 >= 0) && (exp10 < 23) {
1404 result = Bellerophon_MultiplyAndTest(exactInt, exp10, 0);
1405 } else { // exp10 < 0 || exp10 > 23)
1406 result = Bellerophon_MultiplyAndTest(exactInt, exp10, 3);
1407 }
1408 } else {
1409 if ( (exp10 >= 0) && (exp10 < 23) {
1410 result = Bellerophon_MultiplyAndTest(exactInt, exp10, 1);
1411 } else { // exp10 < 0 || exp10 > 23)
1412 result = Bellerophon_MultiplyAndTest(exactInt, exp10, 4);
1413 }
1414 }
1415
1416 *value = negate ? -result : result;
1417 return true;
1418 }
1419
1420 MathUtils::Bellerophon_MultiplyAndTest(BigInteger* f, int32_t exp10, int32_t slop)
1421 {
1422 double x = f->doubleValueOf();
1423 double y = quickPowTen(exp10);
1424 double estimate = x*y;
1425 uint64_t mantissaEstimate = frexp(estimate,&e);
1426 int64_t lowBits = (int64_t)(mantissaEstimate % 2048); // two^p-n = 2^(64-53) = 2^11 = 2048
1427
1428 // check if slop is large enough to make a difference when rounding to 53 bits
1429 int64_t diff = lowBits - 1024;
1430 if ( (diff >= 0 && diff <= slop) ||
1431 (diff < 0 && -diff <= slop) )
1432 {
1433 findClosestFloat(f,exp10,estimate);
1434 }
1435
1436 return estimate;
1437 }
1438
1439 MathUtils::findClosestFloat(BigInteger* f, int32_t e, double estimate)
1440 {
1441 int32_t k;
1442 uint64_t m = frexp(estimate,&k);
1443 BigInteger *m = BigInteger::newFromUint64(m);
1444 BigInteger* compX;
1445 BigInteger* compY;
1446
1447 if (e >= 0) {
1448 if (k >= 0) {
1449 compX = f->multBy(quickPow2(e));
1450 compY = m->multBy(pow(beta,k));
1451 } else {
1452 compX = f->multBy(quickPow2(e))->multBy( pow(beta,-k) );
1453 compY = m;
1454 }
1455 } else { // e < 0
1456 if (k >= 0) {
1457 compX = f;
1458 compY = m->multBy(pow(beta,k))->multBy(quickPow2(-e));
1459 } else {
1460 compX = f->multBy(pow(beta,-k));
1461 compY = m->multBy(quickPow2(-e)) );
1462 }
1463 }
1464 // compare
1465 BigInteger* diff = compX->subtract(compY); // need to modify BigInteger to support negative values
1466 // rest of compare would go here
1467 }
1468
1469 */
1470
1471 //
1472 // Random number generator
1473 //
1474
1475 /*-------------------------------------------------------------------------*/
1476 /* Coefficients used in the pure random number generator. */
1477#define c315731L 15731L
1478#define c2789221L 789221L
1479#define c11376312589L 1376312589L
1480
1481 /* Read-only, XOR masks for random # generation. */
1482 static const uint32_t Random_Xor_Masks[31] =
1483 {
1484 /* First mask is for generating sequences of 3 (2^2 - 1) random numbers,
1485 / Second mask is for generating sequences of 7 (2^3 - 1) random numbers,
1486 / etc. Numbers to the right are number of bits of random number sequence.
1487 / It generates 2^n - 1 numbers. */
1488 0x00000003L, 0x00000006L, 0x0000000CL, 0x00000014L, /* 2,3,4,5 */
1489 0x00000030L, 0x00000060L, 0x000000B8L, 0x00000110L, /* 6,7,8,9 */
1490 0x00000240L, 0x00000500L, 0x00000CA0L, 0x00001B00L, /* 10,11,12,13 */
1491 0x00003500L, 0x00006000L, 0x0000B400L, 0x00012000L, /* 14,15,16,17 */
1492 0x00020400L, 0x00072000L, 0x00090000L, 0x00140000L, /* 18,19,20,21 */
1493 0x00300000L, 0x00400000L, 0x00D80000L, 0x01200000L, /* 22,23,24,25 */
1494 0x03880000L, 0x07200000L, 0x09000000L, 0x14000000L, /* 26,27,28,29 */
1495 0x32800000L, 0x48000000L, 0xA3000000L /* 30,31,32 */
1496
1497 }; /* Randomizer Xor mask values CRK 10-29-90 */
1498 /*-------------------------------------------------------------------------*/
1499
1500 //TRandomFast gRandomFast = { 0, 0, 0}; Initialized in global.cpp
1501
1502 /*-*-------------------------------------------------------------------------
1503 / Function
1504 / RandomFastInit
1505 /
1506 / Purpose
1507 / This initializes the fast random number generator.
1508 /
1509 / Notes
1510 / The fast random number generator has some unique qualities:
1511 /
1512 / 1) The same sequence repeats every 2^n-1 generations.
1513 / 2) Each number can be generated extremely rapidly.
1514 /
1515 / Entry
1516 / pRandomFast = Pointer to data structure for current state of a fast
1517 / pseudorandom number generator.
1518 /--------------------------------------------*/
1519 void MathUtils::RandomFastInit(pTRandomFast pRandomFast)
1520 {
1521 int32_t n = 31; // Changed from 32 to 31 per Prince (you mean "The Artist"?)
1522
1523 pRandomFast->uValue = (uint32_t)(VMPI_getPerformanceCounter());
1524
1525 /* Figure out the sequence length (2^n - 1). */
1526 pRandomFast->uSequenceLength = (1L << n) - 1L;
1527
1528 /* Gather the XOR mask value. */
1529 pRandomFast->uXorMask = Random_Xor_Masks[n - 2];
1530 }
1531 /*-------------------------------------------------------------------------*/
1532
1533 /*-*-------------------------------------------------------------------------
1534 / Function
1535 / imRandomFastNext
1536 /
1537 / Purpose
1538 / This generates a new pseudorandom number using the fast random number
1539 / generator.
1540 /
1541 / Notes
1542 / The fast random number generator must be initialized before use.
1543 /
1544 / This is implemented as a macro for the best possible performance.
1545 /
1546 / Entry
1547 / pRandomFast = Pointer to data structure for current state of a fast
1548 / pseudorandom number generator.
1549 /
1550 / Exit
1551 / Returns the next pseudorandom number in the sequence.
1552 /--------------------------------------------*/
1553#define RandomFastNext(_pRandomFast)( ((_pRandomFast)->uValue & 1L) ? ((_pRandomFast)->
uValue = ((_pRandomFast)->uValue >> 1) ^ (_pRandomFast
)->uXorMask) : ((_pRandomFast)->uValue = ((_pRandomFast
)->uValue >> 1)) ) \
1554( \
1555 ((_pRandomFast)->uValue & 1L) \
1556 ? ((_pRandomFast)->uValue = ((_pRandomFast)->uValue >> 1) ^ (_pRandomFast)->uXorMask) \
1557 : ((_pRandomFast)->uValue = ((_pRandomFast)->uValue >> 1)) \
1558)
1559 /*-------------------------------------------------------------------------*/
1560
1561 /*-*-------------------------------------------------------------------------
1562 / Function
1563 / RandomPureHasher
1564 /
1565 / Purpose
1566 / This generates a pseudorandom number from a given seed using the high
1567 / quality random number generator.
1568 /
1569 / Notes
1570 / The caller must pass in a seed value. This value is typically the
1571 / output of the fast random number generator multiplied by a prime
1572 / number, but can be any seed which was not derived from the output of
1573 / this function.
1574 /
1575 / Please note that this random number generator is really as hasher.
1576 / It will not work well if you pass its own output in as the next input.
1577 /
1578 / A given input seed always generates the same pseudorandom output result.
1579 /
1580 / The feature of this hasher is that the value returned from one seed
1581 / to the next is highly uncorrelated to the seed value. High quality
1582 / multidimensional random numbers can be generated by using a sum of
1583 / prime multiples of the seed values.
1584 /
1585 / For example, to generate a vector given three seed values sx, sy and sz,
1586 / use something like the following:
1587 / result = imRandomPureHasher(59*sx + 67*sy + 71*sz);
1588 /
1589 / Entry
1590 / iSeed = Starting seed value.
1591 /
1592 / Exit
1593 / Returns the next pseudorandom value in the sequence.
1594 /--------------------------------------------*/
1595 int32_t MathUtils::RandomPureHasher(int32_t iSeed)
1596 {
1597 int32_t iResult;
1598
1599 /* Adapted from "A Recursive Implementation of the Perlin Noise Function,"
1600 / by Greg Ward, Graphics Gems II, p. 396. */
1601
1602 /* First, hash s as a preventive measure. This helps ensure the first
1603 / few numbers are more random when used in conjunction with the fast
1604 / random number generator, which starts off with the first 10 or so
1605 / numbers all zero in the lowe bits. */
1606 iSeed = ((iSeed << 13) ^ iSeed) - (iSeed >> 21);
1607
1608 /* Next, use a third order odd polynomial, better than linear. */
1609 iResult = (iSeed*(iSeed*iSeed*c315731L + c2789221L) + c11376312589L) & kRandomPureMax0x7FFFFFFFL;
1610
1611 /* DJ14dec94 -- The above wonderful expression always returns odd
1612 / numbers, and this is easy to prove. So we add the seed back to
1613 / the result which again randomizes bit zero. */
1614 iResult += iSeed;
1615
1616 /* DJ17nov95 -- The above always returns values that are NEVER divisible
1617 / evenly by four, so do additional hashing. */
1618 iResult = ((iResult << 13) ^ iResult) - (iResult >> 21);
1619
1620 return iResult;
1621 }
1622 /* ------------------------------------------------------------------------------ */
1623 int32_t MathUtils::GenerateRandomNumber(pTRandomFast pRandomFast)
1624 {
1625 // Fill out gRandomFast if it is uninitialized.
1626 if (pRandomFast->uValue == 0) {
1627 RandomFastInit(pRandomFast);
1628 }
1629
1630 long aNum = RandomFastNext(pRandomFast)( ((pRandomFast)->uValue & 1L) ? ((pRandomFast)->uValue
= ((pRandomFast)->uValue >> 1) ^ (pRandomFast)->
uXorMask) : ((pRandomFast)->uValue = ((pRandomFast)->uValue
>> 1)) );
1631
1632 /* Use the pure random number generator to hash the result. */
1633 aNum = RandomPureHasher(aNum * 71L);
1634
1635 return aNum & kRandomPureMax0x7FFFFFFFL;
1636 }
1637
1638 int32_t MathUtils::Random(int32_t range, pTRandomFast pRandomFast)
1639 {
1640 if (range > 0) {
1641 return GenerateRandomNumber(pRandomFast) % range;
1642 } else {
1643 return 0;
1644 }
1645 }
1646
1647 void MathUtils::initRandom(TRandomFast *seed)
1648 {
1649 RandomFastInit(seed);
1650 }
1651
1652 double MathUtils::random(TRandomFast *seed)
1653 {
1654 return (double)GenerateRandomNumber(seed) / ((double)kRandomPureMax0x7FFFFFFFL+1.0);
1655 }
1656
1657 /*-------------------------------------------------------------------------*/
1658 /* Begining D2A class implemention, used for converting double values to */
1659 /* strings accurately. */
1660
1661 // This algorithm for printing floating point numbers is based on the paper
1662 // "Printing floating-point numbers quickly and accurately" by Robert G Burger
1663 // and R. Kent Dybvig http://www.cs.indiana.edu/~burger/FP-Printing-PLDI96.pdf
1664 // See that paper for details.
1665 //
1666 // The algorithm generates the shortest, correctly rounded output string that
1667 // converts to the same number when read back in. This implementation assumes
1668 // IEEE unbiased rounding, that the input is a 64 bit (base 2) floating point
1669 // number composed of a sign bit, a hidden bit (valued 1) and 52 explict bits of
1670 // mantissa data, and an 11 bit exponent (in base 2). It also assumes the output
1671 // desired is base 10.
1672
1673
1674 // ------------------------------------------------------
1675 // Support functions
1676
1677 // For results larger than 10^22, error creeps into the double estimate. Use BigIntegers to avoid error
1678 static inline void quickBigPowTen(int32_t exp, BigInteger* result)
1679 {
1680 if (exp < 22 && exp > 0) {
1681 result->setFromDouble(kPowersOfTen[exp]);
1682 } else if (exp > 0) {// double starts to include roundoff error after 1e22, need to compute exactly
1683 result->setFromDouble(kPowersOfTen[21]);
1684 exp -= 21;
1685 while (exp-- > 0)
1686 result->multBy((int32_t)10);
1687 } else { // we won't get here because we only deal in positive exponents
1688 result->setFromDouble(MathUtils::pow(10,exp)); // but just in case
1689 }
1690 }
1691
1692 // Use left shifts to compute powers of 2 < 63
1693 static inline double quickPowTwo(int32_t exp)
1694 {
1695 static uint64_t one = 1; // max we can shift is 64bits
1696 if (exp < 64 && exp > 0)
1697 return (double)(one << exp);
1698 else
1699 return MathUtils::pow(2,exp);
1700 }
1701
1702 // nextDigit() returns the next relevant digit in the number being converted, else -1 if
1703 // there are no more relevant digits.
1704 int32_t D2A::nextDigit()
1705 {
1706 if (finished)
1707 return -1;
1708
1709 boolbool withinLowEndRoundRange;
1710 boolbool withinHighEndRoundRange;
1711 int32_t quotient;
1712
1713 if ( bFastEstimateOk )
1714 {
1715 quotient = (int32_t)(dr / ds);
1716 double mod = MathUtils::mod(dr,ds);
1717 dr = mod;
1718
1719 // check if remaining ratio r/s is within the range of floats which would round to the value we have output
1720 // so far when read in from a string.
1721 withinLowEndRoundRange = (lowOk ? (dr <= dMMinus) : (dr < dMMinus));
1722 withinHighEndRoundRange = (highOk ? (dr+dMPlus >= ds) : (dr+dMPlus > ds));
1723 }
1724 else
1725 {
1726 BigInteger bigQuotient;
1727 bigQuotient.setFromInteger(0);
1728 r.divBy(&s, &bigQuotient); // r = r %s, bigQuotient = r / s.
1729 quotient = (int32_t)(bigQuotient.wordBuffer[0]); // todo: optimize away need for BigInteger result? We know it should always be a single digit
1730 // r <= mMinus : r < rMinus
1731 withinLowEndRoundRange = (lowOk ? (r.compare(&mMinus) != 1) : (r.compare(&mMinus) == -1));
1732 // r+mPlus >= s : r+mPlus > s
1733 withinHighEndRoundRange = (highOk ? (r.compareOffset(&s,&mPlus) != -1) : (r.compareOffset(&s,&mPlus) == 1));
1734 }
1735
1736
1737 if (quotient < 0 || quotient > 9)
1738 {
1739 AvmAssert(quotient >= 0)do { } while (0);
1740 AvmAssert(quotient < 10)do { } while (0);
1741 quotient = 0;
1742 }
1743 if (!withinLowEndRoundRange)
1744 {
1745 if (!withinHighEndRoundRange) // if not within either error range, set up to generate the next digit.
1746 {
1747 if ( bFastEstimateOk )
1748 {
1749 dr *= 10;
1750 dMPlus *= 10;
1751 dMMinus *= 10;
1752 }
1753 else
1754 {
1755 r.multBy((int32_t)10);
1756 mPlus.multBy((int32_t)10);
1757 mMinus.multBy((int32_t)10);
1758 }
1759 }
1760 else
1761 {
1762 quotient++;
1763 finished = truetrue;
1764 }
1765 }
1766 else if (!withinHighEndRoundRange)
1767 {
1768 finished = truetrue;
1769 }
1770 else
1771 {
1772 if ( (bFastEstimateOk && (dr*2 < ds)) ||
1773 (!bFastEstimateOk && (r.compareOffset(&s,&r) == -1)) ) // if (r*2 < s) todo: faster to do lshift and compare?
1774 {
1775 finished = truetrue;
1776 }
1777 else
1778 {
1779 quotient++;
1780 finished = truetrue;
1781 }
1782 }
1783 return quotient;
1784 }
1785
1786
1787 int32_t D2A::fixup_ExponentEstimate(int32_t exponentEstimate)
1788 {
1789 int32_t correctedEstimate;
1790 if (bFastEstimateOk)
1791 {
1792 if (highOk ? (dr+dMPlus) >= ds : dr+dMPlus > ds)
1793 {
1794 correctedEstimate = exponentEstimate+1;
1795 }
1796 else
1797 {
1798 dr *= 10;
1799 dMPlus *= 10;
1800 dMMinus *= 10;
1801 correctedEstimate = exponentEstimate;
1802 }
1803 }
1804 else
1805 {
1806 if (highOk ? (r.compareOffset(&s,&mPlus) != -1) : (r.compareOffset(&s,&mPlus) == 1))
1807 {
1808 correctedEstimate = exponentEstimate+1;
1809 }
1810 else
1811 {
1812 r.multBy((int32_t)10);
1813 mPlus.multBy((int32_t)10);
1814 mMinus.multBy((int32_t)10);
1815 correctedEstimate = exponentEstimate;
1816 }
1817 }
1818 return correctedEstimate;
1819 }
1820
1821 int32_t D2A::scale()
1822 {
1823 // estimate base10 exponent:
1824 int32_t base2Exponent = e + mantissaPrec-1;
1825 int32_t exponentEstimate = (int32_t)MathUtils::ceil( (base2Exponent * kLog2_10) - 0.0000000001);
1826
1827
1828 if (bFastEstimateOk)
1829 {
1830 double scale = quickPowTen( (exponentEstimate > 0) ? exponentEstimate : -exponentEstimate);
1831
1832 if (exponentEstimate >= 0)
1833 {
1834 ds *= scale;
1835 return fixup_ExponentEstimate(exponentEstimate);
1836 }
1837 else
1838 {
1839 dr *= scale;
1840 dMPlus *= scale;
1841 dMMinus *= scale;
1842 return fixup_ExponentEstimate(exponentEstimate);
1843 }
1844 }
1845 else
1846 {
1847 BigInteger scale;
1848 scale.setFromInteger(0);
1849 quickBigPowTen( (exponentEstimate > 0) ? exponentEstimate : -exponentEstimate, &scale);
1850
1851 if (exponentEstimate >= 0)
1852 {
1853 s.multBy(&scale);
1854 return fixup_ExponentEstimate(exponentEstimate);
1855 }
1856 else
1857 {
1858 r.multBy(&scale);
1859 mPlus.multBy(&scale);
1860 mMinus.multBy(&scale);
1861
1862 return fixup_ExponentEstimate(exponentEstimate);
1863 }
1864 }
1865 }
1866
1867 D2A::~D2A()
1868 {
1869 }
1870
1871
1872 D2A::D2A(double avalue, int32_t mode, int32_t minDigits)
1873 : value(avalue), finished(falsefalse), bFastEstimateOk(falsefalse), minPrecision(minDigits)
1874 {
1875 // break double down into integer mantissa (max size unsigned 53 bits) and integer base 2
1876 // expondent (max size 11 signed bits).
1877 mantissa = MathUtils::frexp(value,&e); // value = mantissa*2^e, mantissa and e are integers.
1878
1879 if (mode != MathUtils::DTOSTR_NORMAL)
1880 {
1881 lowOk = highOk = truetrue;
1882 }
1883 else
1884 {
1885 boolbool round = !(mantissa & 0x1); // if mantissa is even, need to round.
1886 lowOk = highOk = round; // IEEE standard unbiased rounding.
1887 }
1888
1889 // calculate #of bits of precision needed to encode this mantissa (max is 53 bits)
1890 long base2Precision = 53;
1891 mantissaPrec = base2Precision; // double has 53 bits of precision in the mantissa, but first bit is assumed
1892 while ( mantissaPrec != 0 && ((mantissa >> --mantissaPrec & 1) == 0) )
1893 ; // skip leading zeros
1894 mantissaPrec++;
1895
1896 int32_t absE = (e > 0 ? e : -e);
1897 if ((absE + mantissaPrec-1) < 50) // we get error prone when there is more than 15 base 10 digits of precision. (15 / kLog2_10) == 49.828921423310435218054791442341)
1898 bFastEstimateOk = truetrue;
1899
1900 // Represent this double as a ratio of two integers. Use infinitely precise integers if bFastEstimate is false.
1901 // mPlus and mMinus represent the range of doubles around this value which would
1902 // round to this same value when converted from string back to number via atod().
1903 if (bFastEstimateOk)
1904 {
1905 if (e >= 0)
1906 {
1907 if (mantissa != 0x10000000000000LL) // == 4503599627370496 == pow(2, base2Precision-1)) == 2^52 (i.e. just under becoming +1 to e and rolling over to 0)
1908 {
1909 double be = quickPowTwo(e); // 2^e
1910
1911 dr = (double)mantissa*be*2;
1912 ds = 2;
1913 dMPlus = be;
1914 dMMinus = be;
1915 }
1916 else
1917 {
1918 double be = quickPowTwo(e);; // 2^e
1919 double be1 = be*2; // 2^(e+1)
1920
1921 dr = (double)mantissa*be1*2;
1922 ds = 4; // 2*2;
1923 dMPlus = be1;
1924 dMMinus = be;
1925 }
1926 }
1927 else if ( mantissa != MathUtils::pow(2, base2Precision-1) )
1928 {
1929 dr = (double)mantissa*2.0;
1930 ds = quickPowTwo(1-e); // pow(2,-e)*2;
1931 dMPlus = 1;
1932 dMMinus = 1;
1933 }
1934 else
1935 {
1936 dr = (double)mantissa*4.0;
1937 ds = quickPowTwo(2-e); // pow(2, 1-e)*2;
1938 dMPlus = 2;
1939 dMMinus = 1;
1940 }
1941
1942 // adjust stopping conditions if a particular number of digits are requested
1943 if (mode != MathUtils::DTOSTR_NORMAL)
1944 {
1945 double fixedPrecisionPowTen = quickPowTen( minPrecision );
1946 ds *= fixedPrecisionPowTen;
1947 dr *= fixedPrecisionPowTen;
1948 }
1949 }
1950 else
1951 {
1952 if (e >= 0)
1953 {
1954 if (mantissa != 0x10000000000000LL) // == 4503599627370496 == pow(2, base2Precision-1)) == 2^53 (i.e. just under becoming +1 to e and rolling over to 0)
1955 {
1956 BigInteger be;
1957 be.setFromInteger(1);
1958 be.lshiftBy(e); // == pow(2,e);
1959
1960 r.setFromDouble(value);
1961 r.lshiftBy(1); // r = mantissa*be*2
1962
1963 s.setFromInteger(2);
1964
1965 mPlus.setFromBigInteger(&be,0,be.numWords); // == be;
1966 mMinus.setFromBigInteger(&be,0,be.numWords); // == be;
1967 }
1968 else
1969 {
1970 BigInteger be;
1971 be.setFromInteger(1);
1972 be.lshiftBy(e); // be = pow(2,e);
1973
1974 BigInteger be1;
1975 be1.setFromInteger(0);
1976 be.lshift(1,&be1); // be1 == be*2;
1977
1978 r.setFromDouble(value*4.0); // r == mantissa*be1*2;
1979 s.setFromInteger(4); // 2*2
1980 mPlus.setFromBigInteger(&be1,0,be1.numWords); // == be1;
1981 mMinus.setFromBigInteger(&be,0,be.numWords); // == be;
1982 }
1983 }
1984 else if ( mantissa != MathUtils::pow(2, base2Precision-1) )
1985 {
1986 r.setFromDouble( (double)(mantissa*2) );
1987 s.setFromInteger(2);
1988 s.lshiftBy(-e); // s = pow(2,-e)*2;
1989 mPlus.setFromInteger(1);
1990 mMinus.setFromInteger(1);
1991 }
1992 else
1993 {
1994 r.setFromDouble( (double)(mantissa*4) );
1995 s.setFromInteger(2);
1996 s.lshiftBy(1-e); // s == pow(2, 1-e)*2;
1997 mPlus.setFromInteger(2);
1998 mMinus.setFromInteger(1);
1999 }
2000
2001 // adjust stopping conditions if a particular number of digits are requested
2002 // Note: this is a cheap approximation of the correct adjustment. It will likely
2003 // lead to occasional off by one errors in the final digit generated for toPrecision() when
2004 // a truncation might be required. This mode is actually less accurate than normal mode anyway.
2005 // Only normal mode is gauranteed to generate a string which will result in the same value when converted back
2006 // from string (base 10) to double (base 2). Extra digits generated by toPrecision will often be junk
2007 // (but will display a value more similar to what you would see in a code debugger, where 0.012 can
2008 // appear as 0.012999999999999999999 or such apparent nonsense)
2009 if (mode != MathUtils::DTOSTR_NORMAL)
2010 {
2011 BigInteger bFixedPrecisionPowTen;
2012 bFixedPrecisionPowTen.setFromInteger(0);
2013 quickBigPowTen( minPrecision, &bFixedPrecisionPowTen );
2014 s.multBy(&bFixedPrecisionPowTen);
2015 r.multBy(&bFixedPrecisionPowTen);
2016 }
2017 }
2018 base10Exp = scale();
2019 }
2020
2021
2022 /*
2023 A2D::A2D(String source,int32_t radix)
2024 {
2025 bitsPerChar = 0;
2026 while ( (radix >>= 1) != 1 )
2027 bitsPerChar++;
2028
2029 source = source;
2030 radix = radix;
2031 next = 0;
2032 end = string.length();
2033 finished = false;
2034
2035 // skip leading 0's
2036 while( nextDigit() == 0);
2037 next--;
2038 }
2039
2040 inline int32_t A2D::parseIntDigit(wchar ch)
2041 {
2042 if (ch >= '0' && ch <= '9') {
2043 return (ch - '0');
2044 } else if (ch >= 'a' && ch <= 'z') {
2045 return (ch - 'a' + 10);
2046 } else if (ch >= 'A' && ch <= 'Z') {
2047 return (ch - 'A' + 10);
2048 } else {
2049 return -1;
2050 }
2051 }
2052 int32_t A2D::nextDigit()
2053 {
2054 if (finished)
2055 return -1;
2056
2057 wc = source[next++];
2058 int32_t result = parseIntDigit(wc);
2059 if (next == end || result == -1)
2060 finished = true;
2061 return result;
2062 }
2063 */
2064}