74 #ifdef WORDS_BIGENDIAN
79 #ifndef HAVE_LONG_LONG
82 #define Omit_Private_Memory
199 typedef unsigned Long ULong;
204 #define Bug(x) {fprintf(stderr, "%s\n", x); exit(1);}
217 #undef Avoid_Underflow
226 #ifndef NO_INFNAN_CHECK
240 #define DBL_MAX_10_EXP 308
241 #define DBL_MAX_EXP 1024
247 #define DBL_MAX_10_EXP 75
248 #define DBL_MAX_EXP 63
250 #define DBL_MAX 7.2370055773322621e+75
255 #define DBL_MAX_10_EXP 38
256 #define DBL_MAX_EXP 127
258 #define DBL_MAX 1.7014118346046923e+38
262 #define LONG_MAX 2147483647
286 #if defined(IEEE_8087) + defined(IEEE_MC68k) + defined(VAX) + defined(IBM) != 1
287 Exactly one of IEEE_8087, IEEE_MC68k, VAX, or IBM should be defined.
291 typedef union {
double d; ULong L[2]; }
U;
296 #define word0(x) ((ULong *)&x)[1]
297 #define word1(x) ((ULong *)&x)[0]
299 #define word0(x) ((ULong *)&x)[0]
300 #define word1(x) ((ULong *)&x)[1]
304 #define word0(x) ((U*)&x)->L[1]
305 #define word1(x) ((U*)&x)->L[0]
307 #define word0(x) ((U*)&x)->L[0]
308 #define word1(x) ((U*)&x)->L[1]
310 #define dval(x) ((U*)&x)->d
317 #if defined(IEEE_8087) + defined(VAX)
318 #define Storeinc(a,b,c) (((unsigned short *)a)[1] = (unsigned short)b, \
319 ((unsigned short *)a)[0] = (unsigned short)c, a++)
321 #define Storeinc(a,b,c) (((unsigned short *)a)[0] = (unsigned short)b, \
322 ((unsigned short *)a)[1] = (unsigned short)c, a++)
333 #define Exp_shift1 20
334 #define Exp_msk1 0x100000
335 #define Exp_msk11 0x100000
336 #define Exp_mask 0x7ff00000
340 #define Exp_1 0x3ff00000
341 #define Exp_11 0x3ff00000
343 #define Frac_mask 0xfffff
344 #define Frac_mask1 0xfffff
347 #define Bndry_mask 0xfffff
348 #define Bndry_mask1 0xfffff
350 #define Sign_bit 0x80000000
356 #ifndef NO_IEEE_Scale
357 #define Avoid_Underflow
359 #undef Sudden_Underflow
365 #define Flt_Rounds FLT_ROUNDS
371 #ifdef Honor_FLT_ROUNDS
372 #define Rounding rounding
373 #undef Check_FLT_ROUNDS
374 #define Check_FLT_ROUNDS
376 #define Rounding Flt_Rounds
380 #undef Check_FLT_ROUNDS
381 #undef Honor_FLT_ROUNDS
383 #undef Sudden_Underflow
384 #define Sudden_Underflow
389 #define Exp_shift1 24
390 #define Exp_msk1 0x1000000
391 #define Exp_msk11 0x1000000
392 #define Exp_mask 0x7f000000
395 #define Exp_1 0x41000000
396 #define Exp_11 0x41000000
398 #define Frac_mask 0xffffff
399 #define Frac_mask1 0xffffff
402 #define Bndry_mask 0xefffff
403 #define Bndry_mask1 0xffffff
405 #define Sign_bit 0x80000000
407 #define Tiny0 0x100000
416 #define Exp_msk1 0x80
417 #define Exp_msk11 0x800000
418 #define Exp_mask 0x7f80
421 #define Exp_1 0x40800000
422 #define Exp_11 0x4080
424 #define Frac_mask 0x7fffff
425 #define Frac_mask1 0xffff007f
428 #define Bndry_mask 0xffff007f
429 #define Bndry_mask1 0xffff007f
431 #define Sign_bit 0x8000
445 #define rounded_product(a,b) a = rnd_prod(a, b)
446 #define rounded_quotient(a,b) a = rnd_quot(a, b)
448 extern double rnd_prod(), rnd_quot();
450 extern double rnd_prod(
double,
double), rnd_quot(
double,
double);
453 #define rounded_product(a,b) a *= b
454 #define rounded_quotient(a,b) a /= b
457 #define Big0 (Frac_mask1 | Exp_msk1*(DBL_MAX_EXP+Bias-1))
458 #define Big1 0xffffffff
465 #define FFFFFFFF ((((unsigned long)0xffff)<<16)|(unsigned long)0xffff)
467 #define FFFFFFFF 0xffffffffUL
482 #define Llong long long
485 #define ULLong unsigned Llong
489 #ifndef MULTIPLE_THREADS
490 #define ACQUIRE_DTOA_LOCK(n)
491 #define FREE_DTOA_LOCK(n)
497 extern "C" double sb_strtod(
const char *s00,
char **se);
503 int k, maxwds, sign, wds;
522 len = (
sizeof(
Bigint) + (x-1)*
sizeof(ULong) +
sizeof(
double) - 1)
527 rv->sign = rv->wds = 0;
542 #define Bcopy(x,y) memcpy((char *)&x->sign, (char *)&y->sign, \
543 y->wds*sizeof(Long) + 2*sizeof(int))
548 (b, m, a)
Bigint *b;
int m, a;
571 y = *x * (ULLong)m + carry;
577 y = (xi & 0xffff) * m + carry;
578 z = (xi >> 16) * m + (y >> 16);
580 *x++ = (z << 16) + (y & 0xffff);
590 if (wds >= b->maxwds) {
605 (s, nd0, nd, y9) CONST
char *s;
int nd0, nd; ULong y9;
607 (CONST
char *s,
int nd0,
int nd, ULong y9)
615 for(k = 0, y = 1; x > y; y <<= 1, k++) ;
622 b->x[0] = y9 & 0xffff;
623 b->wds = (b->x[1] = y9 >> 16) ? 2 : 1;
629 do b = multadd(b, 10, *s++ -
'0');
636 b = multadd(b, 10, *s++ -
'0');
643 (x)
register ULong x;
650 if (!(x & 0xffff0000)) {
654 if (!(x & 0xff000000)) {
658 if (!(x & 0xf0000000)) {
662 if (!(x & 0xc0000000)) {
666 if (!(x & 0x80000000)) {
668 if (!(x & 0x40000000))
683 register ULong x = *y;
748 ULong *x, *xa, *xae, *xb, *xbe, *xc, *xc0;
759 if (a->wds < b->wds) {
771 for(x = c->x, xa = x + wc; x < xa; x++)
779 for(; xb < xbe; xc0++) {
785 z = *x++ * (ULLong)y + *xc + carry;
787 *xc++ = z & FFFFFFFF;
795 for(; xb < xbe; xb++, xc0++) {
796 if (y = *xb & 0xffff) {
801 z = (*x & 0xffff) * y + (*xc & 0xffff) + carry;
803 z2 = (*x++ >> 16) * y + (*xc >> 16) + carry;
816 z = (*x & 0xffff) * y + (*xc >> 16) + carry;
819 z2 = (*x++ >> 16) * y + (*xc & 0xffff) + carry;
827 for(; xb < xbe; xc0++) {
833 z = *x++ * y + *xc + carry;
843 for(xc0 = c->x, xc = xc0 + wc; wc > 0 && !*--xc; --wc) ;
858 static int CONST p05[3] = { 5, 25, 125 };
861 b = multadd(b, p05[i-1], 0);
893 ULong *x, *x1, *xe, z;
902 for(i = b->maxwds; n1 > i; i <<= 1)
906 for(i = 0; i < n; i++)
927 *x1++ = *x << k & 0xffff | z;
951 ULong *xa, *xa0, *xb, *xb0;
957 if (i > 1 && !a->x[i-1])
958 Bug(
"cmp called with a->x[a->wds-1] == 0");
959 if (j > 1 && !b->x[j-1])
960 Bug(
"cmp called with b->x[b->wds-1] == 0");
970 return *xa < *xb ? -1 : 1;
987 ULong *xa, *xae, *xb, *xbe, *xc;
1024 y = (ULLong)*xa++ - *xb++ - borrow;
1025 borrow = y >> 32 & (ULong)1;
1026 *xc++ = y & FFFFFFFF;
1031 borrow = y >> 32 & (ULong)1;
1032 *xc++ = y & FFFFFFFF;
1037 y = (*xa & 0xffff) - (*xb & 0xffff) - borrow;
1038 borrow = (y & 0x10000) >> 16;
1039 z = (*xa++ >> 16) - (*xb++ >> 16) - borrow;
1040 borrow = (z & 0x10000) >> 16;
1045 y = (*xa & 0xffff) - borrow;
1046 borrow = (y & 0x10000) >> 16;
1047 z = (*xa++ >> 16) - borrow;
1048 borrow = (z & 0x10000) >> 16;
1053 y = *xa++ - *xb++ - borrow;
1054 borrow = (y & 0x10000) >> 16;
1060 borrow = (y & 0x10000) >> 16;
1082 L = (word0(x) & Exp_mask) - (P-1)*Exp_msk1;
1083 #ifndef Avoid_Underflow
1084 #ifndef Sudden_Underflow
1093 #ifndef Avoid_Underflow
1094 #ifndef Sudden_Underflow
1097 L = -L >> Exp_shift;
1098 if (L < Exp_shift) {
1099 word0(a) = 0x80000 >> L;
1105 word1(a) = L >= 31 ? 1 : 1 << 31 - L;
1116 (a, e)
Bigint *a;
int *e;
1121 ULong *xa, *xa0, w, y, z;
1135 if (!y) Bug(
"zero y in b2d");
1141 d0 = Exp_1 | y >> (Ebits - k);
1142 w = xa > xa0 ? *--xa : 0;
1143 d1 = y << ((32-Ebits) + k) | w >> (Ebits - k);
1146 z = xa > xa0 ? *--xa : 0;
1148 d0 = Exp_1 | y << k | z >> (32 - k);
1149 y = xa > xa0 ? *--xa : 0;
1150 d1 = z << k | y >> (32 - k);
1157 if (k < Ebits + 16) {
1158 z = xa > xa0 ? *--xa : 0;
1159 d0 = Exp_1 | y << k - Ebits | z >> Ebits + 16 - k;
1160 w = xa > xa0 ? *--xa : 0;
1161 y = xa > xa0 ? *--xa : 0;
1162 d1 = z << k + 16 - Ebits | w << k - Ebits | y >> 16 + Ebits - k;
1165 z = xa > xa0 ? *--xa : 0;
1166 w = xa > xa0 ? *--xa : 0;
1168 d0 = Exp_1 | y << k + 16 | z << k | w >> 16 - k;
1169 y = xa > xa0 ? *--xa : 0;
1170 d1 = w << k + 16 | y << k;
1174 word0(d) = d0 >> 16 | d0 << 16;
1175 word1(d) = d1 >> 16 | d1 << 16;
1186 (d, e, bits)
double d;
int *e, *bits;
1188 (
double _d,
int *e,
int *bits)
1195 #ifndef Sudden_Underflow
1200 d0 = word0(d) >> 16 | word0(d) << 16;
1201 d1 = word1(d) >> 16 | word1(d) << 16;
1217 #ifdef Sudden_Underflow
1218 de = (int)(d0 >> Exp_shift);
1223 if ((de = (
int)(d0 >> Exp_shift)))
1228 if ((k = lo0bits(&y))) {
1229 x[0] = y | z << (32 - k);
1234 #ifndef Sudden_Underflow
1237 b->wds = (x[1] = z) ? 2 : 1;
1242 Bug(
"Zero passed to d2b");
1246 #ifndef Sudden_Underflow
1254 if (k = lo0bits(&y))
1256 x[0] = y | z << 32 - k & 0xffff;
1257 x[1] = z >> k - 16 & 0xffff;
1263 x[1] = y >> 16 | z << 16 - k & 0xffff;
1264 x[2] = z >> k & 0xffff;
1279 Bug(
"Zero passed to d2b");
1297 #ifndef Sudden_Underflow
1301 *e = (de - Bias - (P-1) << 2) + k;
1302 *bits = 4*P + 8 - k - hi0bits(word0(d) & Frac_mask);
1304 *e = de - Bias - (P-1) + k;
1307 #ifndef Sudden_Underflow
1310 *e = de - Bias - (P-1) + 1 + k;
1312 *bits = 32*i - hi0bits(x[i-1]);
1314 *bits = (i+2)*16 - hi0bits(x[i]);
1334 dval(da) = b2d(a, &ka);
1335 dval(db) = b2d(b, &kb);
1337 k = ka - kb + 32*(a->wds - b->wds);
1339 k = ka - kb + 16*(a->wds - b->wds);
1343 word0(da) += (k >> 2)*Exp_msk1;
1349 word0(db) += (k >> 2)*Exp_msk1;
1355 word0(da) += k*Exp_msk1;
1358 word0(db) += k*Exp_msk1;
1361 return dval(da) / dval(db);
1366 1e0, 1e1, 1e2, 1e3, 1e4, 1e5, 1e6, 1e7, 1e8, 1e9,
1367 1e10, 1e11, 1e12, 1e13, 1e14, 1e15, 1e16, 1e17, 1e18, 1e19,
1376 bigtens[] = { 1e16, 1e32, 1e64, 1e128, 1e256 };
1377 static CONST
double tinytens[] = { 1e-16, 1e-32, 1e-64, 1e-128,
1378 #ifdef Avoid_Underflow
1379 9007199254740992.*9007199254740992.e-256
1387 #define Scale_Bit 0x10
1391 bigtens[] = { 1e16, 1e32, 1e64 };
1392 static CONST
double tinytens[] = { 1e-16, 1e-32, 1e-64 };
1395 bigtens[] = { 1e16, 1e32 };
1396 static CONST
double tinytens[] = { 1e-16, 1e-32 };
1404 #define NAN_WORD0 0x7ff80000
1414 (sp, t)
char **sp, *t;
1416 (CONST
char **sp,
char *t)
1420 CONST
char *s = *sp;
1423 if ((c = *++s) >=
'A' && c <=
'Z')
1436 (rvp, sp)
double *rvp; CONST
char **sp;
1438 (
U *rvp, CONST
char **sp)
1443 int havedig, udx0, xshift;
1446 havedig = xshift = 0;
1450 while((c = *(CONST
unsigned char*)(s+1)) && c <=
' ')
1452 if (s[1] ==
'0' && (s[2] ==
'x' || s[2] ==
'X'))
1454 while((c = *(CONST
unsigned char*)++s)) {
1455 if (c >=
'0' && c <=
'9')
1457 else if (c >=
'a' && c <=
'f')
1459 else if (c >=
'A' && c <=
'F')
1461 else if (c <=
' ') {
1462 if (udx0 && havedig) {
1468 #ifdef GDTOA_NON_PEDANTIC_NANCHECK
1469 else if ( c ==
')' && havedig) {
1482 }
while((c = *++s));
1493 x[0] = (x[0] << 4) | (x[1] >> 28);
1494 x[1] = (x[1] << 4) | c;
1496 if ((x[0] &= 0xfffff) || x[1]) {
1497 word0(*rvp) = Exp_mask | x[0];
1507 (s00, se) CONST
char *s00;
char **se;
1509 (CONST
char *s00,
char **se)
1512 #ifdef Avoid_Underflow
1515 int bb2, bb5, bbe, bd2, bd5, bbbits, bs2, c, dsign,
1516 e, e1, esign, i, j, k, nd, nd0, nf, nz, nz0, sign;
1517 CONST
char *s, *s0, *s1;
1522 Bigint *bb, *bb1, *bd, *bd0, *bs, *delta;
1524 int inexact, oldinexact;
1526 #ifdef Honor_FLT_ROUNDS
1533 sign = nz0 = nz = 0;
1535 for(s = s00;;s++)
switch(*s) {
1558 while(*++s ==
'0') ;
1564 for(nd = nf = 0; (c = *s) >=
'0' && c <=
'9'; nd++, s++)
1571 s1 = localeconv()->decimal_point;
1592 for(; c ==
'0'; c = *++s)
1594 if (c >
'0' && c <=
'9') {
1602 for(; c >=
'0' && c <=
'9'; c = *++s) {
1607 for(i = 1; i < nz; i++)
1610 else if (nd <= DBL_DIG + 1)
1614 else if (nd <= DBL_DIG + 1)
1622 if (c ==
'e' || c ==
'E') {
1623 if (!nd && !nz && !nz0) {
1634 if (c >=
'0' && c <=
'9') {
1637 if (c >
'0' && c <=
'9') {
1640 while((c = *++s) >=
'0' && c <=
'9')
1642 if (s - s1 > 8 || L > 19999)
1665 if (match(&s,
"nf")) {
1667 if (!match(&s,
"inity"))
1669 word0(rv) = 0x7ff00000;
1676 if (match(&s,
"an")) {
1677 word0(rv) = NAN_WORD0;
1678 word1(rv) = NAN_WORD1;
1702 k = nd < DBL_DIG + 1 ? nd : DBL_DIG + 1;
1707 oldinexact = get_inexact();
1709 dval(rv) = tens[k - 9] * dval(rv) + z;
1713 #ifndef RND_PRODQUOT
1714 #ifndef Honor_FLT_ROUNDS
1722 if (e <= Ten_pmax) {
1724 goto vax_ovfl_check;
1726 #ifdef Honor_FLT_ROUNDS
1733 rounded_product(dval(rv), tens[e]);
1738 if (e <= Ten_pmax + i) {
1742 #ifdef Honor_FLT_ROUNDS
1750 dval(rv) *= tens[i];
1756 word0(rv) -= P*Exp_msk1;
1757 rounded_product(dval(rv), tens[e]);
1758 if ((word0(rv) & Exp_mask)
1759 > Exp_msk1*(DBL_MAX_EXP+Bias-1-P))
1761 word0(rv) += P*Exp_msk1;
1763 rounded_product(dval(rv), tens[e]);
1768 #ifndef Inaccurate_Divide
1769 else if (e >= -Ten_pmax) {
1770 #ifdef Honor_FLT_ROUNDS
1777 rounded_quotient(dval(rv), tens[-e]);
1788 oldinexact = get_inexact();
1790 #ifdef Avoid_Underflow
1793 #ifdef Honor_FLT_ROUNDS
1794 if ((rounding = Flt_Rounds) >= 2) {
1796 rounding = rounding == 2 ? 0 : 2;
1808 dval(rv) *= tens[i];
1810 if (e1 > DBL_MAX_10_EXP) {
1817 #ifdef Honor_FLT_ROUNDS
1825 word0(rv) = Exp_mask;
1829 word0(rv) = Exp_mask;
1835 dval(rv0) *= dval(rv0);
1846 for(j = 0; e1 > 1; j++, e1 >>= 1)
1848 dval(rv) *= bigtens[j];
1850 word0(rv) -= P*Exp_msk1;
1851 dval(rv) *= bigtens[j];
1852 if ((z = word0(rv) & Exp_mask)
1853 > Exp_msk1*(DBL_MAX_EXP+Bias-P))
1855 if (z > Exp_msk1*(DBL_MAX_EXP+Bias-1-P)) {
1862 word0(rv) += P*Exp_msk1;
1868 dval(rv) /= tens[i];
1870 if (e1 >= 1 << n_bigtens)
1872 #ifdef Avoid_Underflow
1875 for(j = 0; e1 > 0; j++, e1 >>= 1)
1877 dval(rv) *= tinytens[j];
1878 if (scale && (j = 2*P + 1 - ((word0(rv) & Exp_mask)
1879 >> Exp_shift)) > 0) {
1884 word0(rv) = (P+2)*Exp_msk1;
1886 word0(rv) &= 0xffffffff << (j-32);
1889 word1(rv) &= 0xffffffff << j;
1892 for(j = 0; e1 > 1; j++, e1 >>= 1)
1894 dval(rv) *= tinytens[j];
1896 dval(rv0) = dval(rv);
1897 dval(rv) *= tinytens[j];
1899 dval(rv) = 2.*dval(rv0);
1900 dval(rv) *= tinytens[j];
1912 #ifndef Avoid_Underflow
1927 bd0 = s2b(s0, nd0, nd, y);
1930 bd = Balloc(bd0->k);
1932 bb = d2b(dval(rv), &bbe, &bbbits);
1948 #ifdef Honor_FLT_ROUNDS
1952 #ifdef Avoid_Underflow
1960 #ifdef Sudden_Underflow
1962 j = 1 + 4*P - 3 - bbbits + ((bbe + bbbits - 1) & 3);
1977 #ifdef Avoid_Underflow
1980 i = bb2 < bd2 ? bb2 : bd2;
1989 bs = pow5mult(bs, bb5);
1995 bb = lshift(bb, bb2);
1997 bd = pow5mult(bd, bd5);
1999 bd = lshift(bd, bd2);
2001 bs = lshift(bs, bs2);
2002 delta = diff(bb, bd);
2003 dsign = delta->sign;
2006 #ifdef Honor_FLT_ROUNDS
2007 if (rounding != 1) {
2010 if (!delta->x[0] && delta->wds <= 1) {
2026 && !(word0(rv) & Frac_mask)) {
2027 y = word0(rv) & Exp_mask;
2028 #ifdef Avoid_Underflow
2029 if (!scale || y > 2*P*Exp_msk1)
2034 delta = lshift(delta,Log2P);
2035 if (cmp(delta, bs) <= 0)
2040 #ifdef Avoid_Underflow
2041 if (scale && (y = word0(rv) & Exp_mask)
2043 word0(adj) += (2*P+1)*Exp_msk1 - y;
2045 #ifdef Sudden_Underflow
2046 if ((word0(rv) & Exp_mask) <=
2048 word0(rv) += P*Exp_msk1;
2049 dval(rv) += adj*ulp(dval(rv));
2050 word0(rv) -= P*Exp_msk1;
2055 dval(rv) += adj*ulp(dval(rv));
2059 adj = ratio(delta, bs);
2062 if (adj <= 0x7ffffffe) {
2066 if (!((rounding>>1) ^ dsign))
2071 #ifdef Avoid_Underflow
2072 if (scale && (y = word0(rv) & Exp_mask) <= 2*P*Exp_msk1)
2073 word0(adj) += (2*P+1)*Exp_msk1 - y;
2075 #ifdef Sudden_Underflow
2076 if ((word0(rv) & Exp_mask) <= P*Exp_msk1) {
2077 word0(rv) += P*Exp_msk1;
2078 adj *= ulp(dval(rv));
2083 word0(rv) -= P*Exp_msk1;
2088 adj *= ulp(dval(rv));
2101 if (dsign || word1(rv) || word0(rv) & Bndry_mask
2103 #ifdef Avoid_Underflow
2104 || (word0(rv) & Exp_mask) <= (2*P+1)*Exp_msk1
2106 || (word0(rv) & Exp_mask) <= Exp_msk1
2111 if (!delta->x[0] && delta->wds <= 1)
2116 if (!delta->x[0] && delta->wds <= 1) {
2123 delta = lshift(delta,Log2P);
2124 if (cmp(delta, bs) > 0)
2131 if ((word0(rv) & Bndry_mask1) == Bndry_mask1
2133 #ifdef Avoid_Underflow
2134 (scale && (y = word0(rv) & Exp_mask) <= 2*P*Exp_msk1)
2135 ? (0xffffffff & (0xffffffff << (2*P+1-(y>>Exp_shift)))) :
2139 word0(rv) = (word0(rv) & Exp_mask)
2146 #ifdef Avoid_Underflow
2152 else if (!(word0(rv) & Bndry_mask) && !word1(rv)) {
2155 #ifdef Sudden_Underflow
2156 L = word0(rv) & Exp_mask;
2160 #ifdef Avoid_Underflow
2161 if (L <= (scale ? (2*P+1)*Exp_msk1 : Exp_msk1))
2169 #ifdef Avoid_Underflow
2171 L = word0(rv) & Exp_mask;
2172 if (L <= (2*P+1)*Exp_msk1) {
2173 if (L > (P+2)*Exp_msk1)
2182 L = (word0(rv) & Exp_mask) - Exp_msk1;
2184 word0(rv) = L | Bndry_mask1;
2185 word1(rv) = 0xffffffff;
2192 #ifndef ROUND_BIASED
2193 if (!(word1(rv) & LSB))
2197 dval(rv) += ulp(dval(rv));
2198 #ifndef ROUND_BIASED
2200 dval(rv) -= ulp(dval(rv));
2201 #ifndef Sudden_Underflow
2206 #ifdef Avoid_Underflow
2212 if ((aadj = ratio(delta, bs)) <= 2.) {
2214 aadj = dval(aadj1) = 1.;
2215 else if (word1(rv) || word0(rv) & Bndry_mask) {
2216 #ifndef Sudden_Underflow
2217 if (word1(rv) == Tiny1 && !word0(rv))
2227 if (aadj < 2./FLT_RADIX)
2228 aadj = 1./FLT_RADIX;
2231 dval(aadj1) = -aadj;
2236 dval(aadj1) = dsign ? aadj : -aadj;
2237 #ifdef Check_FLT_ROUNDS
2247 if (Flt_Rounds == 0)
2251 y = word0(rv) & Exp_mask;
2255 if (y == Exp_msk1*(DBL_MAX_EXP+Bias-1)) {
2256 dval(rv0) = dval(rv);
2257 word0(rv) -= P*Exp_msk1;
2258 adj = dval(aadj1) * ulp(dval(rv));
2260 if ((word0(rv) & Exp_mask) >=
2261 Exp_msk1*(DBL_MAX_EXP+Bias-P)) {
2262 if (word0(rv0) == Big0 && word1(rv0) == Big1)
2269 word0(rv) += P*Exp_msk1;
2272 #ifdef Avoid_Underflow
2273 if (scale && y <= 2*P*Exp_msk1) {
2274 if (aadj <= 0x7fffffff) {
2275 if ((z = (uint32)aadj) <= 0)
2278 dval(aadj1) = dsign ? aadj : -aadj;
2280 word0(aadj1) += (2*P+1)*Exp_msk1 - y;
2282 adj = dval(aadj1) * ulp(dval(rv));
2285 #ifdef Sudden_Underflow
2286 if ((word0(rv) & Exp_mask) <= P*Exp_msk1) {
2287 dval(rv0) = dval(rv);
2288 word0(rv) += P*Exp_msk1;
2289 adj = aadj1 * ulp(dval(rv));
2292 if ((word0(rv) & Exp_mask) < P*Exp_msk1)
2294 if ((word0(rv) & Exp_mask) <= P*Exp_msk1)
2297 if (word0(rv0) == Tiny0
2298 && word1(rv0) == Tiny1)
2305 word0(rv) -= P*Exp_msk1;
2308 adj = aadj1 * ulp(dval(rv));
2319 if (y <= (P-1)*Exp_msk1 && aadj > 1.) {
2320 aadj1 = (double)(
int)(aadj + 0.5);
2324 adj = aadj1 * ulp(dval(rv));
2329 z = word0(rv) & Exp_mask;
2331 #ifdef Avoid_Underflow
2339 if (dsign || word1(rv) || word0(rv) & Bndry_mask) {
2340 if (aadj < .4999999 || aadj > .5000001)
2343 else if (aadj < .4999999/FLT_RADIX)
2356 word0(rv0) = Exp_1 + (70 << Exp_shift);
2361 else if (!oldinexact)
2364 #ifdef Avoid_Underflow
2366 word0(rv0) = Exp_1 - 2*P*Exp_msk1;
2368 dval(rv) *= dval(rv0);
2371 if (word0(rv) == 0 && word1(rv) == 0)
2377 if (inexact && !(word0(rv) & Exp_mask)) {
2380 dval(rv0) *= dval(rv0);
2392 return sign ? -dval(rv) : dval(rv);
Sphinx's memory allocation/deallocation routines.
Basic type definitions used in Sphinx.
SPHINXBASE_EXPORT void ckd_free(void *ptr)
Test and free a 1-D array.
#define ckd_malloc(sz)
Macro for ckd_malloc
Union to extract the bytes of a double.