diff options
Diffstat (limited to 'src')
| -rw-r--r-- | src/dump.c | 22 | ||||
| -rw-r--r-- | src/fmt_fp.c | 367 | ||||
| -rw-r--r-- | src/numeric.c | 149 |
3 files changed, 387 insertions, 151 deletions
diff --git a/src/dump.c b/src/dump.c index 24ddd90b8..3113a71d3 100644 --- a/src/dump.c +++ b/src/dump.c @@ -18,6 +18,12 @@ #ifdef ENABLE_STDIO +#ifdef MRB_USE_FLOAT +#define MRB_FLOAT_FMT "%.9e" +#else +#define MRB_FLOAT_FMT "%.16e" +#endif + static size_t get_irep_record_size_1(mrb_state *mrb, mrb_irep *irep); #if UINT32_MAX > SIZE_MAX @@ -111,7 +117,6 @@ get_pool_block_size(mrb_state *mrb, mrb_irep *irep) size_t size = 0; size_t pool_no; mrb_value str; - char buf[32]; size += sizeof(uint32_t); /* plen */ size += irep->plen * (sizeof(uint8_t) + sizeof(uint16_t)); /* len(n) */ @@ -130,9 +135,9 @@ get_pool_block_size(mrb_state *mrb, mrb_irep *irep) break; case MRB_TT_FLOAT: + str = mrb_float_to_str(mrb, irep->pool[pool_no], MRB_FLOAT_FMT); { - int len; - len = mrb_float_to_str(buf, mrb_float(irep->pool[pool_no])); + mrb_int len = RSTRING_LEN(str); mrb_assert_int_fit(mrb_int, len, size_t, SIZE_MAX); size += (size_t)len; } @@ -163,7 +168,6 @@ write_pool_block(mrb_state *mrb, mrb_irep *irep, uint8_t *buf) uint16_t len; mrb_value str; const char *char_ptr; - char char_buf[30]; cur += uint32_to_bin(irep->plen, cur); /* number of pool */ @@ -186,13 +190,15 @@ write_pool_block(mrb_state *mrb, mrb_irep *irep, uint8_t *buf) case MRB_TT_FLOAT: cur += uint8_to_bin(IREP_TT_FLOAT, cur); /* data type */ + str = mrb_float_to_str(mrb, irep->pool[pool_no], MRB_FLOAT_FMT); + char_ptr = RSTRING_PTR(str); { - int tlen; - tlen = mrb_float_to_str(char_buf, mrb_float(irep->pool[pool_no])); - mrb_assert_int_fit(int, tlen, uint16_t, UINT16_MAX); + mrb_int tlen; + + tlen = RSTRING_LEN(str); + mrb_assert_int_fit(mrb_int, tlen, uint16_t, UINT16_MAX); len = (uint16_t)tlen; } - char_ptr = &char_buf[0]; break; case MRB_TT_STRING: diff --git a/src/fmt_fp.c b/src/fmt_fp.c new file mode 100644 index 000000000..9345c31d1 --- /dev/null +++ b/src/fmt_fp.c @@ -0,0 +1,367 @@ +/* + +Most code in this file originates from musl (src/stdio/vfprintf.c) +which, just like mruby itself, is licensed under the MIT license. + +Copyright © 2005-2014 Rich Felker, et al. + +Permission is hereby granted, free of charge, to any person obtaining +a copy of this software and associated documentation files (the +"Software"), to deal in the Software without restriction, including +without limitation the rights to use, copy, modify, merge, publish, +distribute, sublicense, and/or sell copies of the Software, and to +permit persons to whom the Software is furnished to do so, subject to +the following conditions: + +The above copyright notice and this permission notice shall be +included in all copies or substantial portions of the Software. + +THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. +IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY +CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, +TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE +SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +*/ + +#include <limits.h> +#include <string.h> +#include <stdint.h> +#include <math.h> +#include <float.h> +#include <ctype.h> + +#include "mruby.h" +#include "mruby/string.h" + +struct fmt_args { + mrb_state *mrb; + mrb_value str; +}; + +#define MAX(a,b) ((a)>(b) ? (a) : (b)) +#define MIN(a,b) ((a)<(b) ? (a) : (b)) + +/* Convenient bit representation for modifier flags, which all fall + * within 31 codepoints of the space character. */ + +#define ALT_FORM (1U<<('#'-' ')) +#define ZERO_PAD (1U<<('0'-' ')) +#define LEFT_ADJ (1U<<('-'-' ')) +#define PAD_POS (1U<<(' '-' ')) +#define MARK_POS (1U<<('+'-' ')) + +static void +out(struct fmt_args *f, const char *s, size_t l) +{ + mrb_str_cat(f->mrb, f->str, s, l); +} + +static void +pad(struct fmt_args *f, char c, int w, int l, int fl) +{ + char pad[256]; + if (fl & (LEFT_ADJ | ZERO_PAD) || l >= w) return; + l = w - l; + memset(pad, c, l>sizeof pad ? sizeof pad : l); + for (; l >= sizeof pad; l -= sizeof pad) + out(f, pad, sizeof pad); + out(f, pad, l); +} + +static const char xdigits[16] = { + '0', '1', '2', '3', '4', '5', '6', '7', '8', '9', 'A', 'B', 'C', 'D', 'E', 'F' +}; + +static char* +fmt_u(uint32_t x, char *s) +{ + for (; x; x /= 10) *--s = '0' + x % 10; + return s; +} + +/* Do not override this check. The floating point printing code below + * depends on the float.h constants being right. If they are wrong, it + * may overflow the stack. */ +#if LDBL_MANT_DIG == 53 +typedef char compiler_defines_long_double_incorrectly[9-(int)sizeof(long double)]; +#endif + +static int +fmt_fp(struct fmt_args *f, long double y, int w, int p, int fl, int t) +{ + uint32_t big[(LDBL_MANT_DIG+28)/29 + 1 // mantissa expansion + + (LDBL_MAX_EXP+LDBL_MANT_DIG+28+8)/9]; // exponent expansion + uint32_t *a, *d, *r, *z; + uint32_t i; + int e2=0, e, j, l; + char buf[9+LDBL_MANT_DIG/4], *s; + const char *prefix="-0X+0X 0X-0x+0x 0x"; + int pl; + char ebuf0[3*sizeof(int)], *ebuf=&ebuf0[3*sizeof(int)], *estr; + + pl=1; + if (signbit(y)) { + y=-y; + } else if (fl & MARK_POS) { + prefix+=3; + } else if (fl & PAD_POS) { + prefix+=6; + } else prefix++, pl=0; + + if (!isfinite(y)) { + char *ss = (t&32)?"inf":"INF"; + if (y!=y) ss=(t&32)?"nan":"NAN"; + pad(f, ' ', w, 3+pl, fl&~ZERO_PAD); + out(f, prefix, pl); + out(f, ss, 3); + pad(f, ' ', w, 3+pl, fl^LEFT_ADJ); + return MAX(w, 3+pl); + } + + y = frexpl(y, &e2) * 2; + if (y) e2--; + + if ((t|32)=='a') { + long double round = 8.0; + int re; + + if (t&32) prefix += 9; + pl += 2; + + if (p<0 || p>=LDBL_MANT_DIG/4-1) re=0; + else re=LDBL_MANT_DIG/4-1-p; + + if (re) { + while (re--) round*=16; + if (*prefix=='-') { + y=-y; + y-=round; + y+=round; + y=-y; + } + else { + y+=round; + y-=round; + } + } + + estr=fmt_u(e2<0 ? -e2 : e2, ebuf); + if (estr==ebuf) *--estr='0'; + *--estr = (e2<0 ? '-' : '+'); + *--estr = t+('p'-'a'); + + s=buf; + do { + int x=y; + *s++=xdigits[x]|(t&32); + y=16*(y-x); + if (s-buf==1 && (y||p>0||(fl&ALT_FORM))) *s++='.'; + } while (y); + + if (p && s-buf-2 < p) + l = (p+2) + (ebuf-estr); + else + l = (s-buf) + (ebuf-estr); + + pad(f, ' ', w, pl+l, fl); + out(f, prefix, pl); + pad(f, '0', w, pl+l, fl^ZERO_PAD); + out(f, buf, s-buf); + pad(f, '0', l-(ebuf-estr)-(s-buf), 0, 0); + out(f, estr, ebuf-estr); + pad(f, ' ', w, pl+l, fl^LEFT_ADJ); + return MAX(w, pl+l); + } + if (p<0) p=6; + + if (y) y *= 268435456.0, e2-=28; + + if (e2<0) a=r=z=big; + else a=r=z=big+sizeof(big)/sizeof(*big) - LDBL_MANT_DIG - 1; + + do { + *z = y; + y = 1000000000*(y-*z++); + } while (y); + + while (e2>0) { + uint32_t carry=0; + int sh=MIN(29,e2); + for (d=z-1; d>=a; d--) { + uint64_t x = ((uint64_t)*d<<sh)+carry; + *d = x % 1000000000; + carry = x / 1000000000; + } + if (carry) *--a = carry; + while (z>a && !z[-1]) z--; + e2-=sh; + } + while (e2<0) { + uint32_t carry=0, *b; + int sh=MIN(9,-e2), need=1+(p+LDBL_MANT_DIG/3+8)/9; + for (d=a; d<z; d++) { + uint32_t rm = *d & (1<<sh)-1; + *d = (*d>>sh) + carry; + carry = (1000000000>>sh) * rm; + } + if (!*a) a++; + if (carry) *z++ = carry; + /* Avoid (slow!) computation past requested precision */ + b = (t|32)=='f' ? r : a; + if (z-b > need) z = b+need; + e2+=sh; + } + + if (a<z) for (i=10, e=9*(r-a); *a>=i; i*=10, e++); + else e=0; + + /* Perform rounding: j is precision after the radix (possibly neg) */ + j = p - ((t|32)!='f')*e - ((t|32)=='g' && p); + if (j < 9*(z-r-1)) { + uint32_t x; + /* We avoid C's broken division of negative numbers */ + d = r + 1 + ((j+9*LDBL_MAX_EXP)/9 - LDBL_MAX_EXP); + j += 9*LDBL_MAX_EXP; + j %= 9; + for (i=10, j++; j<9; i*=10, j++); + x = *d % i; + /* Are there any significant digits past j? */ + if (x || d+1!=z) { + long double round = 2/LDBL_EPSILON; + long double small; + if (*d/i & 1) round += 2; + if (x<i/2) small=0.5; + else if (x==i/2 && d+1==z) small=1.0; + else small=1.5; + if (pl && *prefix=='-') round*=-1, small*=-1; + *d -= x; + /* Decide whether to round by probing round+small */ + if (round+small != round) { + *d = *d + i; + while (*d > 999999999) { + *d--=0; + if (d<a) *--a=0; + (*d)++; + } + for (i=10, e=9*(r-a); *a>=i; i*=10, e++); + } + } + if (z>d+1) z=d+1; + } + for (; z>a && !z[-1]; z--); + + if ((t|32)=='g') { + if (!p) p++; + if (p>e && e>=-4) { + t--; + p-=e+1; + } + else { + t-=2; + p--; + } + if (!(fl&ALT_FORM)) { + /* Count trailing zeros in last place */ + if (z>a && z[-1]) for (i=10, j=0; z[-1]%i==0; i*=10, j++); + else j=9; + if ((t|32)=='f') + p = MIN(p,MAX(0,9*(z-r-1)-j)); + else + p = MIN(p,MAX(0,9*(z-r-1)+e-j)); + } + } + l = 1 + p + (p || (fl&ALT_FORM)); + if ((t|32)=='f') { + if (e>0) l+=e; + } + else { + estr=fmt_u(e<0 ? -e : e, ebuf); + while(ebuf-estr<2) *--estr='0'; + *--estr = (e<0 ? '-' : '+'); + *--estr = t; + l += ebuf-estr; + } + + pad(f, ' ', w, pl+l, fl); + out(f, prefix, pl); + pad(f, '0', w, pl+l, fl^ZERO_PAD); + + if ((t|32)=='f') { + if (a>r) a=r; + for (d=a; d<=r; d++) { + char *ss = fmt_u(*d, buf+9); + if (d!=a) while (ss>buf) *--ss='0'; + else if (ss==buf+9) *--ss='0'; + out(f, ss, buf+9-ss); + } + if (p || (fl&ALT_FORM)) out(f, ".", 1); + for (; d<z && p>0; d++, p-=9) { + char *ss = fmt_u(*d, buf+9); + while (ss>buf) *--ss='0'; + out(f, ss, MIN(9,p)); + } + pad(f, '0', p+9, 9, 0); + } + else { + if (z<=a) z=a+1; + for (d=a; d<z && p>=0; d++) { + char *ss = fmt_u(*d, buf+9); + if (ss==buf+9) *--ss='0'; + if (d!=a) while (ss>buf) *--ss='0'; + else { + out(f, ss++, 1); + if (p>0||(fl&ALT_FORM)) out(f, ".", 1); + } + out(f, ss, MIN(buf+9-ss, p)); + p -= buf+9-ss; + } + pad(f, '0', p+18, 18, 0); + out(f, estr, ebuf-estr); + } + + pad(f, ' ', w, pl+l, fl^LEFT_ADJ); + + return MAX(w, pl+l); +} + +static int +fmt_core(struct fmt_args *f, const char *fmt, mrb_float flo) +{ + int p; + + if (*fmt != '%') { + return -1; + } + ++fmt; + + if (*fmt == '.') { + ++fmt; + for (p = 0; ISDIGIT(*fmt); ++fmt) { + p = 10 * p + (*fmt - '0'); + } + } + else { + p = -1; + } + + switch (*fmt) { + case 'e': case 'f': case 'g': case 'a': + case 'E': case 'F': case 'G': case 'A': + return fmt_fp(f, flo, 0, p, 0, *fmt); + default: + return -1; + } +} + +mrb_value +mrb_float_to_str(mrb_state *mrb, mrb_value flo, const char *fmt) +{ + struct fmt_args f = { mrb, mrb_str_buf_new(mrb, 24) }; + if (fmt_core(&f, fmt, mrb_float(flo)) < 0) { + mrb_raise(mrb, E_ARGUMENT_ERROR, "invalid format string"); + } + return f.str; +} diff --git a/src/numeric.c b/src/numeric.c index e01e661a5..8b6ec4c88 100644 --- a/src/numeric.c +++ b/src/numeric.c @@ -18,13 +18,9 @@ #define floor(f) floorf(f) #define ceil(f) ceilf(f) #define fmod(x,y) fmodf(x,y) -#define FLO_MAX_DIGITS 7 -#define FLO_MAX_SIGN_LENGTH 3 -#define FLO_EPSILON FLT_EPSILON +#define MRB_FLO_TO_STR_FMT "%.7g" #else -#define FLO_MAX_DIGITS 14 -#define FLO_MAX_SIGN_LENGTH 10 -#define FLO_EPSILON DBL_EPSILON +#define MRB_FLO_TO_STR_FMT "%.14g" #endif MRB_API mrb_float @@ -107,142 +103,6 @@ num_div(mrb_state *mrb, mrb_value x) * representation. */ -static mrb_value -mrb_flo_to_str(mrb_state *mrb, mrb_float flo) -{ - double n = (double)flo; - int max_digits = FLO_MAX_DIGITS; - - if (isnan(n)) { - return mrb_str_new_lit(mrb, "NaN"); - } - else if (isinf(n)) { - if (n < 0) { - return mrb_str_new_lit(mrb, "-inf"); - } - else { - return mrb_str_new_lit(mrb, "inf"); - } - } - else { - int digit; - int m = 0; - int exp; - mrb_bool e = FALSE; - char s[48]; - char *c = &s[0]; - int length = 0; - - if (signbit(n)) { - n = -n; - *(c++) = '-'; - } - - if (n != 0.0) { - if (n > 1.0) { - exp = (int)floor(log10(n)); - } - else { - exp = (int)-ceil(-log10(n)); - } - } - else { - exp = 0; - } - - /* preserve significands */ - if (exp < 0) { - int i, beg = -1, end = 0; - double f = n; - double fd = 0; - for (i = 0; i < FLO_MAX_DIGITS; ++i) { - f = (f - fd) * 10.0; - fd = floor(f + FLO_EPSILON); - if (fd != 0) { - if (beg < 0) beg = i; - end = i + 1; - } - } - if (beg >= 0) length = end - beg; - if (length > FLO_MAX_SIGN_LENGTH) length = FLO_MAX_SIGN_LENGTH; - } - - if (abs(exp) + length >= FLO_MAX_DIGITS) { - /* exponent representation */ - e = TRUE; - n = n / pow(10.0, exp); - if (isinf(n)) { - if (s < c) { /* s[0] == '-' */ - return mrb_str_new_lit(mrb, "-0.0"); - } - else { - return mrb_str_new_lit(mrb, "0.0"); - } - } - } - else { - /* un-exponent (normal) representation */ - if (exp > 0) { - m = exp; - } - } - - /* puts digits */ - while (max_digits >= 0) { - double weight = (m < 0) ? 0.0 : pow(10.0, m); - double fdigit = (m < 0) ? n * 10.0 : n / weight; - - if (fdigit < 0) fdigit = n = 0; - if (m < -1 && fdigit < FLO_EPSILON) { - if (e || exp > 0 || m <= -abs(exp)) { - break; - } - } - digit = (int)floor(fdigit + FLO_EPSILON); - if (m == 0 && digit > 9) { - n /= 10.0; - exp++; - continue; - } - *(c++) = '0' + digit; - n = (m < 0) ? n * 10.0 - digit : n - (digit * weight); - max_digits--; - if (m-- == 0) { - *(c++) = '.'; - } - } - if (c[-1] == '0') { - while (&s[0] < c && c[-1] == '0') { - c--; - } - c++; - } - - if (e) { - *(c++) = 'e'; - if (exp > 0) { - *(c++) = '+'; - } - else { - *(c++) = '-'; - exp = -exp; - } - - if (exp >= 100) { - *(c++) = '0' + exp / 100; - exp -= exp / 100 * 100; - } - - *(c++) = '0' + exp / 10; - *(c++) = '0' + exp % 10; - } - - *c = '\0'; - - return mrb_str_new(mrb, &s[0], c - &s[0]); - } -} - /* 15.2.9.3.16(x) */ /* * call-seq: @@ -257,7 +117,10 @@ mrb_flo_to_str(mrb_state *mrb, mrb_float flo) static mrb_value flo_to_s(mrb_state *mrb, mrb_value flt) { - return mrb_flo_to_str(mrb, mrb_float(flt)); + if (isnan(mrb_float(flt))) { + return mrb_str_new_lit(mrb, "NaN"); + } + return mrb_float_to_str(mrb, flt, MRB_FLO_TO_STR_FMT); } /* 15.2.9.3.2 */ |