2008年02月15日

浮動小数点数をできるだけ短い10進数に可逆変換する

はてなブックマークに登録

以前もこのBlogで紹介したswfmillは、swf内の情報をXML形式に変換したりその逆変換ができたりするのですが、swf内の実数をXML内では10進表記している関係で、逆変換しても値が元に戻らない場合がありました。しかし、10進数から2進数の変換では、1/5が2進数で表現できないために正確な変換ができないのに対して、2進数から10進数の変換なら誤差のない正確な変換ができるはずです。とりあえず10進数の有効桁数を倍精度浮動小数点の有効桁数より大きくしてみて、パッチをswfmillのMLに送ってみたところ、「netlibにあるdtoa.cg_fmt.cを使うといいよ」というアドバイスをもらいました。

そもそも、IEEE 754 倍精度浮動小数点での仮数は53桁であり、10の16乗 > 2の53乗 です。2進から10進への変換と10進から2進への変換できちんと丸めを行っているのであれば、10進数では16桁前後の有効桁数があれば可逆変換ができるはずです。dtoa.cは実際にそれを実装したものでした。

dtoa.cには、strtod()とdtoa()という関数が実装されています。strtod()が正確に近似方向への丸めを行う文字列からdoubleへの変換で、昨年Ruby 1.8に取り込まれたものです。dtoa()は逆にdoubleから文字列への変換で、引数のmodeに0を渡すことで、strtod()で正確に元の値に戻せる最低限の桁数で10進化してくれます。ただし、dtoa()が作成する文字列には数字しか含まれておらず、小数点の位置は別に返します。g_fmt.cのg_fmt()というラッパ関数を使うことで、dtoa()の結果を元に通常表記と指数表記のうち短くなる表現を選んで文字列にしてくれます。

実際にg_fmt()とstrtod()を利用したサンプルプログラムを置いておきます。下の実行例では、sprintfに%.18lgとした場合に正確に元の値に戻っていますが、g_fmt()を使うと一文字少ない10進数でも元の値に戻っているのが判ります。

#include <cstdio>
#include <cstdlib>

using namespace std;

extern "C" void g_fmt(char *, double d);

// doubleを分解して表示する.
void dumpdouble(double d)
{
    unsigned long long val = *(unsigned long long*)&d;
    unsigned long long fraction = val & ((1ULL << 52) - 1);
    val >>= 52;
    int exp = val & ((1 << 11) - 1);
    exp -= 1023;
    val >>= 11;
    int sign = val ? '-' : '+';

    printf("sign: %c , exp: %5d , 1.", sign, exp);
    for (int i = 51; i >= 0; --i) {
        int c = fraction & (1ULL << i) ? '1' : '0';
        putchar(c);
    }
    putchar('\n');
}

void test(double d)
{
    puts("--- from ---");
    dumpdouble(d);
    puts("---  to  ---");

    char buff[1024];
    double e;

    sprintf(buff, "%lg", d);
    printf("sprintf(%%lg) = %s\n", buff);
    e = strtod(buff, 0);
    dumpdouble(e);

    sprintf(buff, "%.18lg", d);
    printf("sprintf(%%.18lg) = %s\n", buff);
    e = strtod(buff, 0);
    dumpdouble(e);

    g_fmt(buff, d);
    printf("g_fmt() = %s\n", buff);
    e = strtod(buff, 0);
    dumpdouble(e);
}

int main()
{
    unsigned long long val = 0x400FFFFFF0000000ULL;
    double d = *(double*)&val;

    test(d);

    return 0;
}
--- from ---
sign: + , exp:     1 , 1.1111111111111111111111110000000000000000000000000000
---  to  ---
sprintf(%lg) = 4
sign: + , exp:     2 , 1.0000000000000000000000000000000000000000000000000000
sprintf(%.18lg) = 3.99999988079071045
sign: + , exp:     1 , 1.1111111111111111111111110000000000000000000000000000
g_fmt() = 3.9999998807907104
sign: + , exp:     1 , 1.1111111111111111111111110000000000000000000000000000

人間が読み書きできる設定ファイルに実数を保存するときなど、実数を10進数で表現する機会には試してみてください。


@methane
klab_gijutsu2 at 15:32│Comments(1)TrackBack(0)

トラックバックURL

この記事へのコメント

1. Posted by 通りすがり   2008年02月22日 11:37
いつも楽しく読ませてもらってます。

netlib懐かしいなー。
学生の時めちゃくちゃお世話になった。
自分の場合は、fortranだったけどw

久しぶりに漁ってみようかしら・・・

この記事にコメントする

名前:
URL:
  情報を記憶: 評価: 顔   
 
 
 
Blog内検索
Archives
このブログについて
DSASとは、KLab が構築し運用しているコンテンツサービス用のLinuxベースのインフラです。現在5ヶ所のデータセンタにて構築し、運用していますが、我々はDSASをより使いやすく、より安全に、そしてより省力で運用できることを目指して、日々改良に勤しんでいます。
このブログでは、そんな DSAS で使っている技術の紹介や、実験してみた結果の報告、トラブルに巻き込まれた時の経験談など、広く深く、色々な話題を織りまぜて紹介していきたいと思います。
最新コメント