]>
git.proxmox.com Git - mirror_edk2.git/blob - StdLib/LibC/gdtoa/strtodg.c
1 /* $NetBSD: strtodg.c,v 1.5.14.1 2008/04/08 21:10:55 jdc Exp $ */
3 /****************************************************************
5 The author of this software is David M. Gay.
7 Copyright (C) 1998-2001 by Lucent Technologies
10 Permission to use, copy, modify, and distribute this software and
11 its documentation for any purpose and without fee is hereby
12 granted, provided that the above copyright notice appear in all
13 copies and that both that the copyright notice and this
14 permission notice and warranty disclaimer appear in supporting
15 documentation, and that the name of Lucent or any of its entities
16 not be used in advertising or publicity pertaining to
17 distribution of the software without specific, written prior
20 LUCENT DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE,
21 INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS.
22 IN NO EVENT SHALL LUCENT OR ANY OF ITS ENTITIES BE LIABLE FOR ANY
23 SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
24 WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER
25 IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION,
26 ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF
29 ****************************************************************/
31 /* Please send bug reports to David M. Gay (dmg at acm dot org,
32 * with " at " changed at "@" and " dot " changed to "."). */
33 #include <LibConfig.h>
42 // Disable warnings about assignment within conditional expressions.
43 #pragma warning ( disable : 4706 )
47 fivesbits
[] = { 0, 3, 5, 7, 10, 12, 14, 17, 19, 21,
48 24, 26, 28, 31, 33, 35, 38, 40, 42, 45,
57 increment(b
) Bigint
*b
;
72 if (*x
< (ULong
)0xffffffffL
) {
89 if (b
->wds
>= b
->maxwds
) {
104 decrement(b
) Bigint
*b
;
128 borrow
= (y
& 0x10000) >> 16;
130 } while(borrow
&& x
< xe
);
132 return STRTOG_Inexlo
;
137 all_on(b
, n
) CONST Bigint
*b
; int n
;
139 all_on(CONST Bigint
*b
, int n
)
145 xe
= x
+ ((unsigned int)n
>> kshift
);
147 if ((*x
++ & ALL_ON
) != ALL_ON
)
150 return ((*x
| (ALL_ON
<< n
)) & ALL_ON
) == ALL_ON
;
156 set_ones(b
, n
) Bigint
*b
; int n
;
158 set_ones(Bigint
*b
, int n
)
164 k
= (unsigned int)(n
+ ((1 << kshift
) - 1)) >> kshift
;
171 k
= (unsigned int)n
>> kshift
;
180 x
[-1] >>= ULbits
- n
;
187 (d
, fpi
, expt
, bits
, exact
, rd
, irv
)
188 double d
; CONST FPI
*fpi
; Long
*expt
; ULong
*bits
; int exact
, rd
, *irv
;
190 (double d
, CONST FPI
*fpi
, Long
*expt
, ULong
*bits
, int exact
, int rd
, int *irv
)
194 ULong carry
, inex
, lostbits
;
195 int bdif
, e
, j
, k
, k1
, nb
, rv
;
198 b
= d2b(d
, &e
, &bdif
);
199 bdif
-= nb
= fpi
->nbits
;
208 #ifndef IMPRECISE_INEXACT
225 default: /* round near */
234 if (b
->x
[(unsigned int)k
>>kshift
] & ((ULong
)1 << (k
& kmask
)))
238 /* "break" cases: round up 1 bit, then truncate; bdif > 0 */
243 if ( (lostbits
= any_on(b
, bdif
)) !=0)
244 inex
= STRTOG_Inexlo
;
247 inex
= STRTOG_Inexhi
;
249 if ( (j
= nb
& kmask
) !=0)
251 if (hi0bits(b
->x
[b
->wds
- 1]) != j
) {
253 lostbits
= b
->x
[0] & 1;
260 b
= lshift(b
, -bdif
);
264 if (k
> nb
|| fpi
->sudden_underflow
) {
266 *irv
= STRTOG_Underflow
| STRTOG_Inexlo
;
270 if (k1
> 0 && !lostbits
)
271 lostbits
= any_on(b
, k1
);
272 if (!lostbits
&& !exact
)
275 carry
= b
->x
[(unsigned int)k1
>>kshift
] &
276 (ULong
)(1 << ((unsigned int)k1
& kmask
));
278 *irv
= STRTOG_Denormal
;
281 inex
= STRTOG_Inexhi
| STRTOG_Underflow
;
284 inex
= STRTOG_Inexlo
| STRTOG_Underflow
;
287 else if (e
> fpi
->emax
) {
289 *irv
= STRTOG_Infinite
| STRTOG_Overflow
| STRTOG_Inexhi
;
296 copybits(bits
, nb
, b
);
307 mantbits(d
) double d
;
314 L
= word1(d
) << 16 | word1(d
) >> 16;
317 if ( (L
= word1(d
)) !=0)
319 return P
- lo0bits(&L
);
321 L
= word0(d
) << 16 | word0(d
) >> 16 | Exp_msk11
;
323 L
= word0(d
) | Exp_msk1
;
325 return P
- 32 - lo0bits(&L
);
332 (s00
, se
, fpi
, expt
, bits
)
333 CONST
char *s00
; char **se
; CONST FPI
*fpi
; Long
*expt
; ULong
*bits
;
335 (CONST
char *s00
, char **se
, CONST FPI
*fpi
, Long
*expt
, ULong
*bits
)
338 int abe
, abits
, asub
;
339 int bb0
, bb2
, bb5
, bbe
, bd2
, bd5
, bbbits
, bs2
, c
, decpt
, denorm
;
340 int dsign
, e
, e1
, e2
, emin
, esign
, finished
, i
, inex
, irv
;
341 int j
, k
, nbits
, nd
, nd0
, nf
, nz
, nz0
, rd
, rvbits
, rve
, rve1
, sign
;
342 int sudden_underflow
= 0; /* pacify gcc */
343 CONST
char *s
, *s0
, *s1
;
344 double adj
, adj0
, rv
, tol
;
347 Bigint
*ab
, *bb
, *bb1
, *bd
, *bd0
, *bs
, *delta
, *rvb
, *rvb0
;
349 e2
= 0; /* XXX gcc */
352 denorm
= sign
= nz0
= nz
= 0;
356 for(s
= s00
;;s
++) switch(*s
) {
366 irv
= STRTOG_NoNumber
;
385 irv
= gethex(&s
, fpi
, expt
, &rvb
, sign
);
386 if (irv
== STRTOG_NoNumber
) {
398 sudden_underflow
= fpi
->sudden_underflow
;
401 for(decpt
= nd
= nf
= 0; (c
= *s
) >= '0' && c
<= '9'; nd
++, s
++)
408 if (c
== *localeconv()->decimal_point
)
416 for(; c
== '0'; c
= *++s
)
418 if (c
> '0' && c
<= '9') {
426 for(; c
>= '0' && c
<= '9'; c
= *++s
) {
431 for(i
= 1; i
< nz
; i
++)
434 else if (nd
<= DBL_DIG
+ 1)
438 else if (nd
<= DBL_DIG
+ 1)
446 if (c
== 'e' || c
== 'E') {
447 if (!nd
&& !nz
&& !nz0
) {
448 irv
= STRTOG_NoNumber
;
461 if (c
>= '0' && c
<= '9') {
464 if (c
> '0' && c
<= '9') {
467 while((c
= *++s
) >= '0' && c
<= '9')
469 if (s
- s1
> 8 || L
> 19999)
470 /* Avoid confusion from exponents
471 * so large that e might overflow.
473 e
= 19999; /* safe for 16 bit ints */
488 /* Check for Nan and Infinity */
493 if (match(&s
,"nf")) {
495 if (!match(&s
,"inity"))
497 irv
= STRTOG_Infinite
;
503 if (match(&s
, "an")) {
505 *expt
= fpi
->emax
+ 1;
508 irv
= hexnan(&s
, fpi
, bits
);
513 #endif /* INFNAN_CHECK */
514 irv
= STRTOG_NoNumber
;
523 switch(fpi
->rounding
& 3) {
534 /* Now we have nd0 digits, starting at s0, followed by a
535 * decimal point, followed by nd-nd0 digits. The number we're
536 * after is the integer represented by those digits times
541 k
= nd
< DBL_DIG
+ 1 ? nd
: DBL_DIG
+ 1;
542 dval(rv
) = (double)y
;
544 dval(rv
) = tens
[k
- 9] * dval(rv
) + z
;
546 if (nbits
<= P
&& nd
<= DBL_DIG
) {
548 if (rvOK(dval(rv
), fpi
, expt
, bits
, 1, rd
, &irv
))
556 i
= fivesbits
[e
] + mantbits(dval(rv
)) <= P
;
557 /* rv = */ rounded_product(dval(rv
), tens
[e
]);
558 if (rvOK(dval(rv
), fpi
, expt
, bits
, i
, rd
, &irv
))
565 if (e
<= Ten_pmax
+ i
) {
566 /* A fancier test would sometimes let us do
567 * this for larger i values.
573 /* VAX exponent range is so narrow we must
574 * worry about overflow here...
577 dval(adj
) = dval(rv
);
578 word0(adj
) -= P
*Exp_msk1
;
579 /* adj = */ rounded_product(dval(adj
), tens
[e2
]);
580 if ((word0(adj
) & Exp_mask
)
581 > Exp_msk1
*(DBL_MAX_EXP
+Bias
-1-P
))
583 word0(adj
) += P
*Exp_msk1
;
584 dval(rv
) = dval(adj
);
586 /* rv = */ rounded_product(dval(rv
), tens
[e2
]);
588 if (rvOK(dval(rv
), fpi
, expt
, bits
, 0, rd
, &irv
))
593 #ifndef Inaccurate_Divide
594 else if (e
>= -Ten_pmax
) {
595 /* rv = */ rounded_quotient(dval(rv
), tens
[-e
]);
596 if (rvOK(dval(rv
), fpi
, expt
, bits
, 0, rd
, &irv
))
605 /* Get starting approximation = rv * 10**e1 */
609 if ( (i
= e1
& 15) !=0)
612 e1
= (unsigned int)e1
>> 4;
613 while(e1
>= (1 << (n_bigtens
-1))) {
614 e2
+= ((word0(rv
) & Exp_mask
)
615 >> Exp_shift1
) - Bias
;
616 word0(rv
) &= ~Exp_mask
;
617 word0(rv
) |= Bias
<< Exp_shift1
;
618 dval(rv
) *= bigtens
[n_bigtens
-1];
619 e1
-= 1 << (n_bigtens
-1);
621 e2
+= ((word0(rv
) & Exp_mask
) >> Exp_shift1
) - Bias
;
622 word0(rv
) &= ~Exp_mask
;
623 word0(rv
) |= Bias
<< Exp_shift1
;
624 for(j
= 0; e1
> 0; j
++, e1
= (unsigned int)e1
>> 1)
626 dval(rv
) *= bigtens
[j
];
631 if ( (i
= e1
& 15) !=0)
634 e1
= (unsigned int)e1
>> 4;
635 while(e1
>= (1 << (n_bigtens
-1))) {
636 e2
+= ((word0(rv
) & Exp_mask
)
637 >> Exp_shift1
) - Bias
;
638 word0(rv
) &= ~Exp_mask
;
639 word0(rv
) |= Bias
<< Exp_shift1
;
640 dval(rv
) *= tinytens
[n_bigtens
-1];
641 e1
-= 1 << (n_bigtens
-1);
643 e2
+= ((word0(rv
) & Exp_mask
) >> Exp_shift1
) - Bias
;
644 word0(rv
) &= ~Exp_mask
;
645 word0(rv
) |= Bias
<< Exp_shift1
;
646 for(j
= 0; e1
> 0; j
++, e1
= (unsigned int)e1
>> 1)
648 dval(rv
) *= tinytens
[j
];
652 /* e2 is a correction to the (base 2) exponent of the return
653 * value, reflecting adjustments above to avoid overflow in the
654 * native arithmetic. For native IBM (base 16) arithmetic, we
655 * must multiply e2 by 4 to change from base 16 to 2.
659 rvb
= d2b(dval(rv
), &rve
, &rvbits
); /* rv = rvb * 2^rve */
661 return STRTOG_NoMemory
;
663 if ((j
= rvbits
- nbits
) > 0) {
668 bb0
= 0; /* trailing zero bits in rvb */
669 e2
= rve
+ rvbits
- nbits
;
670 if (e2
> fpi
->emax
+ 1)
672 rve1
= rve
+ rvbits
- nbits
;
673 if (e2
< (emin
= fpi
->emin
)) {
677 rvb
= lshift(rvb
, j
);
688 irv
= STRTOG_Underflow
| STRTOG_Inexlo
;
691 rvb
->x
[0] = rvb
->wds
= rvbits
= 1;
697 if (sudden_underflow
&& e2
+ 1 < emin
)
701 /* Now the hard part -- adjusting rv to the correct value.*/
703 /* Put digits into bd: true value = bd * 10^e */
705 bd0
= s2b(s0
, nd0
, nd
, y
);
710 return STRTOG_NoMemory
;
714 return STRTOG_NoMemory
;
716 bbbits
= rvbits
- bb0
;
720 return STRTOG_NoMemory
;
735 j
= nbits
+ 1 - bbbits
;
736 i
= bbe
+ bbbits
- nbits
;
737 if (i
< emin
) /* denormal */
741 i
= bb2
< bd2
? bb2
: bd2
;
750 bs
= pow5mult(bs
, bb5
);
752 return STRTOG_NoMemory
;
755 return STRTOG_NoMemory
;
761 bb
= lshift(bb
, bb2
);
763 return STRTOG_NoMemory
;
768 bd
= pow5mult(bd
, bd5
);
770 return STRTOG_NoMemory
;
773 bd
= lshift(bd
, bd2
);
775 return STRTOG_NoMemory
;
778 bs
= lshift(bs
, bs2
);
780 return STRTOG_NoMemory
;
783 inex
= STRTOG_Inexhi
;
784 delta
= diff(bb
, bd
);
786 return STRTOG_NoMemory
;
787 if (delta
->wds
<= 1 && !delta
->x
[0])
790 delta
->sign
= finished
= 0;
795 if ( (finished
= dsign
^ (rd
&1)) !=0) {
797 irv
|= STRTOG_Inexhi
;
800 irv
|= STRTOG_Inexlo
;
803 for(i
= 0, j
= nbits
; j
>= ULbits
;
805 if (rvb
->x
[i
] & ALL_ON
)
808 if (j
> 1 && lo0bits(rvb
->x
+ i
) < j
- 1)
811 rvb
= set_ones(rvb
, rvbits
= nbits
);
813 return STRTOG_NoMemory
;
816 irv
|= dsign
? STRTOG_Inexlo
: STRTOG_Inexhi
;
820 /* Error is less than half an ulp -- check for
821 * special case of mantissa a power of two.
824 ? STRTOG_Normal
| STRTOG_Inexlo
825 : STRTOG_Normal
| STRTOG_Inexhi
;
826 if (dsign
|| bbbits
> 1 || denorm
|| rve1
== emin
)
828 delta
= lshift(delta
,1);
830 return STRTOG_NoMemory
;
831 if (cmp(delta
, bs
) > 0) {
832 irv
= STRTOG_Normal
| STRTOG_Inexlo
;
838 /* exactly half-way between */
840 if (denorm
&& all_on(rvb
, rvbits
)) {
841 /*boundary case -- increment exponent*/
844 rve
= emin
+ nbits
- (rvbits
= 1);
845 irv
= STRTOG_Normal
| STRTOG_Inexhi
;
849 irv
= STRTOG_Normal
| STRTOG_Inexlo
;
851 else if (bbbits
== 1) {
854 /* boundary case -- decrement exponent */
856 irv
= STRTOG_Normal
| STRTOG_Inexhi
;
857 if (rvb
->wds
== 1 && rvb
->x
[0] == 1)
858 sudden_underflow
= 1;
862 rvb
= set_ones(rvb
, rvbits
= nbits
);
864 return STRTOG_NoMemory
;
868 irv
= STRTOG_Normal
| STRTOG_Inexhi
;
869 if ((bbbits
< nbits
&& !denorm
) || !(rvb
->x
[0] & 1))
872 rvb
= increment(rvb
);
874 return STRTOG_NoMemory
;
875 if ( (j
= rvbits
& kmask
) !=0)
877 if (hi0bits(rvb
->x
[(unsigned int)(rvb
->wds
- 1)
881 irv
= STRTOG_Normal
| STRTOG_Inexhi
;
887 irv
= STRTOG_Normal
| STRTOG_Inexlo
;
891 if ((dval(adj
) = ratio(delta
, bs
)) <= 2.) {
893 inex
= STRTOG_Inexlo
;
896 inex
= STRTOG_Inexhi
;
898 else if (denorm
&& bbbits
<= 1) {
902 irv
= STRTOG_Underflow
| STRTOG_Inexlo
;
905 adj0
= dval(adj
) = 1.;
908 adj0
= dval(adj
) *= 0.5;
911 inex
= STRTOG_Inexlo
;
913 if (dval(adj
) < 2147483647.) {
922 if (asub
&& adj0
> 0.)
926 if (!asub
&& adj0
> 0.) {
929 inex
= STRTOG_Inexact
- inex
;
932 dval(adj
) = (double)L
;
937 /* adj *= ulp(dval(rv)); */
938 /* if (asub) rv -= adj; else rv += adj; */
940 if (!denorm
&& rvbits
< nbits
) {
941 rvb
= lshift(rvb
, j
= nbits
- rvbits
);
943 return STRTOG_NoMemory
;
947 ab
= d2b(dval(adj
), &abe
, &abits
);
949 return STRTOG_NoMemory
;
953 ab
= lshift(ab
, abe
);
957 j
= hi0bits(rvb
->x
[rvb
->wds
-1]);
960 return STRTOG_NoMemory
;
964 else if (rvb
->wds
<= k
965 || hi0bits( rvb
->x
[k
]) >
966 hi0bits(rvb0
->x
[k
])) {
967 /* unlikely; can only have lost 1 high bit */
973 rvb
= lshift(rvb
, 1);
975 return STRTOG_NoMemory
;
985 return STRTOG_NoMemory
;
988 || hi0bits(rvb
->x
[k
]) < hi0bits(rvb0
->x
[k
])) {
990 if (++rvbits
== nbits
)
1008 /* Can we stop now? */
1009 tol
= dval(adj
) * 5e-16; /* > max rel error */
1010 dval(adj
) = adj0
- .5;
1011 if (dval(adj
) < -tol
) {
1017 else if (dval(adj
) > tol
&& adj0
< 1. - tol
) {
1022 bb0
= denorm
? 0 : trailz(rvb
);
1028 if (!denorm
&& (j
= nbits
- rvbits
)) {
1030 rvb
= lshift(rvb
, j
);
1041 if (rve
> fpi
->emax
) {
1044 irv
= STRTOG_Infinite
| STRTOG_Overflow
| STRTOG_Inexhi
;
1051 *expt
= fpi
->emax
+ 1;
1055 if (sudden_underflow
) {
1057 irv
= STRTOG_Underflow
| STRTOG_Inexlo
;
1060 irv
= (irv
& ~STRTOG_Retmask
) |
1061 (rvb
->wds
> 0 ? STRTOG_Denormal
: STRTOG_Zero
);
1062 if (irv
& STRTOG_Inexact
)
1063 irv
|= STRTOG_Underflow
;
1071 copybits(bits
, nbits
, rvb
);