Deleted Added
sdiff udiff text old ( 187808 ) new ( 219557 )
full compact
strtod.c (187808) strtod.c (219557)
1/****************************************************************
2
3The author of this software is David M. Gay.
4
5Copyright (C) 1998-2001 by Lucent Technologies
6All Rights Reserved
7
8Permission to use, copy, modify, and distribute this software and

--- 15 unchanged lines hidden (view full) ---

24ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF
25THIS SOFTWARE.
26
27****************************************************************/
28
29/* Please send bug reports to David M. Gay (dmg at acm dot org,
30 * with " at " changed at "@" and " dot " changed to "."). */
31
1/****************************************************************
2
3The author of this software is David M. Gay.
4
5Copyright (C) 1998-2001 by Lucent Technologies
6All Rights Reserved
7
8Permission to use, copy, modify, and distribute this software and

--- 15 unchanged lines hidden (view full) ---

24ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF
25THIS SOFTWARE.
26
27****************************************************************/
28
29/* Please send bug reports to David M. Gay (dmg at acm dot org,
30 * with " at " changed at "@" and " dot " changed to "."). */
31
32/* $FreeBSD: head/contrib/gdtoa/strtod.c 187808 2009-01-28 04:36:34Z das $ */
32/* $FreeBSD: head/contrib/gdtoa/strtod.c 219557 2011-03-12 07:03:06Z das $ */
33
34#include "gdtoaimp.h"
35#ifndef NO_FENV_H
36#include <fenv.h>
37#endif
38
39#ifdef USE_LOCALE
40#include "locale.h"

--- 13 unchanged lines hidden (view full) ---

54
55#ifdef Honor_FLT_ROUNDS
56#undef Check_FLT_ROUNDS
57#define Check_FLT_ROUNDS
58#else
59#define Rounding Flt_Rounds
60#endif
61
33
34#include "gdtoaimp.h"
35#ifndef NO_FENV_H
36#include <fenv.h>
37#endif
38
39#ifdef USE_LOCALE
40#include "locale.h"

--- 13 unchanged lines hidden (view full) ---

54
55#ifdef Honor_FLT_ROUNDS
56#undef Check_FLT_ROUNDS
57#define Check_FLT_ROUNDS
58#else
59#define Rounding Flt_Rounds
60#endif
61
62#ifdef Avoid_Underflow /*{*/
63 static double
64sulp
65#ifdef KR_headers
66 (x, scale) U *x; int scale;
67#else
68 (U *x, int scale)
69#endif
70{
71 U u;
72 double rv;
73 int i;
74
75 rv = ulp(x);
76 if (!scale || (i = 2*P + 1 - ((word0(x) & Exp_mask) >> Exp_shift)) <= 0)
77 return rv; /* Is there an example where i <= 0 ? */
78 word0(&u) = Exp_1 + (i << Exp_shift);
79 word1(&u) = 0;
80 return rv * u.d;
81 }
82#endif /*}*/
83
62 double
63strtod
64#ifdef KR_headers
65 (s00, se) CONST char *s00; char **se;
66#else
67 (CONST char *s00, char **se)
68#endif
69{
70#ifdef Avoid_Underflow
71 int scale;
72#endif
73 int bb2, bb5, bbe, bd2, bd5, bbbits, bs2, c, decpt, dsign,
74 e, e1, esign, i, j, k, nd, nd0, nf, nz, nz0, sign;
75 CONST char *s, *s0, *s1;
84 double
85strtod
86#ifdef KR_headers
87 (s00, se) CONST char *s00; char **se;
88#else
89 (CONST char *s00, char **se)
90#endif
91{
92#ifdef Avoid_Underflow
93 int scale;
94#endif
95 int bb2, bb5, bbe, bd2, bd5, bbbits, bs2, c, decpt, dsign,
96 e, e1, esign, i, j, k, nd, nd0, nf, nz, nz0, sign;
97 CONST char *s, *s0, *s1;
76 double aadj, aadj1, adj, rv, rv0;
98 double aadj;
77 Long L;
99 Long L;
100 U adj, aadj1, rv, rv0;
78 ULong y, z;
79 Bigint *bb, *bb1, *bd, *bd0, *bs, *delta;
101 ULong y, z;
102 Bigint *bb, *bb1, *bd, *bd0, *bs, *delta;
103#ifdef Avoid_Underflow
104 ULong Lsb, Lsb1;
105#endif
80#ifdef SET_INEXACT
81 int inexact, oldinexact;
82#endif
83#ifdef USE_LOCALE /*{{*/
84#ifdef NO_LOCALE_CACHE
85 char *decimalpoint = localeconv()->decimal_point;
86 int dplen = strlen(decimalpoint);
87#else
88 char *decimalpoint;
89 static char *decimalpoint_cache;
90 static int dplen;
91 if (!(s0 = decimalpoint_cache)) {
92 s0 = localeconv()->decimal_point;
106#ifdef SET_INEXACT
107 int inexact, oldinexact;
108#endif
109#ifdef USE_LOCALE /*{{*/
110#ifdef NO_LOCALE_CACHE
111 char *decimalpoint = localeconv()->decimal_point;
112 int dplen = strlen(decimalpoint);
113#else
114 char *decimalpoint;
115 static char *decimalpoint_cache;
116 static int dplen;
117 if (!(s0 = decimalpoint_cache)) {
118 s0 = localeconv()->decimal_point;
93 if ((decimalpoint_cache = (char*)malloc(strlen(s0) + 1))) {
119 if ((decimalpoint_cache = (char*)MALLOC(strlen(s0) + 1))) {
94 strcpy(decimalpoint_cache, s0);
95 s0 = decimalpoint_cache;
96 }
97 dplen = strlen(s0);
98 }
99 decimalpoint = (char*)s0;
100#endif /*NO_LOCALE_CACHE*/
101#else /*USE_LOCALE}{*/

--- 10 unchanged lines hidden (view full) ---

112 case FE_TOWARDZERO: Rounding = 0; break;
113 case FE_UPWARD: Rounding = 2; break;
114 case FE_DOWNWARD: Rounding = 3;
115 }
116#endif /*}}*/
117#endif /*}*/
118
119 sign = nz0 = nz = decpt = 0;
120 strcpy(decimalpoint_cache, s0);
121 s0 = decimalpoint_cache;
122 }
123 dplen = strlen(s0);
124 }
125 decimalpoint = (char*)s0;
126#endif /*NO_LOCALE_CACHE*/
127#else /*USE_LOCALE}{*/

--- 10 unchanged lines hidden (view full) ---

138 case FE_TOWARDZERO: Rounding = 0; break;
139 case FE_UPWARD: Rounding = 2; break;
140 case FE_DOWNWARD: Rounding = 3;
141 }
142#endif /*}}*/
143#endif /*}*/
144
145 sign = nz0 = nz = decpt = 0;
120 dval(rv) = 0.;
146 dval(&rv) = 0.;
121 for(s = s00;;s++) switch(*s) {
122 case '-':
123 sign = 1;
124 /* no break */
125 case '+':
126 if (*++s)
127 goto break2;
128 /* no break */

--- 15 unchanged lines hidden (view full) ---

144 {
145 static FPI fpi = { 53, 1-1023-53+1, 2046-1023-53+1, 1, SI };
146 Long exp;
147 ULong bits[2];
148 switch(s[1]) {
149 case 'x':
150 case 'X':
151 {
147 for(s = s00;;s++) switch(*s) {
148 case '-':
149 sign = 1;
150 /* no break */
151 case '+':
152 if (*++s)
153 goto break2;
154 /* no break */

--- 15 unchanged lines hidden (view full) ---

170 {
171 static FPI fpi = { 53, 1-1023-53+1, 2046-1023-53+1, 1, SI };
172 Long exp;
173 ULong bits[2];
174 switch(s[1]) {
175 case 'x':
176 case 'X':
177 {
152#if defined(FE_DOWNWARD) && defined(FE_TONEAREST) && defined(FE_TOWARDZERO) && defined(FE_UPWARD) /*{{*/
178#ifdef Honor_FLT_ROUNDS
153 FPI fpi1 = fpi;
179 FPI fpi1 = fpi;
154#ifdef Honor_FLT_ROUNDS /*{{*/
155 fpi1.rounding = Rounding;
180 fpi1.rounding = Rounding;
156#else /*}{*/
157 switch(fegetround()) {
158 case FE_TOWARDZERO: fpi1.rounding = 0; break;
159 case FE_UPWARD: fpi1.rounding = 2; break;
160 case FE_DOWNWARD: fpi1.rounding = 3;
161 }
162#endif /*}}*/
163#else /*}{*/
181#else
164#define fpi1 fpi
182#define fpi1 fpi
165#endif /*}}*/
183#endif
166 switch((i = gethex(&s, &fpi1, &exp, &bb, sign)) & STRTOG_Retmask) {
167 case STRTOG_NoNumber:
168 s = s00;
169 sign = 0;
170 case STRTOG_Zero:
171 break;
172 default:
173 if (bb) {

--- 108 unchanged lines hidden (view full) ---

282 if (!decpt)
283 switch(c) {
284 case 'i':
285 case 'I':
286 if (match(&s,"nf")) {
287 --s;
288 if (!match(&s,"inity"))
289 ++s;
184 switch((i = gethex(&s, &fpi1, &exp, &bb, sign)) & STRTOG_Retmask) {
185 case STRTOG_NoNumber:
186 s = s00;
187 sign = 0;
188 case STRTOG_Zero:
189 break;
190 default:
191 if (bb) {

--- 108 unchanged lines hidden (view full) ---

300 if (!decpt)
301 switch(c) {
302 case 'i':
303 case 'I':
304 if (match(&s,"nf")) {
305 --s;
306 if (!match(&s,"inity"))
307 ++s;
290 word0(rv) = 0x7ff00000;
291 word1(rv) = 0;
308 word0(&rv) = 0x7ff00000;
309 word1(&rv) = 0;
292 goto ret;
293 }
294 break;
295 case 'n':
296 case 'N':
297 if (match(&s, "an")) {
298#ifndef No_Hex_NaN
299 if (*s == '(' /*)*/
300 && hexnan(&s, &fpinan, bits)
301 == STRTOG_NaNbits) {
310 goto ret;
311 }
312 break;
313 case 'n':
314 case 'N':
315 if (match(&s, "an")) {
316#ifndef No_Hex_NaN
317 if (*s == '(' /*)*/
318 && hexnan(&s, &fpinan, bits)
319 == STRTOG_NaNbits) {
302 word0(rv) = 0x7ff80000 | bits[1];
303 word1(rv) = bits[0];
320 word0(&rv) = 0x7ff80000 | bits[1];
321 word1(&rv) = bits[0];
304 }
305 else {
306#endif
322 }
323 else {
324#endif
307 word0(rv) = NAN_WORD0;
308 word1(rv) = NAN_WORD1;
325 word0(&rv) = NAN_WORD0;
326 word1(&rv) = NAN_WORD1;
309#ifndef No_Hex_NaN
310 }
311#endif
312 goto ret;
313 }
314 }
315#endif /* INFNAN_CHECK */
316 ret0:

--- 7 unchanged lines hidden (view full) ---

324 /* Now we have nd0 digits, starting at s0, followed by a
325 * decimal point, followed by nd-nd0 digits. The number we're
326 * after is the integer represented by those digits times
327 * 10**e */
328
329 if (!nd0)
330 nd0 = nd;
331 k = nd < DBL_DIG + 1 ? nd : DBL_DIG + 1;
327#ifndef No_Hex_NaN
328 }
329#endif
330 goto ret;
331 }
332 }
333#endif /* INFNAN_CHECK */
334 ret0:

--- 7 unchanged lines hidden (view full) ---

342 /* Now we have nd0 digits, starting at s0, followed by a
343 * decimal point, followed by nd-nd0 digits. The number we're
344 * after is the integer represented by those digits times
345 * 10**e */
346
347 if (!nd0)
348 nd0 = nd;
349 k = nd < DBL_DIG + 1 ? nd : DBL_DIG + 1;
332 dval(rv) = y;
350 dval(&rv) = y;
333 if (k > 9) {
334#ifdef SET_INEXACT
335 if (k > DBL_DIG)
336 oldinexact = get_inexact();
337#endif
351 if (k > 9) {
352#ifdef SET_INEXACT
353 if (k > DBL_DIG)
354 oldinexact = get_inexact();
355#endif
338 dval(rv) = tens[k - 9] * dval(rv) + z;
356 dval(&rv) = tens[k - 9] * dval(&rv) + z;
339 }
340 bd0 = 0;
341 if (nd <= DBL_DIG
342#ifndef RND_PRODQUOT
343#ifndef Honor_FLT_ROUNDS
344 && Flt_Rounds == 1
345#endif
346#endif
347 ) {
348 if (!e)
349 goto ret;
357 }
358 bd0 = 0;
359 if (nd <= DBL_DIG
360#ifndef RND_PRODQUOT
361#ifndef Honor_FLT_ROUNDS
362 && Flt_Rounds == 1
363#endif
364#endif
365 ) {
366 if (!e)
367 goto ret;
368#ifndef ROUND_BIASED_without_Round_Up
350 if (e > 0) {
351 if (e <= Ten_pmax) {
352#ifdef VAX
353 goto vax_ovfl_check;
354#else
355#ifdef Honor_FLT_ROUNDS
356 /* round correctly FLT_ROUNDS = 2 or 3 */
357 if (sign) {
369 if (e > 0) {
370 if (e <= Ten_pmax) {
371#ifdef VAX
372 goto vax_ovfl_check;
373#else
374#ifdef Honor_FLT_ROUNDS
375 /* round correctly FLT_ROUNDS = 2 or 3 */
376 if (sign) {
358 rv = -rv;
377 rv.d = -rv.d;
359 sign = 0;
360 }
361#endif
378 sign = 0;
379 }
380#endif
362 /* rv = */ rounded_product(dval(rv), tens[e]);
381 /* rv = */ rounded_product(dval(&rv), tens[e]);
363 goto ret;
364#endif
365 }
366 i = DBL_DIG - nd;
367 if (e <= Ten_pmax + i) {
368 /* A fancier test would sometimes let us do
369 * this for larger i values.
370 */
371#ifdef Honor_FLT_ROUNDS
372 /* round correctly FLT_ROUNDS = 2 or 3 */
373 if (sign) {
382 goto ret;
383#endif
384 }
385 i = DBL_DIG - nd;
386 if (e <= Ten_pmax + i) {
387 /* A fancier test would sometimes let us do
388 * this for larger i values.
389 */
390#ifdef Honor_FLT_ROUNDS
391 /* round correctly FLT_ROUNDS = 2 or 3 */
392 if (sign) {
374 rv = -rv;
393 rv.d = -rv.d;
375 sign = 0;
376 }
377#endif
378 e -= i;
394 sign = 0;
395 }
396#endif
397 e -= i;
379 dval(rv) *= tens[i];
398 dval(&rv) *= tens[i];
380#ifdef VAX
381 /* VAX exponent range is so narrow we must
382 * worry about overflow here...
383 */
384 vax_ovfl_check:
399#ifdef VAX
400 /* VAX exponent range is so narrow we must
401 * worry about overflow here...
402 */
403 vax_ovfl_check:
385 word0(rv) -= P*Exp_msk1;
386 /* rv = */ rounded_product(dval(rv), tens[e]);
387 if ((word0(rv) & Exp_mask)
404 word0(&rv) -= P*Exp_msk1;
405 /* rv = */ rounded_product(dval(&rv), tens[e]);
406 if ((word0(&rv) & Exp_mask)
388 > Exp_msk1*(DBL_MAX_EXP+Bias-1-P))
389 goto ovfl;
407 > Exp_msk1*(DBL_MAX_EXP+Bias-1-P))
408 goto ovfl;
390 word0(rv) += P*Exp_msk1;
409 word0(&rv) += P*Exp_msk1;
391#else
410#else
392 /* rv = */ rounded_product(dval(rv), tens[e]);
411 /* rv = */ rounded_product(dval(&rv), tens[e]);
393#endif
394 goto ret;
395 }
396 }
397#ifndef Inaccurate_Divide
398 else if (e >= -Ten_pmax) {
399#ifdef Honor_FLT_ROUNDS
400 /* round correctly FLT_ROUNDS = 2 or 3 */
401 if (sign) {
412#endif
413 goto ret;
414 }
415 }
416#ifndef Inaccurate_Divide
417 else if (e >= -Ten_pmax) {
418#ifdef Honor_FLT_ROUNDS
419 /* round correctly FLT_ROUNDS = 2 or 3 */
420 if (sign) {
402 rv = -rv;
421 rv.d = -rv.d;
403 sign = 0;
404 }
405#endif
422 sign = 0;
423 }
424#endif
406 /* rv = */ rounded_quotient(dval(rv), tens[-e]);
425 /* rv = */ rounded_quotient(dval(&rv), tens[-e]);
407 goto ret;
408 }
409#endif
426 goto ret;
427 }
428#endif
429#endif /* ROUND_BIASED_without_Round_Up */
410 }
411 e1 += nd - k;
412
413#ifdef IEEE_Arith
414#ifdef SET_INEXACT
415 inexact = 1;
416 if (k <= DBL_DIG)
417 oldinexact = get_inexact();

--- 11 unchanged lines hidden (view full) ---

429 }
430#endif
431#endif /*IEEE_Arith*/
432
433 /* Get starting approximation = rv * 10**e1 */
434
435 if (e1 > 0) {
436 if ( (i = e1 & 15) !=0)
430 }
431 e1 += nd - k;
432
433#ifdef IEEE_Arith
434#ifdef SET_INEXACT
435 inexact = 1;
436 if (k <= DBL_DIG)
437 oldinexact = get_inexact();

--- 11 unchanged lines hidden (view full) ---

449 }
450#endif
451#endif /*IEEE_Arith*/
452
453 /* Get starting approximation = rv * 10**e1 */
454
455 if (e1 > 0) {
456 if ( (i = e1 & 15) !=0)
437 dval(rv) *= tens[i];
457 dval(&rv) *= tens[i];
438 if (e1 &= ~15) {
439 if (e1 > DBL_MAX_10_EXP) {
440 ovfl:
458 if (e1 &= ~15) {
459 if (e1 > DBL_MAX_10_EXP) {
460 ovfl:
441#ifndef NO_ERRNO
442 errno = ERANGE;
443#endif
444 /* Can't trust HUGE_VAL */
445#ifdef IEEE_Arith
446#ifdef Honor_FLT_ROUNDS
447 switch(Rounding) {
448 case 0: /* toward 0 */
449 case 3: /* toward -infinity */
461 /* Can't trust HUGE_VAL */
462#ifdef IEEE_Arith
463#ifdef Honor_FLT_ROUNDS
464 switch(Rounding) {
465 case 0: /* toward 0 */
466 case 3: /* toward -infinity */
450 word0(rv) = Big0;
451 word1(rv) = Big1;
467 word0(&rv) = Big0;
468 word1(&rv) = Big1;
452 break;
453 default:
469 break;
470 default:
454 word0(rv) = Exp_mask;
455 word1(rv) = 0;
471 word0(&rv) = Exp_mask;
472 word1(&rv) = 0;
456 }
457#else /*Honor_FLT_ROUNDS*/
473 }
474#else /*Honor_FLT_ROUNDS*/
458 word0(rv) = Exp_mask;
459 word1(rv) = 0;
475 word0(&rv) = Exp_mask;
476 word1(&rv) = 0;
460#endif /*Honor_FLT_ROUNDS*/
461#ifdef SET_INEXACT
462 /* set overflow bit */
477#endif /*Honor_FLT_ROUNDS*/
478#ifdef SET_INEXACT
479 /* set overflow bit */
463 dval(rv0) = 1e300;
464 dval(rv0) *= dval(rv0);
480 dval(&rv0) = 1e300;
481 dval(&rv0) *= dval(&rv0);
465#endif
466#else /*IEEE_Arith*/
482#endif
483#else /*IEEE_Arith*/
467 word0(rv) = Big0;
468 word1(rv) = Big1;
484 word0(&rv) = Big0;
485 word1(&rv) = Big1;
469#endif /*IEEE_Arith*/
486#endif /*IEEE_Arith*/
470 if (bd0)
471 goto retfree;
487 range_err:
488 if (bd0) {
489 Bfree(bb);
490 Bfree(bd);
491 Bfree(bs);
492 Bfree(bd0);
493 Bfree(delta);
494 }
495#ifndef NO_ERRNO
496 errno = ERANGE;
497#endif
472 goto ret;
473 }
474 e1 >>= 4;
475 for(j = 0; e1 > 1; j++, e1 >>= 1)
476 if (e1 & 1)
498 goto ret;
499 }
500 e1 >>= 4;
501 for(j = 0; e1 > 1; j++, e1 >>= 1)
502 if (e1 & 1)
477 dval(rv) *= bigtens[j];
503 dval(&rv) *= bigtens[j];
478 /* The last multiplication could overflow. */
504 /* The last multiplication could overflow. */
479 word0(rv) -= P*Exp_msk1;
480 dval(rv) *= bigtens[j];
481 if ((z = word0(rv) & Exp_mask)
505 word0(&rv) -= P*Exp_msk1;
506 dval(&rv) *= bigtens[j];
507 if ((z = word0(&rv) & Exp_mask)
482 > Exp_msk1*(DBL_MAX_EXP+Bias-P))
483 goto ovfl;
484 if (z > Exp_msk1*(DBL_MAX_EXP+Bias-1-P)) {
485 /* set to largest number */
486 /* (Can't trust DBL_MAX) */
508 > Exp_msk1*(DBL_MAX_EXP+Bias-P))
509 goto ovfl;
510 if (z > Exp_msk1*(DBL_MAX_EXP+Bias-1-P)) {
511 /* set to largest number */
512 /* (Can't trust DBL_MAX) */
487 word0(rv) = Big0;
488 word1(rv) = Big1;
513 word0(&rv) = Big0;
514 word1(&rv) = Big1;
489 }
490 else
515 }
516 else
491 word0(rv) += P*Exp_msk1;
517 word0(&rv) += P*Exp_msk1;
492 }
493 }
494 else if (e1 < 0) {
495 e1 = -e1;
496 if ( (i = e1 & 15) !=0)
518 }
519 }
520 else if (e1 < 0) {
521 e1 = -e1;
522 if ( (i = e1 & 15) !=0)
497 dval(rv) /= tens[i];
523 dval(&rv) /= tens[i];
498 if (e1 >>= 4) {
499 if (e1 >= 1 << n_bigtens)
500 goto undfl;
501#ifdef Avoid_Underflow
502 if (e1 & Scale_Bit)
503 scale = 2*P;
504 for(j = 0; e1 > 0; j++, e1 >>= 1)
505 if (e1 & 1)
524 if (e1 >>= 4) {
525 if (e1 >= 1 << n_bigtens)
526 goto undfl;
527#ifdef Avoid_Underflow
528 if (e1 & Scale_Bit)
529 scale = 2*P;
530 for(j = 0; e1 > 0; j++, e1 >>= 1)
531 if (e1 & 1)
506 dval(rv) *= tinytens[j];
507 if (scale && (j = 2*P + 1 - ((word0(rv) & Exp_mask)
532 dval(&rv) *= tinytens[j];
533 if (scale && (j = 2*P + 1 - ((word0(&rv) & Exp_mask)
508 >> Exp_shift)) > 0) {
509 /* scaled rv is denormal; zap j low bits */
510 if (j >= 32) {
534 >> Exp_shift)) > 0) {
535 /* scaled rv is denormal; zap j low bits */
536 if (j >= 32) {
511 word1(rv) = 0;
537 word1(&rv) = 0;
512 if (j >= 53)
538 if (j >= 53)
513 word0(rv) = (P+2)*Exp_msk1;
539 word0(&rv) = (P+2)*Exp_msk1;
514 else
540 else
515 word0(rv) &= 0xffffffff << j-32;
541 word0(&rv) &= 0xffffffff << (j-32);
516 }
517 else
542 }
543 else
518 word1(rv) &= 0xffffffff << j;
544 word1(&rv) &= 0xffffffff << j;
519 }
520#else
521 for(j = 0; e1 > 1; j++, e1 >>= 1)
522 if (e1 & 1)
545 }
546#else
547 for(j = 0; e1 > 1; j++, e1 >>= 1)
548 if (e1 & 1)
523 dval(rv) *= tinytens[j];
549 dval(&rv) *= tinytens[j];
524 /* The last multiplication could underflow. */
550 /* The last multiplication could underflow. */
525 dval(rv0) = dval(rv);
526 dval(rv) *= tinytens[j];
527 if (!dval(rv)) {
528 dval(rv) = 2.*dval(rv0);
529 dval(rv) *= tinytens[j];
551 dval(&rv0) = dval(&rv);
552 dval(&rv) *= tinytens[j];
553 if (!dval(&rv)) {
554 dval(&rv) = 2.*dval(&rv0);
555 dval(&rv) *= tinytens[j];
530#endif
556#endif
531 if (!dval(rv)) {
557 if (!dval(&rv)) {
532 undfl:
558 undfl:
533 dval(rv) = 0.;
534#ifndef NO_ERRNO
535 errno = ERANGE;
536#endif
537 if (bd0)
538 goto retfree;
539 goto ret;
559 dval(&rv) = 0.;
560 goto range_err;
540 }
541#ifndef Avoid_Underflow
561 }
562#ifndef Avoid_Underflow
542 word0(rv) = Tiny0;
543 word1(rv) = Tiny1;
563 word0(&rv) = Tiny0;
564 word1(&rv) = Tiny1;
544 /* The refinement below will clean
545 * this approximation up.
546 */
547 }
548#endif
549 }
550 }
551
552 /* Now the hard part -- adjusting rv to the correct value.*/
553
554 /* Put digits into bd: true value = bd * 10^e */
555
556 bd0 = s2b(s0, nd0, nd, y, dplen);
557
558 for(;;) {
559 bd = Balloc(bd0->k);
560 Bcopy(bd, bd0);
565 /* The refinement below will clean
566 * this approximation up.
567 */
568 }
569#endif
570 }
571 }
572
573 /* Now the hard part -- adjusting rv to the correct value.*/
574
575 /* Put digits into bd: true value = bd * 10^e */
576
577 bd0 = s2b(s0, nd0, nd, y, dplen);
578
579 for(;;) {
580 bd = Balloc(bd0->k);
581 Bcopy(bd, bd0);
561 bb = d2b(dval(rv), &bbe, &bbbits); /* rv = bb * 2^bbe */
582 bb = d2b(dval(&rv), &bbe, &bbbits); /* rv = bb * 2^bbe */
562 bs = i2b(1);
563
564 if (e >= 0) {
565 bb2 = bb5 = 0;
566 bd2 = bd5 = e;
567 }
568 else {
569 bb2 = bb5 = -e;

--- 4 unchanged lines hidden (view full) ---

574 else
575 bd2 -= bbe;
576 bs2 = bb2;
577#ifdef Honor_FLT_ROUNDS
578 if (Rounding != 1)
579 bs2++;
580#endif
581#ifdef Avoid_Underflow
583 bs = i2b(1);
584
585 if (e >= 0) {
586 bb2 = bb5 = 0;
587 bd2 = bd5 = e;
588 }
589 else {
590 bb2 = bb5 = -e;

--- 4 unchanged lines hidden (view full) ---

595 else
596 bd2 -= bbe;
597 bs2 = bb2;
598#ifdef Honor_FLT_ROUNDS
599 if (Rounding != 1)
600 bs2++;
601#endif
602#ifdef Avoid_Underflow
603 Lsb = LSB;
604 Lsb1 = 0;
582 j = bbe - scale;
583 i = j + bbbits - 1; /* logb(rv) */
605 j = bbe - scale;
606 i = j + bbbits - 1; /* logb(rv) */
584 if (i < Emin) /* denormal */
585 j += P - Emin;
586 else
587 j = P + 1 - bbbits;
607 j = P + 1 - bbbits;
608 if (i < Emin) { /* denormal */
609 i = Emin - i;
610 j -= i;
611 if (i < 32)
612 Lsb <<= i;
613 else
614 Lsb1 = Lsb << (i-32);
615 }
588#else /*Avoid_Underflow*/
589#ifdef Sudden_Underflow
590#ifdef IBM
591 j = 1 + 4*P - 3 - bbbits + ((bbe + bbbits - 1) & 3);
592#else
593 j = P + 1 - bbbits;
594#endif
595#else /*Sudden_Underflow*/
596 j = bbe;
616#else /*Avoid_Underflow*/
617#ifdef Sudden_Underflow
618#ifdef IBM
619 j = 1 + 4*P - 3 - bbbits + ((bbe + bbbits - 1) & 3);
620#else
621 j = P + 1 - bbbits;
622#endif
623#else /*Sudden_Underflow*/
624 j = bbe;
597 i = j + bbbits - 1; /* logb(rv) */
625 i = j + bbbits - 1; /* logb(&rv) */
598 if (i < Emin) /* denormal */
599 j += P - Emin;
600 else
601 j = P + 1 - bbbits;
602#endif /*Sudden_Underflow*/
603#endif /*Avoid_Underflow*/
604 bb2 += j;
605 bd2 += j;

--- 34 unchanged lines hidden (view full) ---

640 /* exact */
641#ifdef SET_INEXACT
642 inexact = 0;
643#endif
644 break;
645 }
646 if (Rounding) {
647 if (dsign) {
626 if (i < Emin) /* denormal */
627 j += P - Emin;
628 else
629 j = P + 1 - bbbits;
630#endif /*Sudden_Underflow*/
631#endif /*Avoid_Underflow*/
632 bb2 += j;
633 bd2 += j;

--- 34 unchanged lines hidden (view full) ---

668 /* exact */
669#ifdef SET_INEXACT
670 inexact = 0;
671#endif
672 break;
673 }
674 if (Rounding) {
675 if (dsign) {
648 adj = 1.;
676 dval(&adj) = 1.;
649 goto apply_adj;
650 }
651 }
652 else if (!dsign) {
677 goto apply_adj;
678 }
679 }
680 else if (!dsign) {
653 adj = -1.;
654 if (!word1(rv)
655 && !(word0(rv) & Frac_mask)) {
656 y = word0(rv) & Exp_mask;
681 dval(&adj) = -1.;
682 if (!word1(&rv)
683 && !(word0(&rv) & Frac_mask)) {
684 y = word0(&rv) & Exp_mask;
657#ifdef Avoid_Underflow
658 if (!scale || y > 2*P*Exp_msk1)
659#else
660 if (y)
661#endif
662 {
663 delta = lshift(delta,Log2P);
664 if (cmp(delta, bs) <= 0)
685#ifdef Avoid_Underflow
686 if (!scale || y > 2*P*Exp_msk1)
687#else
688 if (y)
689#endif
690 {
691 delta = lshift(delta,Log2P);
692 if (cmp(delta, bs) <= 0)
665 adj = -0.5;
693 dval(&adj) = -0.5;
666 }
667 }
668 apply_adj:
669#ifdef Avoid_Underflow
694 }
695 }
696 apply_adj:
697#ifdef Avoid_Underflow
670 if (scale && (y = word0(rv) & Exp_mask)
698 if (scale && (y = word0(&rv) & Exp_mask)
671 <= 2*P*Exp_msk1)
699 <= 2*P*Exp_msk1)
672 word0(adj) += (2*P+1)*Exp_msk1 - y;
700 word0(&adj) += (2*P+1)*Exp_msk1 - y;
673#else
674#ifdef Sudden_Underflow
701#else
702#ifdef Sudden_Underflow
675 if ((word0(rv) & Exp_mask) <=
703 if ((word0(&rv) & Exp_mask) <=
676 P*Exp_msk1) {
704 P*Exp_msk1) {
677 word0(rv) += P*Exp_msk1;
678 dval(rv) += adj*ulp(dval(rv));
679 word0(rv) -= P*Exp_msk1;
705 word0(&rv) += P*Exp_msk1;
706 dval(&rv) += adj*ulp(&rv);
707 word0(&rv) -= P*Exp_msk1;
680 }
681 else
682#endif /*Sudden_Underflow*/
683#endif /*Avoid_Underflow*/
708 }
709 else
710#endif /*Sudden_Underflow*/
711#endif /*Avoid_Underflow*/
684 dval(rv) += adj*ulp(dval(rv));
712 dval(&rv) += adj.d*ulp(&rv);
685 }
686 break;
687 }
713 }
714 break;
715 }
688 adj = ratio(delta, bs);
689 if (adj < 1.)
690 adj = 1.;
691 if (adj <= 0x7ffffffe) {
692 /* adj = Rounding ? ceil(adj) : floor(adj); */
693 y = adj;
694 if (y != adj) {
716 dval(&adj) = ratio(delta, bs);
717 if (adj.d < 1.)
718 dval(&adj) = 1.;
719 if (adj.d <= 0x7ffffffe) {
720 /* dval(&adj) = Rounding ? ceil(&adj) : floor(&adj); */
721 y = adj.d;
722 if (y != adj.d) {
695 if (!((Rounding>>1) ^ dsign))
696 y++;
723 if (!((Rounding>>1) ^ dsign))
724 y++;
697 adj = y;
725 dval(&adj) = y;
698 }
699 }
700#ifdef Avoid_Underflow
726 }
727 }
728#ifdef Avoid_Underflow
701 if (scale && (y = word0(rv) & Exp_mask) <= 2*P*Exp_msk1)
702 word0(adj) += (2*P+1)*Exp_msk1 - y;
729 if (scale && (y = word0(&rv) & Exp_mask) <= 2*P*Exp_msk1)
730 word0(&adj) += (2*P+1)*Exp_msk1 - y;
703#else
704#ifdef Sudden_Underflow
731#else
732#ifdef Sudden_Underflow
705 if ((word0(rv) & Exp_mask) <= P*Exp_msk1) {
706 word0(rv) += P*Exp_msk1;
707 adj *= ulp(dval(rv));
733 if ((word0(&rv) & Exp_mask) <= P*Exp_msk1) {
734 word0(&rv) += P*Exp_msk1;
735 dval(&adj) *= ulp(&rv);
708 if (dsign)
736 if (dsign)
709 dval(rv) += adj;
737 dval(&rv) += adj;
710 else
738 else
711 dval(rv) -= adj;
712 word0(rv) -= P*Exp_msk1;
739 dval(&rv) -= adj;
740 word0(&rv) -= P*Exp_msk1;
713 goto cont;
714 }
715#endif /*Sudden_Underflow*/
716#endif /*Avoid_Underflow*/
741 goto cont;
742 }
743#endif /*Sudden_Underflow*/
744#endif /*Avoid_Underflow*/
717 adj *= ulp(dval(rv));
745 dval(&adj) *= ulp(&rv);
718 if (dsign) {
746 if (dsign) {
719 if (word0(rv) == Big0 && word1(rv) == Big1)
747 if (word0(&rv) == Big0 && word1(&rv) == Big1)
720 goto ovfl;
748 goto ovfl;
721 dval(rv) += adj;
749 dval(&rv) += adj.d;
722 }
723 else
750 }
751 else
724 dval(rv) -= adj;
752 dval(&rv) -= adj.d;
725 goto cont;
726 }
727#endif /*Honor_FLT_ROUNDS*/
728
729 if (i < 0) {
730 /* Error is less than half an ulp -- check for
731 * special case of mantissa a power of two.
732 */
753 goto cont;
754 }
755#endif /*Honor_FLT_ROUNDS*/
756
757 if (i < 0) {
758 /* Error is less than half an ulp -- check for
759 * special case of mantissa a power of two.
760 */
733 if (dsign || word1(rv) || word0(rv) & Bndry_mask
761 if (dsign || word1(&rv) || word0(&rv) & Bndry_mask
734#ifdef IEEE_Arith
735#ifdef Avoid_Underflow
762#ifdef IEEE_Arith
763#ifdef Avoid_Underflow
736 || (word0(rv) & Exp_mask) <= (2*P+1)*Exp_msk1
764 || (word0(&rv) & Exp_mask) <= (2*P+1)*Exp_msk1
737#else
765#else
738 || (word0(rv) & Exp_mask) <= Exp_msk1
766 || (word0(&rv) & Exp_mask) <= Exp_msk1
739#endif
740#endif
741 ) {
742#ifdef SET_INEXACT
743 if (!delta->x[0] && delta->wds <= 1)
744 inexact = 0;
745#endif
746 break;

--- 8 unchanged lines hidden (view full) ---

755 delta = lshift(delta,Log2P);
756 if (cmp(delta, bs) > 0)
757 goto drop_down;
758 break;
759 }
760 if (i == 0) {
761 /* exactly half-way between */
762 if (dsign) {
767#endif
768#endif
769 ) {
770#ifdef SET_INEXACT
771 if (!delta->x[0] && delta->wds <= 1)
772 inexact = 0;
773#endif
774 break;

--- 8 unchanged lines hidden (view full) ---

783 delta = lshift(delta,Log2P);
784 if (cmp(delta, bs) > 0)
785 goto drop_down;
786 break;
787 }
788 if (i == 0) {
789 /* exactly half-way between */
790 if (dsign) {
763 if ((word0(rv) & Bndry_mask1) == Bndry_mask1
764 && word1(rv) == (
791 if ((word0(&rv) & Bndry_mask1) == Bndry_mask1
792 && word1(&rv) == (
765#ifdef Avoid_Underflow
793#ifdef Avoid_Underflow
766 (scale && (y = word0(rv) & Exp_mask) <= 2*P*Exp_msk1)
794 (scale && (y = word0(&rv) & Exp_mask) <= 2*P*Exp_msk1)
767 ? (0xffffffff & (0xffffffff << (2*P+1-(y>>Exp_shift)))) :
768#endif
769 0xffffffff)) {
770 /*boundary case -- increment exponent*/
795 ? (0xffffffff & (0xffffffff << (2*P+1-(y>>Exp_shift)))) :
796#endif
797 0xffffffff)) {
798 /*boundary case -- increment exponent*/
771 word0(rv) = (word0(rv) & Exp_mask)
799 if (word0(&rv) == Big0 && word1(&rv) == Big1)
800 goto ovfl;
801 word0(&rv) = (word0(&rv) & Exp_mask)
772 + Exp_msk1
773#ifdef IBM
774 | Exp_msk1 >> 4
775#endif
776 ;
802 + Exp_msk1
803#ifdef IBM
804 | Exp_msk1 >> 4
805#endif
806 ;
777 word1(rv) = 0;
807 word1(&rv) = 0;
778#ifdef Avoid_Underflow
779 dsign = 0;
780#endif
781 break;
782 }
783 }
808#ifdef Avoid_Underflow
809 dsign = 0;
810#endif
811 break;
812 }
813 }
784 else if (!(word0(rv) & Bndry_mask) && !word1(rv)) {
814 else if (!(word0(&rv) & Bndry_mask) && !word1(&rv)) {
785 drop_down:
786 /* boundary case -- decrement exponent */
787#ifdef Sudden_Underflow /*{{*/
815 drop_down:
816 /* boundary case -- decrement exponent */
817#ifdef Sudden_Underflow /*{{*/
788 L = word0(rv) & Exp_mask;
818 L = word0(&rv) & Exp_mask;
789#ifdef IBM
790 if (L < Exp_msk1)
791#else
792#ifdef Avoid_Underflow
793 if (L <= (scale ? (2*P+1)*Exp_msk1 : Exp_msk1))
794#else
795 if (L <= Exp_msk1)
796#endif /*Avoid_Underflow*/
797#endif /*IBM*/
798 goto undfl;
799 L -= Exp_msk1;
800#else /*Sudden_Underflow}{*/
801#ifdef Avoid_Underflow
802 if (scale) {
819#ifdef IBM
820 if (L < Exp_msk1)
821#else
822#ifdef Avoid_Underflow
823 if (L <= (scale ? (2*P+1)*Exp_msk1 : Exp_msk1))
824#else
825 if (L <= Exp_msk1)
826#endif /*Avoid_Underflow*/
827#endif /*IBM*/
828 goto undfl;
829 L -= Exp_msk1;
830#else /*Sudden_Underflow}{*/
831#ifdef Avoid_Underflow
832 if (scale) {
803 L = word0(rv) & Exp_mask;
833 L = word0(&rv) & Exp_mask;
804 if (L <= (2*P+1)*Exp_msk1) {
805 if (L > (P+2)*Exp_msk1)
806 /* round even ==> */
807 /* accept rv */
808 break;
809 /* rv = smallest denormal */
810 goto undfl;
811 }
812 }
813#endif /*Avoid_Underflow*/
834 if (L <= (2*P+1)*Exp_msk1) {
835 if (L > (P+2)*Exp_msk1)
836 /* round even ==> */
837 /* accept rv */
838 break;
839 /* rv = smallest denormal */
840 goto undfl;
841 }
842 }
843#endif /*Avoid_Underflow*/
814 L = (word0(rv) & Exp_mask) - Exp_msk1;
844 L = (word0(&rv) & Exp_mask) - Exp_msk1;
815#endif /*Sudden_Underflow}}*/
845#endif /*Sudden_Underflow}}*/
816 word0(rv) = L | Bndry_mask1;
817 word1(rv) = 0xffffffff;
846 word0(&rv) = L | Bndry_mask1;
847 word1(&rv) = 0xffffffff;
818#ifdef IBM
819 goto cont;
820#else
821 break;
822#endif
823 }
824#ifndef ROUND_BIASED
848#ifdef IBM
849 goto cont;
850#else
851 break;
852#endif
853 }
854#ifndef ROUND_BIASED
825 if (!(word1(rv) & LSB))
855#ifdef Avoid_Underflow
856 if (Lsb1) {
857 if (!(word0(&rv) & Lsb1))
858 break;
859 }
860 else if (!(word1(&rv) & Lsb))
826 break;
861 break;
862#else
863 if (!(word1(&rv) & LSB))
864 break;
827#endif
865#endif
866#endif
828 if (dsign)
867 if (dsign)
829 dval(rv) += ulp(dval(rv));
868#ifdef Avoid_Underflow
869 dval(&rv) += sulp(&rv, scale);
870#else
871 dval(&rv) += ulp(&rv);
872#endif
830#ifndef ROUND_BIASED
831 else {
873#ifndef ROUND_BIASED
874 else {
832 dval(rv) -= ulp(dval(rv));
875#ifdef Avoid_Underflow
876 dval(&rv) -= sulp(&rv, scale);
877#else
878 dval(&rv) -= ulp(&rv);
879#endif
833#ifndef Sudden_Underflow
880#ifndef Sudden_Underflow
834 if (!dval(rv))
881 if (!dval(&rv))
835 goto undfl;
836#endif
837 }
838#ifdef Avoid_Underflow
839 dsign = 1 - dsign;
840#endif
841#endif
842 break;
843 }
844 if ((aadj = ratio(delta, bs)) <= 2.) {
845 if (dsign)
882 goto undfl;
883#endif
884 }
885#ifdef Avoid_Underflow
886 dsign = 1 - dsign;
887#endif
888#endif
889 break;
890 }
891 if ((aadj = ratio(delta, bs)) <= 2.) {
892 if (dsign)
846 aadj = aadj1 = 1.;
847 else if (word1(rv) || word0(rv) & Bndry_mask) {
893 aadj = dval(&aadj1) = 1.;
894 else if (word1(&rv) || word0(&rv) & Bndry_mask) {
848#ifndef Sudden_Underflow
895#ifndef Sudden_Underflow
849 if (word1(rv) == Tiny1 && !word0(rv))
896 if (word1(&rv) == Tiny1 && !word0(&rv))
850 goto undfl;
851#endif
852 aadj = 1.;
897 goto undfl;
898#endif
899 aadj = 1.;
853 aadj1 = -1.;
900 dval(&aadj1) = -1.;
854 }
855 else {
856 /* special case -- power of FLT_RADIX to be */
857 /* rounded down... */
858
859 if (aadj < 2./FLT_RADIX)
860 aadj = 1./FLT_RADIX;
861 else
862 aadj *= 0.5;
901 }
902 else {
903 /* special case -- power of FLT_RADIX to be */
904 /* rounded down... */
905
906 if (aadj < 2./FLT_RADIX)
907 aadj = 1./FLT_RADIX;
908 else
909 aadj *= 0.5;
863 aadj1 = -aadj;
910 dval(&aadj1) = -aadj;
864 }
865 }
866 else {
867 aadj *= 0.5;
911 }
912 }
913 else {
914 aadj *= 0.5;
868 aadj1 = dsign ? aadj : -aadj;
915 dval(&aadj1) = dsign ? aadj : -aadj;
869#ifdef Check_FLT_ROUNDS
870 switch(Rounding) {
871 case 2: /* towards +infinity */
916#ifdef Check_FLT_ROUNDS
917 switch(Rounding) {
918 case 2: /* towards +infinity */
872 aadj1 -= 0.5;
919 dval(&aadj1) -= 0.5;
873 break;
874 case 0: /* towards 0 */
875 case 3: /* towards -infinity */
920 break;
921 case 0: /* towards 0 */
922 case 3: /* towards -infinity */
876 aadj1 += 0.5;
923 dval(&aadj1) += 0.5;
877 }
878#else
879 if (Flt_Rounds == 0)
924 }
925#else
926 if (Flt_Rounds == 0)
880 aadj1 += 0.5;
927 dval(&aadj1) += 0.5;
881#endif /*Check_FLT_ROUNDS*/
882 }
928#endif /*Check_FLT_ROUNDS*/
929 }
883 y = word0(rv) & Exp_mask;
930 y = word0(&rv) & Exp_mask;
884
885 /* Check for overflow */
886
887 if (y == Exp_msk1*(DBL_MAX_EXP+Bias-1)) {
931
932 /* Check for overflow */
933
934 if (y == Exp_msk1*(DBL_MAX_EXP+Bias-1)) {
888 dval(rv0) = dval(rv);
889 word0(rv) -= P*Exp_msk1;
890 adj = aadj1 * ulp(dval(rv));
891 dval(rv) += adj;
892 if ((word0(rv) & Exp_mask) >=
935 dval(&rv0) = dval(&rv);
936 word0(&rv) -= P*Exp_msk1;
937 dval(&adj) = dval(&aadj1) * ulp(&rv);
938 dval(&rv) += dval(&adj);
939 if ((word0(&rv) & Exp_mask) >=
893 Exp_msk1*(DBL_MAX_EXP+Bias-P)) {
940 Exp_msk1*(DBL_MAX_EXP+Bias-P)) {
894 if (word0(rv0) == Big0 && word1(rv0) == Big1)
941 if (word0(&rv0) == Big0 && word1(&rv0) == Big1)
895 goto ovfl;
942 goto ovfl;
896 word0(rv) = Big0;
897 word1(rv) = Big1;
943 word0(&rv) = Big0;
944 word1(&rv) = Big1;
898 goto cont;
899 }
900 else
945 goto cont;
946 }
947 else
901 word0(rv) += P*Exp_msk1;
948 word0(&rv) += P*Exp_msk1;
902 }
903 else {
904#ifdef Avoid_Underflow
905 if (scale && y <= 2*P*Exp_msk1) {
906 if (aadj <= 0x7fffffff) {
907 if ((z = aadj) <= 0)
908 z = 1;
909 aadj = z;
949 }
950 else {
951#ifdef Avoid_Underflow
952 if (scale && y <= 2*P*Exp_msk1) {
953 if (aadj <= 0x7fffffff) {
954 if ((z = aadj) <= 0)
955 z = 1;
956 aadj = z;
910 aadj1 = dsign ? aadj : -aadj;
957 dval(&aadj1) = dsign ? aadj : -aadj;
911 }
958 }
912 word0(aadj1) += (2*P+1)*Exp_msk1 - y;
959 word0(&aadj1) += (2*P+1)*Exp_msk1 - y;
913 }
960 }
914 adj = aadj1 * ulp(dval(rv));
915 dval(rv) += adj;
961 dval(&adj) = dval(&aadj1) * ulp(&rv);
962 dval(&rv) += dval(&adj);
916#else
917#ifdef Sudden_Underflow
963#else
964#ifdef Sudden_Underflow
918 if ((word0(rv) & Exp_mask) <= P*Exp_msk1) {
919 dval(rv0) = dval(rv);
920 word0(rv) += P*Exp_msk1;
921 adj = aadj1 * ulp(dval(rv));
922 dval(rv) += adj;
965 if ((word0(&rv) & Exp_mask) <= P*Exp_msk1) {
966 dval(&rv0) = dval(&rv);
967 word0(&rv) += P*Exp_msk1;
968 dval(&adj) = dval(&aadj1) * ulp(&rv);
969 dval(&rv) += adj;
923#ifdef IBM
970#ifdef IBM
924 if ((word0(rv) & Exp_mask) < P*Exp_msk1)
971 if ((word0(&rv) & Exp_mask) < P*Exp_msk1)
925#else
972#else
926 if ((word0(rv) & Exp_mask) <= P*Exp_msk1)
973 if ((word0(&rv) & Exp_mask) <= P*Exp_msk1)
927#endif
928 {
974#endif
975 {
929 if (word0(rv0) == Tiny0
930 && word1(rv0) == Tiny1)
976 if (word0(&rv0) == Tiny0
977 && word1(&rv0) == Tiny1)
931 goto undfl;
978 goto undfl;
932 word0(rv) = Tiny0;
933 word1(rv) = Tiny1;
979 word0(&rv) = Tiny0;
980 word1(&rv) = Tiny1;
934 goto cont;
935 }
936 else
981 goto cont;
982 }
983 else
937 word0(rv) -= P*Exp_msk1;
984 word0(&rv) -= P*Exp_msk1;
938 }
939 else {
985 }
986 else {
940 adj = aadj1 * ulp(dval(rv));
941 dval(rv) += adj;
987 dval(&adj) = dval(&aadj1) * ulp(&rv);
988 dval(&rv) += adj;
942 }
943#else /*Sudden_Underflow*/
989 }
990#else /*Sudden_Underflow*/
944 /* Compute adj so that the IEEE rounding rules will
945 * correctly round rv + adj in some half-way cases.
946 * If rv * ulp(rv) is denormalized (i.e.,
991 /* Compute dval(&adj) so that the IEEE rounding rules will
992 * correctly round rv + dval(&adj) in some half-way cases.
993 * If rv * ulp(&rv) is denormalized (i.e.,
947 * y <= (P-1)*Exp_msk1), we must adjust aadj to avoid
948 * trouble from bits lost to denormalization;
949 * example: 1.2e-307 .
950 */
951 if (y <= (P-1)*Exp_msk1 && aadj > 1.) {
994 * y <= (P-1)*Exp_msk1), we must adjust aadj to avoid
995 * trouble from bits lost to denormalization;
996 * example: 1.2e-307 .
997 */
998 if (y <= (P-1)*Exp_msk1 && aadj > 1.) {
952 aadj1 = (double)(int)(aadj + 0.5);
999 dval(&aadj1) = (double)(int)(aadj + 0.5);
953 if (!dsign)
1000 if (!dsign)
954 aadj1 = -aadj1;
1001 dval(&aadj1) = -dval(&aadj1);
955 }
1002 }
956 adj = aadj1 * ulp(dval(rv));
957 dval(rv) += adj;
1003 dval(&adj) = dval(&aadj1) * ulp(&rv);
1004 dval(&rv) += adj;
958#endif /*Sudden_Underflow*/
959#endif /*Avoid_Underflow*/
960 }
1005#endif /*Sudden_Underflow*/
1006#endif /*Avoid_Underflow*/
1007 }
961 z = word0(rv) & Exp_mask;
1008 z = word0(&rv) & Exp_mask;
962#ifndef SET_INEXACT
963#ifdef Avoid_Underflow
964 if (!scale)
965#endif
966 if (y == z) {
967 /* Can we stop now? */
968 L = (Long)aadj;
969 aadj -= L;
970 /* The tolerances below are conservative. */
1009#ifndef SET_INEXACT
1010#ifdef Avoid_Underflow
1011 if (!scale)
1012#endif
1013 if (y == z) {
1014 /* Can we stop now? */
1015 L = (Long)aadj;
1016 aadj -= L;
1017 /* The tolerances below are conservative. */
971 if (dsign || word1(rv) || word0(rv) & Bndry_mask) {
1018 if (dsign || word1(&rv) || word0(&rv) & Bndry_mask) {
972 if (aadj < .4999999 || aadj > .5000001)
973 break;
974 }
975 else if (aadj < .4999999/FLT_RADIX)
976 break;
977 }
978#endif
979 cont:
980 Bfree(bb);
981 Bfree(bd);
982 Bfree(bs);
983 Bfree(delta);
984 }
1019 if (aadj < .4999999 || aadj > .5000001)
1020 break;
1021 }
1022 else if (aadj < .4999999/FLT_RADIX)
1023 break;
1024 }
1025#endif
1026 cont:
1027 Bfree(bb);
1028 Bfree(bd);
1029 Bfree(bs);
1030 Bfree(delta);
1031 }
1032 Bfree(bb);
1033 Bfree(bd);
1034 Bfree(bs);
1035 Bfree(bd0);
1036 Bfree(delta);
985#ifdef SET_INEXACT
986 if (inexact) {
987 if (!oldinexact) {
1037#ifdef SET_INEXACT
1038 if (inexact) {
1039 if (!oldinexact) {
988 word0(rv0) = Exp_1 + (70 << Exp_shift);
989 word1(rv0) = 0;
990 dval(rv0) += 1.;
1040 word0(&rv0) = Exp_1 + (70 << Exp_shift);
1041 word1(&rv0) = 0;
1042 dval(&rv0) += 1.;
991 }
992 }
993 else if (!oldinexact)
994 clear_inexact();
995#endif
996#ifdef Avoid_Underflow
997 if (scale) {
1043 }
1044 }
1045 else if (!oldinexact)
1046 clear_inexact();
1047#endif
1048#ifdef Avoid_Underflow
1049 if (scale) {
998 word0(rv0) = Exp_1 - 2*P*Exp_msk1;
999 word1(rv0) = 0;
1000 dval(rv) *= dval(rv0);
1050 word0(&rv0) = Exp_1 - 2*P*Exp_msk1;
1051 word1(&rv0) = 0;
1052 dval(&rv) *= dval(&rv0);
1001#ifndef NO_ERRNO
1002 /* try to avoid the bug of testing an 8087 register value */
1003#ifdef IEEE_Arith
1053#ifndef NO_ERRNO
1054 /* try to avoid the bug of testing an 8087 register value */
1055#ifdef IEEE_Arith
1004 if (!(word0(rv) & Exp_mask))
1056 if (!(word0(&rv) & Exp_mask))
1005#else
1057#else
1006 if (word0(rv) == 0 && word1(rv) == 0)
1058 if (word0(&rv) == 0 && word1(&rv) == 0)
1007#endif
1008 errno = ERANGE;
1009#endif
1010 }
1011#endif /* Avoid_Underflow */
1012#ifdef SET_INEXACT
1059#endif
1060 errno = ERANGE;
1061#endif
1062 }
1063#endif /* Avoid_Underflow */
1064#ifdef SET_INEXACT
1013 if (inexact && !(word0(rv) & Exp_mask)) {
1065 if (inexact && !(word0(&rv) & Exp_mask)) {
1014 /* set underflow bit */
1066 /* set underflow bit */
1015 dval(rv0) = 1e-300;
1016 dval(rv0) *= dval(rv0);
1067 dval(&rv0) = 1e-300;
1068 dval(&rv0) *= dval(&rv0);
1017 }
1018#endif
1069 }
1070#endif
1019 retfree:
1020 Bfree(bb);
1021 Bfree(bd);
1022 Bfree(bs);
1023 Bfree(bd0);
1024 Bfree(delta);
1025 ret:
1026 if (se)
1027 *se = (char *)s;
1071 ret:
1072 if (se)
1073 *se = (char *)s;
1028 return sign ? -dval(rv) : dval(rv);
1074 return sign ? -dval(&rv) : dval(&rv);
1029 }
1030
1075 }
1076