shithub: riscv

Info  •  Files  •  Log  •  Branches

ape: bring strtod() in line with plan9's libc version

--- a/sys/src/ape/lib/fmt/strtod.c

+++ b/sys/src/ape/lib/fmt/strtod.c

@@ -27,17 +27,6 @@

 #define ulong _fmtulong

 typedef unsigned long ulong;

-static ulong

-umuldiv(ulong a, ulong b, ulong c)

-{

-	double d;

-

-	d = ((double)a * (double)b) / (double)c;

-	if(d >= 4294967295.)

-		d = 4294967295.;

-	return (ulong)d;

-}

-

 /*

  * This routine will convert to arbitrary precision

  * floating point entirely in multi-precision fixed.

@@ -55,6 +44,7 @@

 {

 	Nbits	= 28,				/* bits safely represented in a ulong */

 	Nmant	= 53,				/* bits of precision required */

+	Bias	= 1022,

 	Prec	= (Nmant+Nbits+1)/Nbits,	/* words of Nbits each to represent mantissa */

 	Sigbit	= 1<<(Prec*Nbits-Nmant),	/* first significant bit of Prec-th word */

 	Ndig	= 1500,

@@ -62,10 +52,6 @@

 	Half	= (ulong)(One>>1),

 	Maxe	= 310,

-	Fsign	= 1<<0,		/* found - */

-	Fesign	= 1<<1,		/* found e- */

-	Fdpoint	= 1<<2,		/* found . */

-

 	S0	= 0,		/* _		_S0	+S1	#S2	.S3 */

 	S1,			/* _+		#S2	.S3 */

 	S2,			/* _+#		#S2	.S4	eS5 */

@@ -74,6 +60,10 @@

 	S5,			/* _+#.#e	+S6	#S7 */

 	S6,			/* _+#.#e+	#S7 */

 	S7,			/* _+#.#e+#	#S7 */

+

+	Fsign	= 1<<0,		/* found - */

+	Fesign	= 1<<1,		/* found e- */

+	Fdpoint	= 1<<2,		/* found . */

 };

 static	int	xcmp(char*, char*);

@@ -81,6 +71,8 @@

 static	void	frnorm(ulong*);

 static	void	divascii(char*, int*, int*, int*);

 static	void	mulascii(char*, int*, int*, int*);

+static	void	divby(char*, int*, int);

+static ulong	umuldiv(ulong, ulong, ulong);

 typedef	struct	Tab	Tab;

 struct	Tab

@@ -93,8 +85,8 @@

 double

 fmtstrtod(const char *as, char **aas)

 {

-	int na, ex, dp, bp, c, i, flag, state;

-	ulong low[Prec], hig[Prec], mid[Prec];

+	int na, ona, ex, dp, bp, c, i, flag, state;

+	ulong low[Prec], hig[Prec], mid[Prec], num, den;

 	double d;

 	char *s, a[Ndig];

@@ -225,7 +217,7 @@

 	if(flag & Fesign)

 		ex = -ex;

 	dp += ex;

-	if(dp < -Maxe){

+	if(dp < -Maxe-Nmant/3){ /* actually -Nmant*log(2)/log(10), but Nmant/3 close enough */

 		errno = ERANGE;

 		goto ret0;	/* underflow by exp */

 	} else

@@ -241,18 +233,33 @@

 		divascii(a, &na, &dp, &bp);

 	while(dp < 0 || a[0] < '5')

 		mulascii(a, &na, &dp, &bp);

+	a[na] = 0;

+	/*

+	 * very small numbers are represented using

+	 * bp = -Bias+1.  adjust accordingly.

+	 */

+	if(bp < -Bias+1){

+		ona = na;

+		divby(a, &na, -bp-Bias+1);

+		if(na < ona){

+			memmove(a+ona-na, a, na);

+			memset(a, '0', ona-na);

+			na = ona;

+		}

+		a[na] = 0;

+		bp = -Bias+1;

+	}

+

 	/* close approx by naive conversion */

-	mid[0] = 0;

-	mid[1] = 1;

-	for(i=0; c=a[i]; i++) {

-		mid[0] = mid[0]*10 + (c-'0');

-		mid[1] = mid[1]*10;

-		if(i >= 8)

-			break;

+	num = 0;

+	den = 1;

+	for(i=0; i<9 && (c=a[i]); i++) {

+		num = num*10 + (c-'0');

+		den *= 10;

 	}

-	low[0] = umuldiv(mid[0], One, mid[1]);

-	hig[0] = umuldiv(mid[0]+1, One, mid[1]);

+	low[0] = umuldiv(num, One, den);

+	hig[0] = umuldiv(num+1, One, den);

 	for(i=1; i<Prec; i++) {

 		low[i] = 0;

 		hig[i] = One-1;

@@ -386,7 +393,7 @@

 }

 static void

-divby(char *a, int *na, int b)

+_divby(char *a, int *na, int b)

 {

 	int n, c;

 	char *p;

@@ -428,6 +435,18 @@

 	*p = 0;

 }

+static void

+divby(char *a, int *na, int b)

+{

+	while(b > 9){

+		_divby(a, na, 9);

+		a[*na] = 0;

+		b -= 9;

+	}

+	if(b > 0)

+		_divby(a, na, b);

+}

+

 static	Tab	tab1[] =

 {

 	 1,  0, "",

@@ -535,4 +554,10 @@

 			return 1;

 	}

 	return 0;

+}

+

+static ulong

+umuldiv(ulong a, ulong b, ulong c)

+{

+	return ((uvlong)a * (uvlong)b) / c;

 }