From 07e77936d535e58b0aeb4f2a11400c1050556739 Mon Sep 17 00:00:00 2001 From: Nicholas Noll Date: Sun, 5 Dec 2021 16:16:21 -0800 Subject: Feat: added math library Used Plan9's libc as starting point. This cleans up dangling references due to loss of libc. --- src/base/string/raw/asfloat.c | 474 ++++++++++++++++++++++++++++++++++++++++++ src/base/string/raw/atof.c | 8 + src/base/string/raw/atoi.c | 4 +- 3 files changed, 484 insertions(+), 2 deletions(-) create mode 100644 src/base/string/raw/asfloat.c create mode 100644 src/base/string/raw/atof.c (limited to 'src/base/string/raw') diff --git a/src/base/string/raw/asfloat.c b/src/base/string/raw/asfloat.c new file mode 100644 index 0000000..7d35626 --- /dev/null +++ b/src/base/string/raw/asfloat.c @@ -0,0 +1,474 @@ +#include +#include + +#define Nbits 28 +#define Nmant 53 +#define Prec ((Nmant+Nbits+1)/Nbits) + +#define Sigbit (1<<(Prec*Nbits-Nmant)) /* first significant bit of Prec-th word */ +#define Ndig 1500 +#define One (ulong)(1<>1) +#define Maxe 310 + +#define Fsign (1<<0) /* found - */ +#define Fesign (1<<1) /* found e- */ +#define Fdpoint (1<<2) /* found . */ + +#define S0 0 /* _ _S0 +S1 #S2 .S3 */ +#define S1 1 /* _+ #S2 .S3 */ +#define S2 2 /* _+# #S2 .S4 eS5 */ +#define S3 3 /* _+. #S4 */ +#define S4 4 /* _+#.# #S4 eS5 */ +#define S5 5 /* _+#.#e +S6 #S7 */ +#define S6 6 /* _+#.#e+ #S7 */ +#define S7 7 /* _+#.#e+# #S7 */ + +typedef struct Tab Tab; +struct Tab +{ + int bp; + int siz; + char *cmp; +}; + +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; +} + +static void +frnorm(ulong *f) +{ + int i, c; + + c = 0; + for(i=Prec-1; i>0; i--) { + f[i] += c; + c = f[i] >> Nbits; + f[i] &= One-1; + } + f[0] += c; +} + +static int +fpcmp(char *a, ulong* f) +{ + ulong tf[Prec]; + int i, d, c; + + for(i=0; i> Nbits) + '0'; + tf[0] &= One-1; + + /* compare next digit */ + c = *a; + if(c == 0) { + if('0' < d) + return -1; + if(tf[0] != 0) + goto cont; + for(i=1; i d) + return +1; + if(c < d) + return -1; + a++; + cont:; +} +} + +static void +divby(char *a, int *na, int b) +{ + int n, c; + char *p; + + p = a; + n = 0; + while(n>>b == 0){ + c = *a++; + if(c == 0) { + while(n) { + c = n*10; + if(c>>b) + break; + n = c; + } + goto xx; + } + n = n*10 + c-'0'; + (*na)--; + } + for(;;){ + c = n>>b; + n -= c<>b; + n -= c<= (int)(arrlen(tab1))) + d = (int)(arrlen(tab1))-1; + t = tab1 + d; + b = t->bp; + if(mem·compare(a, t->siz, t->cmp) > 0) + d--; + *dp -= d; + *bp += b; + divby(a, na, b); +} + +static void +mulby(char *a, char *p, char *q, int b) +{ + int n, c; + + n = 0; + *p = 0; + for(;;) { + q--; + if(q < a) + break; + c = *q - '0'; + c = (c<= (int)(arrlen(tab2))) + d = (int)(arrlen(tab2))-1; + t = tab2 + d; + b = t->bp; + if(mem·compare(a, t->siz, t->cmp) < 0) + d--; + p = a + *na; + *bp -= b; + *dp += d; + *na += d; + mulby(a, p+d, p, b); +} + +static int +cmp(char *a, char *b) +{ + int c1, c2; + + while((c1 = *b++) != '\0') { + c2 = *a++; + if(utf8·isupper(c2)) + c2 = utf8·tolower(c2); + if(c1 != c2) + return 1; + } + return 0; +} + +double +str·asfloat(char *as, char **aas) +{ + int na, ex, dp, bp, c, i, flag, state; + ulong low[Prec], hig[Prec], mid[Prec]; + double d; + char *s, a[Ndig]; + + flag = 0; /* Fsign, Fesign, Fdpoint */ + na = 0; /* number of digits of a[] */ + dp = 0; /* na of decimal point */ + ex = 0; /* exonent */ + + state = S0; + for(s=as;;s++){ + c = *s; + if('0' <= c && c <= '9'){ + switch(state){ + case S0: case S1: case S2: + state = S2; + break; + case S3: case S4: + state = S4; + break; + case S5: case S6: case S7: + state = S7; + ex = ex*10 + (c-'0'); + continue; + } + + if(na == 0 && c == '0'){ + dp--; + continue; + } + if(na < Ndig-50) + a[na++] = c; + continue; + } + switch(c){ + case '\t': case '\n': case '\v': case '\f': case '\r': case ' ': + if(state == S0) + continue; + break; + case '-': + if(state == S0) + flag |= Fsign; + else + flag |= Fesign; + case '+': + if(state == S0) + state = S1; + else + if(state == S5) + state = S6; + else + break; /* syntax */ + continue; + case '.': + flag |= Fdpoint; + dp = na; + if(state == S0 || state == S1){ + state = S3; + continue; + } + if(state == S2){ + state = S4; + continue; + } + break; + case 'e': case 'E': + if(state == S2 || state == S4){ + state = S5; + continue; + } + break; + } + break; + } + + /* clean up return char-pointer */ + switch(state) { + case S0: + if(cmp(s, "nan") == 0){ + if(aas != nil) + *aas = s+3; + goto retnan; + } + case S1: + if(cmp(s, "infinity") == 0){ + if(aas != nil) + *aas = s+8; + goto retinf; + } + if(cmp(s, "inf") == 0){ + if(aas != nil) + *aas = s+3; + goto retinf; + } + case S3: + if(aas != nil) + *aas = as; + goto ret0; /* no digits found */ + case S6: + s--; /* back over +- */ + case S5: + s--; /* back over e */ + break; + } + if(aas != nil) + *aas = s; + + if(flag & Fdpoint) + while(na > 0 && a[na-1] == '0') + na--; + if(na == 0) + goto ret0; /* zero */ + a[na] = 0; + if(!(flag & Fdpoint)) + dp = na; + if(flag & Fesign) + ex = -ex; + dp += ex; + if(dp < -Maxe){ + /* errno = ERANGE; */ + goto ret0; /* underflow by exp */ + } else + if(dp > +Maxe) + goto retinf; /* overflow by exp */ + + /* + * normalize the decimal ascii number + * to range .[5-9][0-9]* e0 + */ + bp = 0; /* binary exponent */ + while(dp > 0) + divascii(a, &na, &dp, &bp); + while(dp < 0 || a[0] < '5') + mulascii(a, &na, &dp, &bp); + + /* close approx by naive conversion */ + mid[0] = 0; + mid[1] = 1; + for(i=0; (c=a[i]) != '\0'; i++) { + mid[0] = mid[0]*10 + (c-'0'); + mid[1] = mid[1]*10; + if(i >= 8) + break; + } + low[0] = umuldiv(mid[0], One, mid[1]); + hig[0] = umuldiv(mid[0]+1, One, mid[1]); + for(i=1; i>= 1; + } + frnorm(mid); + + /* compare */ + c = fpcmp(a, mid); + if(c > 0) { + c = 1; + for(i=0; i= Sigbit/2) { + mid[Prec-1] += Sigbit; + frnorm(mid); + } + goto out; + +ret0: + return 0; + +retnan: + return math·NaN(); + +retinf: + if(flag & Fsign) + return math·Inf(-1); + return math·Inf(+1); + +out: + d = 0; + for(i=0; i +#include + +double +str·atof(char *s) +{ + return str·asfloat(s, nil); +} diff --git a/src/base/string/raw/atoi.c b/src/base/string/raw/atoi.c index 8084e3e..e62ec62 100644 --- a/src/base/string/raw/atoi.c +++ b/src/base/string/raw/atoi.c @@ -1,7 +1,7 @@ -int +long str·atoi(char *s) { - int n = 0; + long n = 0; while(*s) n = 10*n + (*s++ - '0'); -- cgit v1.2.1