#define FDIGIT 30 #define FDEFLT 6 #define NSIGNIF 17 static uvlong uvnan = ((uvlong)0x7FF00000<<32)|0x00000001; static uvlong uvinf = ((uvlong)0x7FF00000<<32)|0x00000000; static uvlong uvneginf = ((uvlong)0xFFF00000<<32)|0x00000000; static char *special[] = { "NaN", "NaN", "+Inf", "+Inf", "-Inf", "-Inf" }; static int isNaN(double val) { union{ uvlong i; double f; }x; x.f = val; return (x.i&uvinf) == uvinf && (x.i&~uvneginf) != 0; } static double NaN(void) { union{ uvlong i; double f; }x; x.i = uvnan; return x.f; } static int isInf(double val, int sign) { union{ uvlong i; double f; }x; x.f = val; if(sign == 0) return x.i == uvinf || x.i == uvneginf; else if(sign == 1) return x.i == uvinf; else return x.i == uvneginf; } static double pows10[] = { 1e0, 1e1, 1e2, 1e3, 1e4, 1e5, 1e6, 1e7, 1e8, 1e9, 1e10, 1e11, 1e12, 1e13, 1e14, 1e15, 1e16, 1e17, 1e18, 1e19, 1e20, 1e21, 1e22, 1e23, 1e24, 1e25, 1e26, 1e27, 1e28, 1e29, 1e30, 1e31, 1e32, 1e33, 1e34, 1e35, 1e36, 1e37, 1e38, 1e39, 1e40, 1e41, 1e42, 1e43, 1e44, 1e45, 1e46, 1e47, 1e48, 1e49, 1e50, 1e51, 1e52, 1e53, 1e54, 1e55, 1e56, 1e57, 1e58, 1e59, 1e60, 1e61, 1e62, 1e63, 1e64, 1e65, 1e66, 1e67, 1e68, 1e69, 1e70, 1e71, 1e72, 1e73, 1e74, 1e75, 1e76, 1e77, 1e78, 1e79, 1e80, 1e81, 1e82, 1e83, 1e84, 1e85, 1e86, 1e87, 1e88, 1e89, 1e90, 1e91, 1e92, 1e93, 1e94, 1e95, 1e96, 1e97, 1e98, 1e99, 1e100, 1e101, 1e102, 1e103, 1e104, 1e105, 1e106, 1e107, 1e108, 1e109, 1e110, 1e111, 1e112, 1e113, 1e114, 1e115, 1e116, 1e117, 1e118, 1e119, 1e120, 1e121, 1e122, 1e123, 1e124, 1e125, 1e126, 1e127, 1e128, 1e129, 1e130, 1e131, 1e132, 1e133, 1e134, 1e135, 1e136, 1e137, 1e138, 1e139, 1e140, 1e141, 1e142, 1e143, 1e144, 1e145, 1e146, 1e147, 1e148, 1e149, 1e150, 1e151, 1e152, 1e153, 1e154, 1e155, 1e156, 1e157, 1e158, 1e159, }; static double fpow10(int n) { double d; int neg; neg = 0; if(n < 0){ neg = 1; n = -n; } if(n NSIGNIF) return 0; for(b = a+n-1; b >= a; b--){ c = *b + 1; if(c <= '9'){ *b = c; return 0; } *b = '0'; } /* * need to overflow adding digit. * shift number down and insert 1 at beginning. * decimal is known to be 0s or we wouldn't * have gotten this far. (e.g., 99999+1 => 00000) */ a[0] = '1'; return 1; } static int sub1(char *a, int n) { int c; char *b; if(n < 0 || n > NSIGNIF) return 0; for(b = a+n-1; b >= a; b--){ c = *b - 1; if(c >= '0'){ if(c == '0' && b == a){ /* * just zeroed the top digit; shift everyone up. * decimal is known to be 9s or we wouldn't * have gotten this far. (e.g., 10000-1 => 09999) */ *b = '9'; return 1; } *b = c; return 0; } *b = '9'; } /* * can't get here. the number a is always normalized * so that it has a nonzero first digit. */ abort(); } // ----------------------------------------------------------------------- // strtod #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(memcmp(a, t->cmp, t->siz) > 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(memcmp(a, t->cmp, t->siz) < 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(isupper(c2)) c2 = tolower(c2); if(c1 != c2) return 1; } return 0; } double fmtstrtod(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 NaN(); retinf: /* Unix strtod requires these. Plan 9 would return Inf(0) or Inf(-1). */ errno = ERANGE; if(flag & Fsign) return -HUGE_VAL; return HUGE_VAL; out: d = 0; for(i=0; i 0) *p++ = se[--i]; *p++ = '\0'; } /* * compute decimal integer m, exp such that: * f = m*10^exp * m is as short as possible with losing exactness * assumes special cases (NaN, +Inf, -Inf) have been handled. */ static void dtoa(double f, char *s, int *exp, int *neg, int *len) { int c, d, e2, e, ee, i, ndigit, oerrno; char buf[NSIGNIF+10]; double g; oerrno = errno; *neg = 0; if(f < 0){ f = -f; *neg = 1; } if(f == 0){ *exp = 0; s[0] = '0'; s[1] = 0; *len = 1; return; } frexp(f, &e2); e = (int)(e2 * .301029995664); g = f * fpow10(-e); while(g < 1) { e--; g = f * fpow10(-e); } while(g >= 10){ e++; g = f * fpow10(-e); } /* convert nsignif digits as a first approximation */ for(i=0; i g) { if(add1(s, NSIGNIF)){ /* gained a digit */ e--; fmtexp(s+NSIGNIF, e, 0); } continue; } if(f < g){ if(sub1(s, NSIGNIF)){ /* lost a digit */ e++; fmtexp(s+NSIGNIF, e, 0); } continue; } break; } /* * bump last few digits down to 0 as we can. */ for(i=NSIGNIF-1; i>=NSIGNIF-3; i--){ c = s[i]; if(c != '0'){ s[i] = '0'; g=fmtstrtod(s, nil); if(g != f){ s[i] = c; break; } } } /* * remove trailing zeros. */ ndigit = NSIGNIF; while(ndigit > 1 && s[ndigit-1] == '0'){ e++; --ndigit; } s[ndigit] = 0; *exp = e; *len = ndigit; errno = oerrno; } static int fmtfloat(fmt·State *io) { char buf[NSIGNIF+10], *dot, *digits, *p, *end, suf[10], *cur; double val; int c, verb, ndot, e, exp, f, ndigits, neg, newndigits; int npad, pt, prec, realverb, sign, nsuf, ucase, n, z1, z2; if(io->flag&fmt·Long) val = va_arg(io->args, long double); else val = va_arg(io->args, double); /* extract formatting flags */ f = io->flag; io->flag = 0; prec = FDEFLT; if(f & fmt·Prec) prec = io->prec; verb = io->verb; ucase = 0; switch(verb) { case 'A': case 'E': case 'F': case 'G': verb += 'a'-'A'; ucase = 1; break; } /* pick off special numbers. */ if(isNaN(val)) { end = special[0+ucase]; special: io->flag = f & (fmt·Width|fmt·Left); return copy(io, end, strlen(end), strlen(end)); } if(isInf(val, 1)) { end = special[2+ucase]; goto special; } if(isInf(val, -1)) { end = special[4+ucase]; goto special; } /* get exact representation. */ digits = buf; dtoa(val, digits, &exp, &neg, &ndigits); /* get locale's decimal point. */ dot = io->decimal; if(dot == nil) dot = "."; ndot = utf8·len(dot); /* * now the formatting fun begins. * compute parameters for actual fmt: * * pad: number of spaces to insert before/after field. * z1: number of zeros to insert before digits * z2: number of zeros to insert after digits * point: number of digits to print before decimal point * ndigits: number of digits to use from digits[] * suf: trailing suffix, like "e-5" */ realverb = verb; switch(verb){ case 'g': /* convert to at most prec significant digits. (prec=0 means 1) */ if(prec == 0) prec = 1; if(ndigits > prec) { if(digits[prec] >= '5' && add1(digits, prec)) exp++; exp += ndigits-prec; ndigits = prec; } /* * extra rules for %g (implemented below): * trailing zeros removed after decimal unless FmtSharp. * decimal point only if digit follows. */ /* fall through to %e */ default: case 'e': /* one significant digit before decimal, no leading zeros. */ pt = 1; z1 = 0; /* * decimal point is after ndigits digits right now. * slide to be after first. */ e = exp + (ndigits-1); /* if this is %g, check exponent and convert prec */ if(realverb == 'g') { if(-4 <= e && e < prec) goto casef; prec--; /* one digit before decimal; rest after */ } /* compute trailing zero padding or truncate digits. */ if(1+prec >= ndigits) z2 = 1+prec - ndigits; else { /* truncate digits */ assert(realverb != 'g'); newndigits = 1+prec; if(digits[newndigits] >= '5' && add1(digits, newndigits)) { /* had 999e4, now have 100e5 */ e++; } ndigits = newndigits; z2 = 0; } fmtexp(suf, e, ucase); nsuf = strlen(suf); break; casef: case 'f': /* determine where digits go with respect to decimal point */ if(ndigits+exp > 0) { pt = ndigits+exp; z1 = 0; } else { pt = 1; z1 = 1 + -(ndigits+exp); } /* * %g specifies prec = number of significant digits * convert to number of digits after decimal point */ if(realverb == 'g') prec += z1 - pt; /* compute trailing zero padding or truncate digits. */ if(pt+prec >= z1+ndigits) z2 = pt+prec - (z1+ndigits); else{ /* truncate digits */ assert(realverb != 'g'); newndigits = pt+prec - z1; if(newndigits < 0){ z1 += newndigits; newndigits = 0; }else if(newndigits == 0){ /* perhaps round up */ if(digits[0] >= '5'){ digits[0] = '1'; newndigits = 1; goto newdigit; } }else if(digits[newndigits] >= '5' && add1(digits, newndigits)){ /* digits was 999, is now 100; make it 1000 */ digits[newndigits++] = '0'; newdigit: /* account for new digit */ if(z1) /* 0.099 => 0.100 or 0.99 => 1.00*/ z1--; else /* 9.99 => 10.00 */ pt++; } z2 = 0; ndigits = newndigits; } nsuf = 0; break; } /* * if %g is given without FmtSharp, remove trailing zeros. * must do after truncation, so that e.g. print %.3g 1.001 * produces 1, not 1.00. sorry, but them's the rules. */ if(realverb == 'g' && !(f & fmt·Sharp)) { if(z1+ndigits+z2 >= pt) { if(z1+ndigits < pt) z2 = pt - (z1+ndigits); else{ z2 = 0; while(z1+ndigits > pt && digits[ndigits-1] == '0') ndigits--; } } } /* * compute width of all digits and decimal point and suffix if any */ n = z1+ndigits+z2; if(n > pt) n += ndot; else if(n == pt){ if(f & fmt·Sharp) n += ndot; else pt++; /* do not print any decimal point */ } n += nsuf; /* * determine sign */ sign = 0; if(neg) sign = '-'; else if(f & fmt·Sign) sign = '+'; else if(f & fmt·Space) sign = ' '; if(sign) n++; /* compute padding */ npad = 0; if((f & fmt·Width) && io->width > n) npad = io->width - n; if(npad && !(f & fmt·Left) && (f & fmt·Zero)){ z1 += npad; pt += npad; npad = 0; } /* format the actual field. too bad about doing this twice. */ if(npad && !(f & fmt·Left) && pad(io, npad < 0)) return -1; cur = io->buffer.cur; end = io->buffer.end; if(sign){ if(cur+1 > end){ if(!(cur=flush(io,cur,1))) return -1; end = io->buffer.end; } *cur++ = sign; } while(z1>0 || ndigits>0 || z2>0){ if(z1 > 0){ z1--; c = '0'; }else if(ndigits > 0){ ndigits--; c = *digits++; }else{ z2--; c = '0'; } if(cur+1 > end){ if(!(cur=flush(io,cur,1))) return -1; end = io->buffer.end; } *cur++ = c; if(--pt == 0) for(p=dot; *p; p++){ if(cur+1 > end){ if(!(cur=flush(io,cur,1))) return -1; end = io->buffer.end; } *cur++ = *p; } } io->n += cur - (char*)io->buffer.cur; io->buffer.cur = cur; if(nsuf && copy(io, suf, nsuf, nsuf) < 0) return -1; if(npad && (f & fmt·Left) && pad(io, npad < 0)) return -1; return 0; }