mrubyの浮動小数点数演算をソフトウェアで置き換える
この記事はmruby Advent Calendar 2017の23日目の記事です
前提環境
この記事ではx86_64アーキテクチャをターゲットにしています。
mrubyはmasterの2017年12月5日のコミット Need to set `ci->proc` when we have `RProc` structure. · mruby/mruby@4a99f49 · GitHub の状態です。
はじめに
一ヶ月ほど前mrubyにMRB_WITHOUT_FLOATという浮動小数点数関連の処理を取り除くオプションが追加されました。
今まではハードウェアでの浮動小数点演算に制限がある場合、自分でmrubyソースコードを改変して浮動小数点演算を取り除く必要がありましたが、このオプションが追加されたため、オプション一つで浮動小数点演算に制限がある環境でもmrubyを簡単に動作できるようになりました。
しかし、今回あえて浮動小数点演算に制限がある環境でMRB_WITHOUT_FLOATを使わないで別のやり方で対応してみました。 浮動小数点数演算にはハードウェアでの実行、ソフトウェアでの実現の2種類があります。 普段、x86_64ではハードウェアで浮動小数点演算を実行します。MRB_WITHOUT_FLOATは浮動小数点演算を排除することでハードウェアでの浮動小数点演算を無くそうということです。この記事ではmrubyの浮動小数点演算をハードウェアの実行ではなく、ソフトウェアで実現することで浮動小数点演算の制限がある環境でも問題なく浮動小数点演算を可能にする試みを紹介します。
gccだけでなんとかしたかった
浮動小数点演算に制限があるため、まずはビルドオプションでなんとかする方法を探ってみます。x86_64ではストリーミングSIMD拡張命令(SSE)を用いて浮動小数点演算を行ないます。なのでこのSSEを使わないオプションとsoftfloatを有効化するオプションをbuild_configに入れてみましょう。
MRuby::Build.new do |conf| toolchain :gcc enable_debug conf.cc.flags << "-mno-sse -mno-sse2 -mno-mmx -mno-3dnow -msoft-float" conf.cc.defines << %w(MRB_INT64) end
これで終わり！解散！とは行かずに
error : SSE register return with SSE disabledとビルドが失敗してしまいます。これはgccがx86_64アーキテクチャがすべてのx86-64プロセッサにSSEがあると仮定して上手くソフトウェアで変換してくれない所にあるようです。
gccだけで何とか出来ませんでした。
コンパイラの代わりに手で置き換える
コンパイラが上手に置き換えをやってくれないのでx86_64でmrubyの浮動小数点演算をソフトウェアで実現するにはmrubyの浮動小数点数周りを自分でどうにか置き換えてあげる必要がありますね。今回はBerkeley SoftFloatというライブラリを使って置き換えて行きたいと思います。
置き換え方法
詳細は公立はこだて未来大学アドベントカレンダー18日目にも書いてありますので基本的な部分に留めます。
float型はBerkeley SoftFloatの float32_tで置き換える
double型はfloat64_tで置き換える
float32_tやfloat64_tといった型は構造体として定義されており、構造体のメンバvで値を16進数で保持する
初期化は構造体メンバの初期化と同様に行うか下記の整数から変換する関数などを用いる 例: float32_t hoge = {0x3F800000};
整数からBerkeley SoftFloatの実数型への変換はi64_to_f64(uint64_t)といった関数を使用する(この場合はlong long int からfloat64_tへの変換)
Berkeley SoftFloatには実数型の四則演算と剰余、平方根関数が用意されているためそれを用いて計算する float64_t同士の加算の例: float64_t a = f64_add(b, c);
関係演算子 ==、<、<=の関数が用意されているのでそれを用いて比較する ==の例: bool hoge = f64_eq(a, b);
実際の置き換え
では、ここからは実際のmrubyのコードとdiffを見ながら実際に置き換えていきましょう。 しかし全部載せると記事がdiffまみれになってしまうので一部抜粋しながら紹介してきます。
parse.y
まずはmrubyの構文解析器から手を出していきましょう。
@@ -4969,12 +4969,12 @@ parser_yylex(parser_state *p) - double d; + float64_t d; char *endp; errno = 0; d = mrb_float_read(tok(p), &endp); - if (d == 0 && endp == tok(p)) { + if (f64_eq(d,i64_to_f64(0)) && endp == tok(p)) { yywarning_s(p, "corrupted float value %s", tok(p)); } else if (errno == ERANGE) {
parse.yの4972行目と4977行目の部分です。parser_yylex関数内のここではmrb_float_to_readに文字列トークンを渡して浮動小数点数に変換してますね。このdouble型変数はfloat64_tに置き換えます。float64_tに置き換えているため、if (d == 0 && endp == tok(p)) {のd ==0の部分も置き換えます。
vm.c
vm内のOP_ADD,OP_SUB,OP_MUL,OP_DIVといったオペコードの処理の部分も置き換えて行きます。
@@ -2198,8 +2200,18 @@ mrb_vm_exec(mrb_state *mrb, struct RProc *proc, mrb_code *pc) } #define TYPES2(a,b) ((((uint16_t)(a))<<8)|(((uint16_t)(b))&0xff)) -#define OP_MATH_BODY(op,v1,v2) do {\ - v1(regs[a]) = v1(regs[a]) op v2(regs[a+1]);\ + +#define OP_MATH_BODY_F64_F64(op,v1,v2) do {\ + v1(regs[a]) = op(v1(regs[a]),v2(regs[a+1]));\ +} while(0) +#define OP_MATH_BODY_I64_F64(op,v1,v2) do {\ + v1(regs[a]) = op(i64_to_f64(v1(regs[a])),v2(regs[a+1]));\ +} while(0) +#define OP_MATH_BODY_F64_I64(op,v1,v2) do {\ + v1(regs[a]) = op(v1(regs[a]),i64_to_f64(v2(regs[a+1])));\ +} while(0) +#define OP_MATH_BODY_I64_I64(op,v1,v2) do {\ + v1(regs[a]) = v1(regs[a]) op v2(regs[a+1]);\ } while(0) @@ -2217,7 +2229,7 @@ mrb_vm_exec(mrb_state *mrb, struct RProc *proc, mrb_code *pc) y = mrb_fixnum(regs_a[1]); if (mrb_int_add_overflow(x, y, &z)) { #ifndef MRB_WITHOUT_FLOAT - SET_FLOAT_VALUE(mrb, regs_a[0], (mrb_float)x + (mrb_float)y); + SET_FLOAT_VALUE(mrb, regs_a[0], f64_add(i64_to_f64(x),i64_to_f64(y))); break; #endif } @@ -2229,7 +2241,7 @@ mrb_vm_exec(mrb_state *mrb, struct RProc *proc, mrb_code *pc) { mrb_int x = mrb_fixnum(regs[a]); mrb_float y = mrb_float(regs[a+1]); - SET_FLOAT_VALUE(mrb, regs[a], (mrb_float)x + y); + SET_FLOAT_VALUE(mrb, regs[a], f64_add(i64_to_f64(x),y)); } break; case TYPES2(MRB_TT_FLOAT,MRB_TT_FIXNUM): @@ -2240,7 +2252,7 @@ mrb_vm_exec(mrb_state *mrb, struct RProc *proc, mrb_code *pc) SET_FLOAT_VALUE(mrb, regs[a], x + y); } #else - OP_MATH_BODY(+,mrb_float,mrb_fixnum); + OP_MATH_BODY_F64_I64(f64_add,mrb_float,mrb_fixnum); #endif break; case TYPES2(MRB_TT_FLOAT,MRB_TT_FLOAT): @@ -2251,7 +2263,7 @@ mrb_vm_exec(mrb_state *mrb, struct RProc *proc, mrb_code *pc) SET_FLOAT_VALUE(mrb, regs[a], x + y); } #else - OP_MATH_BODY(+,mrb_float,mrb_float); + OP_MATH_BODY_F64_F64(f64_add,mrb_float,mrb_float); #endif break; #endif
ここでは、オペランドがどちらも数値だった場合、vmで直接計算をしていますが、浮動小数点数と整数が混合するのでここでは上手く分けてあげなければいけません。置き換え前のOP_MATH_BODYマクロは上手くキャストして計算していましたが、浮動小数点数をソフトウェアで再現するにあたって浮動小数点数の四則演算はライブラリの関数を使ってあげる必要があるため、そのままでは上手く動きません。ですので整数・整数の場合、整数・浮動小数点数の場合、浮動小数点数・整数の場合、浮動小数点数・浮動小数点数の場合の４つで考えて上げる必要があります。また、SET_FLOAT_VALUE(mrb, regs[a], (mrb_float)x + y);といった計算もmrb_floatがfloat64_tに置き換えられているため、このままだとキャストが上手く動きませんのでここも置き換えていきましょう。
次にOP_EQやOP_LT、OP_LE、OP_GT、OP_GEといったオペコードも同様に浮動小数点数と整数が混合するのでここでは上手く分けてあげなければいけません。
@@ -2537,21 +2556,48 @@ mrb_vm_exec(mrb_state *mrb, struct RProc *proc, mrb_code *pc) }\ } while(0) #else -#define OP_CMP(op) do {\ +#define OP_CMP(op1,op2) do {\ int result;\ /* need to check if - is overridden */\ switch (TYPES2(mrb_type(regs[a]),mrb_type(regs[a+1]))) {\ case TYPES2(MRB_TT_FIXNUM,MRB_TT_FIXNUM):\ - result = OP_CMP_BODY(op,mrb_fixnum,mrb_fixnum);\ + result = OP_CMP_BODY_I64_I64(op1,mrb_fixnum,mrb_fixnum);\ + break;\ + case TYPES2(MRB_TT_FIXNUM,MRB_TT_FLOAT):\ + result = OP_CMP_BODY_I64_F64(op2,mrb_fixnum,mrb_float);\ + break;\ + case TYPES2(MRB_TT_FLOAT,MRB_TT_FIXNUM):\ + result = OP_CMP_BODY_F64_I64(op2,mrb_float,mrb_fixnum);\ + break;\ + case TYPES2(MRB_TT_FLOAT,MRB_TT_FLOAT):\ + result = OP_CMP_BODY_F64_F64(op2,mrb_float,mrb_float);\ + break;\ + default:\ + goto L_SEND;\ + }\ + if (result) {\ + SET_TRUE_VALUE(regs[a]);\ + }\ + else {\ + SET_FALSE_VALUE(regs[a]);\ + }\ +} while(0) + +#define OP_CMP_REV(op1,op2) do {\ + int result;\ + /* need to check if - is overridden */\ + switch (TYPES2(mrb_type(regs[a]),mrb_type(regs[a+1]))) {\ + case TYPES2(MRB_TT_FIXNUM,MRB_TT_FIXNUM):\ + result = OP_CMP_BODY_REV_I64_I64(op1,mrb_fixnum,mrb_fixnum);\ break;\ case TYPES2(MRB_TT_FIXNUM,MRB_TT_FLOAT):\ - result = OP_CMP_BODY(op,mrb_fixnum,mrb_float);\ + result = OP_CMP_BODY_REV_I64_F64(op2,mrb_fixnum,mrb_float);\ break;\ case TYPES2(MRB_TT_FLOAT,MRB_TT_FIXNUM):\ - result = OP_CMP_BODY(op,mrb_float,mrb_fixnum);\ + result = OP_CMP_BODY_REV_F64_I64(op2,mrb_float,mrb_fixnum);\ break;\ case TYPES2(MRB_TT_FLOAT,MRB_TT_FLOAT):\ - result = OP_CMP_BODY(op,mrb_float,mrb_float);\ + result = OP_CMP_BODY_REV_F64_F64(op2,mrb_float,mrb_float);\ break;\ default:\ goto L_SEND;\ @@ -2572,36 +2618,36 @@ mrb_vm_exec(mrb_state *mrb, struct RProc *proc, mrb_code *pc) SET_TRUE_VALUE(regs[a]); } else { - OP_CMP(==); + OP_CMP(==,f64_eq); } NEXT; } CASE(OP_LT) { /* A B C R(A) := R(A)<R(A+1) (Syms[B]=:<,C=1)*/ int a = GETARG_A(i); - OP_CMP(<); + OP_CMP(<,f64_lt); NEXT; } CASE(OP_LE) { /* A B C R(A) := R(A)<=R(A+1) (Syms[B]=:<=,C=1)*/ int a = GETARG_A(i); - OP_CMP(<=); + OP_CMP(<=,f64_le); NEXT; } CASE(OP_GT) { /* A B C R(A) := R(A)>R(A+1) (Syms[B]=:>,C=1)*/ int a = GETARG_A(i); - OP_CMP(>); + OP_CMP_REV(>,f64_lt); NEXT; } CASE(OP_GE) { /* A B C R(A) := R(A)>=R(A+1) (Syms[B]=:>=,C=1)*/ int a = GETARG_A(i); - OP_CMP(>=); + OP_CMP_REV(>=,f64_lt); NEXT; }
OP_CMPマクロは型を判定して引数の関係演算子の比較結果をresultに代入していますが、浮動小数点数の比較はBerkeley SoftFloatの関数を使って行ってあげないといけないため置き換える必要があります。なので今回はOP_CMPマクロに渡す引数を2つにして関係演算子とBerkeley SoftFloatのメソッド名を渡すようにしました。また、Berkeley SoftFloatの比較を行う関数はless than（未満）とless than or equal (以下）、equal(一致)のみなので、以上と超については比較させる内容を反転させたりしてうまくやって上げる必要があります。本当はもうちょっといいやり方があったかもしれませんが、今回はOP_CMPとOP＿CMP_REVの2つを用意して対処しました。
残りnumeric.c numeric.cにはNumeric, Integer, Float, Fixnumといった各クラスとメソッドが実装されています。ここでは様々な数値関連メソッドが定義されているため、numeric.cが一番置き換える量が多いです。また、math.hとfloat.hといったヘッダファイルがインクルードされていますが、これらは浮動小数点数の関係で使うことができません。ですので定義されている関数っぽいものををBerkeley SoftFloatを用いて自分で作って上げる必要があります。
static float64_t NAN = {0xFFFFFFFFFFFFFFFF}; static float64_t INFINITY = {0x7FF0000000000000}; static float64_t NEGATIVE_INFINITY = {0xFFF0000000000000}; #define f64_isnan( a ) (((~(a.v) & UINT64_C( 0x7FF0000000000000 )) == 0) && ((a.v) & UINT64_C( 0x000FFFFFFFFFFFFF ))) #define f64_floor (a) (f64_roundToInt(a,softfloat_round_min,0)) #define f64_ceil (a) (f64_roundToInt(a,softfloat_round_max,0)) float64_t f64_pow(float64_t a,float64_t b) { int n = f64_to_i64(b); float64_t buf = {a.v}; if(b.v == 0 ||a.v == 0x3FF0000000000000){ buf.v =0x3FF0000000000000; return buf; } if(a.v == 0){ buf.v = 0; return buf; } for(int m = 1; m < n; m++){ buf = f64_mul(buf,a); } return buf; } int f64_isinf (float64_t f) { if ((((f.v>>52) & 0x07FF) != 0x07FF) || (f.v<<12) != 0)){ return 0; } return 1; } int f64_isfinite(float64_t f) { if(f64_isinf(f) || f64_isnan(f)) return 0; return 1; }
まず、NANやINFINITYといったmath.hにある定数をfloat64_tで定義してあげます。これらの定数はflodivmod(剰余を計算する関数）で必要となります。また、floor(切り捨て)やceil(切り上げ)はBerkeley SoftFloatの64_roundToIntを用いて定義します。isinf(無限大か判定)、isfinite(有限値か判定)に関しては自分で実装します。 これらを実装することでnumeric.cを置き換える準備ができましたので、実際に置き換えていきます。
@@ -233,21 +258,21 @@ flodivmod(mrb_state *mrb, mrb_float x, mrb_float y, mrb_float *divp, mrb_float * mrb_float div; mrb_float mod; - if (y == 0.0) { - if (x > 0.0) div = INFINITY; - else if (x < 0.0) div = -INFINITY; + if (f64_eq(y,i64_to_f64(0))) { + if (f64_lt(i64_to_f64(0),x)) div = INFINITY; + else if (f64_lt(x,i64_to_f64(0))) div = NEGATIVE_INFINITY; else div = NAN; /* x == 0.0 */ mod = NAN; } else { - mod = fmod(x, y); - if (isinf(x) && isfinite(y)) + mod = f64_rem(x, y); + if (f64_isinf(x) && f64_isfinite(y)) div = x; else - div = (x - mod) / y; - if (y*mod < 0) { - mod += y; - div -= 1.0; + div = f64_div(f64_sub(x,mod),y); + if (f64_lt(f64_mul(y,mod),i64_to_f64(0))){ + mod = f64_add(mod,y); + div = f64_sub(div,i64_to_f64(1));
特にこれといって今までの置き換えと変わりませんが、先程定義した定数とかを使用しておきかえています。
fmt_fp.c
最後にfmt_fp.cをを置き換えていきます。置き換えの量はnumeric.cに比べれば少ないですが、こちらのほうが難しいです。ここではmusl libcのsrc/stdio/vfprintf.c内にあるvfprintfをベースに浮動小数点数から文字列に変換する関数があります。これを置き換える必要があるのですが、まずそのためにはfrexp(仮数部と指数部を分割する関数)を実装してあげる必要があります。これを一から自分で作るのは非常にしんどいので、今回はfdlibmのfrexpの実装をベースにBerkeley SoftFloatを使って改変していきます。
#define __HI(x) *(1+(int *)&x.v) #define __LO(x) *(int *)&x.v static float64_t two54 = {0x40000000000000}; /* 0x43500000, 0x00000000 */ float64_t frexp(float64_t x, int *eptr) { int hx, ix, lx; hx = __HI(x); ix = 0x7fffffff&hx; lx = __LO(x); *eptr = 0; if(ix>=0x7ff00000||((ix|lx)==0)) return x; /* 0,inf,nan */ if (ix<0x00100000) { /* subnormal */ x =f64_mul(x,two54); hx = __HI(x); ix = hx&0x7fffffff; *eptr = -54; } *eptr += (ix>>20)-1022; hx = (hx&0x800fffff)|0x3fe00000; __HI(x) = hx; return x; }
さて、ここまで改変してあげたので準備万端と言いたいところなのですがfmt_fpの中にはlong doubleいうちょっとめんどくさい存在が使われています。拡張倍精度浮動小数点フォーマットと呼ばれるものです。ちなみにIEEE754でも拡張倍精度について具体的な形式について定義はしていなかったりします。ただ、幸いにもBerkeley SoftFloatには拡張倍精度浮動小数点フォーマットもextFloat80_tという名前で定義されているのでこれを使っていきましょう。
+ float64_t epsilon={0x3C00000000000000}; +#define extF80_signbit(f)((f.signExp>>15)) +#define extF80_isnan(a) (((~(a.signExp) & UINT64_C(0x7FFF)) == 0) && ((a.signif) & UINT64_C(0xFFFFFFFFFFFFFFFF))) +#define extF80_isinf(a) (((a.signExp & 0x7FFF)==0x7FFF) && (a.signif == 0)) +#define extF80_to_i64(x) extF80_to_i64((x),softfloat_round_min,1) +#define f64_to_i64(x) f64_to_i64((x),softfloat_round_min,1) +int +extF80_isfinite(extFloat80_t f) +{ + if(extF80_isinf(f) || extF80_isnan(f)) return 0; + return 1; +} … -fmt_fp(struct fmt_args *f, long double y, ptrdiff_t p, uint8_t fl, int t) +fmt_fp(struct fmt_args *f, extFloat80_t y, ptrdiff_t p, uint8_t 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 @@ -106,29 +140,29 @@ fmt_fp(struct fmt_args *f, long double y, ptrdiff_t p, uint8_t fl, int t) char ebuf0[3*sizeof(int)], *ebuf=&ebuf0[3*sizeof(int)], *estr; pl=1; - if (signbit(y)) { - y=-y; + if (extF80_signbit(y)) { + y = extF80_mul(y,i64_to_extF80(-1)); } else if (fl & MARK_POS) { prefix+=3; } else if (fl & PAD_POS) { prefix+=6; } else prefix++, pl=0; - if (!isfinite(y)) { + if (!extF80_isfinite(y)) { const char *ss = (t&32)?"inf":"INF"; - if (y!=y) ss=(t&32)?"nan":"NAN"; + if (!extF80_eq(y,y)) ss=(t&32)?"nan":"NAN"; pad(f, ' ', 0, 3+pl, fl&~ZERO_PAD); out(f, prefix, pl); out(f, ss, 3); pad(f, ' ', 0, 3+pl, fl^LEFT_ADJ); return 3+(int)pl; } - y = frexp((double)y, &e2) * 2; - if (y) e2--; + y = f64_to_extF80(f64_mul(frexp(extF80_to_f64(y), &e2),i64_to_f64(2))); + if (!extF80_eq(y,i64_to_extF80(0))) e2--; if ((t|32)=='a') { - long double round = 8.0; + extFloat80_t round = i64_to_extF80(8); ptrdiff_t re; if (t&32) prefix += 9; @@ -138,16 +172,16 @@ fmt_fp(struct fmt_args *f, long double y, ptrdiff_t p, uint8_t fl, int t) else re=LDBL_MANT_DIG/4-1-p; if (re) { - while (re--) round*=16; + while (re--) round = extF80_mul(round,i64_to_extF80(16)); if (*prefix=='-') { - y=-y; - y-=round; - y+=round; - y=-y; + y=extF80_mul(y,i64_to_extF80(-1)); + y=extF80_sub(y,round); + y=extF80_add(y,round); + y=extF80_sub(y,y); } else { - y+=round; - y-=round; + y=extF80_add(y,round); + y=extF80_sub(y,round); } } @@ -158,11 +192,11 @@ fmt_fp(struct fmt_args *f, long double y, ptrdiff_t p, uint8_t fl, int t) s=buf; do { - int x=(int)y; + int x=extF80_to_i64(y); *s++=xdigits[x]|(t&32); - y=16*(y-x); - if (s-buf==1 && (y||p>0||(fl&ALT_FORM))) *s++='.'; - } while (y); + y=extF80_mul(i64_to_extF80(16),extF80_sub(y,i64_to_extF80(x))); + if (s-buf==1 && (!extF80_eq(y,i64_to_extF80(0))||p>0||(fl&ALT_FORM))) *s++='.'; + } while (!extF80_eq(y,i64_to_extF80(0))); if (p && s-buf-2 < p) l = (p+2) + (ebuf-estr); @@ -180,15 +214,15 @@ fmt_fp(struct fmt_args *f, long double y, ptrdiff_t p, uint8_t fl, int t) } if (p<0) p=6; - if (y) y *= 268435456.0, e2-=28; + if (!extF80_eq(y,i64_to_extF80(0))) y =extF80_mul(y,i64_to_extF80(268435456)), e2-=28; if (e2<0) a=r=z=big; else a=r=z=big+sizeof(big)/sizeof(*big) - LDBL_MANT_DIG - 1; do { - *z = (uint32_t)y; - y = 1000000000*(y-*z++); - } while (y); + *z = extF80_to_ui32(y,softfloat_round_min,1); + y = extF80_mul(i64_to_extF80(1000000000),extF80_sub(y,i32_to_extF80(*z++))); + } while (!extF80_eq(y,i64_to_extF80(0))); while (e2>0) { uint32_t carry=0; @@ -233,16 +267,19 @@ fmt_fp(struct fmt_args *f, long double y, ptrdiff_t p, uint8_t fl, int t) 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; + float64_t small1 ={0x3FE0000000000000}; + float64_t small2 ={0xFF8000000000000}; + extFloat80_t round = extF80_div(i64_to_extF80(2),f64_to_extF80(epsilon)); + extFloat80_t small; + if (*d/i & 1) round =extF80_add(round,i64_to_extF80(2)); + if (x<i/2) small=f64_to_extF80(small1); + else if (x==i/2 && d+1==z) small=i64_to_extF80(1.0); + else small=f64_to_extF80(small2); + if (pl && *prefix=='-') round=extF80_mul(round,i64_to_extF80(-1)), small=extF80_mul(small,i64_to_extF80(-1)); *d -= x; /* Decide whether to round by probing round+small */ - if (round+small != round) { + if (!extF80_eq(extF80_add(round,small),round)){ *d = *d + i; while (*d > 999999999) { *d--=0;
まず、extFloat80_tを使うにあたって一つ注意するべきことがあります。extFloat80_tは80bitという都合上、uint16_t signExpとuint64_t signif の2つのメンバで構成されています。また、この構造体extFloat80_tはリトルエンディアンとビッグエンディアンで定義が違うので、ビルド時にオプションでLITTLEENDIANと指定しましょう。(x86_64がリトルエンディアンだから)そうしないとエンディアンで殺されてしまいます。また、isnanやisinfをマクロや関数で作ってあげるときも他の構造体とはメンバが違うことに注意して作ってあげましょう。そこにさえ気にすれば問題なく動きます。以上で置き換えは大体終了です。
おわりに
正直な話、量が大きいnumeric.cよりもfmt_fp.cやvm.cの置き換えのほうがしんどいです。やっぱり、量より内容ですね。今回載せた置き換えは一部ですが、こちらのほうに置き換えた全てが載っています。まあ興味がある人はよっぽど居ないと思いますが一応載せておきます。
github.com この置き換え自体はなかなか使うケースが無いとは思いますが、ハイパーバイザやLinuxカーネルで動かすための一つの対応手段として参考になればと思います。(僕自身はこの置き換えを使ってハイパーバイザにmrubyを突っ込んでいます)。
P.S. 2日遅れになってしまいました、すみません。
