コンピュータでおかしなことになる計算例 (1)

小清水です。初ブログになります。

カタストロフィックに 計算結果がおかしくなる 例を紹介したいと思います。

Excel を使用した例と C言語の例です。

---

お膳立て：(飛ばしてOK)

コンピュータの中で 実数 は一般にどのように表されているのでしょうか？

浮動小数点数 と呼ばれる数で表されることが多いです。

例えば、 はどうでしょうか？

この円周率は コンピュータ内部において、

という 浮動小数点数 で近似されます。(単精度を使用した例)

2進法ではわかりづらいため、この数を 10進法 でかいてみると

3.141592 7410125732 という数になります。

本来の値は 3.141592 65358979 … なので 誤差 が生じたことになります。

この誤差を 丸め誤差 というわけです。

丸め誤差って何なのか、詳しく知りたい方はこちらをどうぞ！

math-koshimizu.hatenablog.jp

---

本題！おかしな例

今回の話は、 コンピュータで近似計算をしているため、 おかしな結果がでるよという話です。

とはいえ、色んな所で使われている Excel でもそんなことできるんでしょうか？

Excel による例

Excel 2011 (Mac) を使って次の式を計算したいと思います：

図の通り、 1.17260394 という計算ができました！

しかしながら実際の -0.827396... です。

もはや近いというレベルではありません。符号も違います。

(追記)

この -0.827396... の導出方法を解説しました。

math-koshimizu.hatenablog.jp

C言語による例

プログラマー向けに C言語によるプログラムも。よく使われる double で計算してみます。

#include <stdio.h>

double rump1(double a, double b)
{
    return 333.75 * b*b*b*b*b*b 
        + a*a * (11.0 * a* a* b* b - b*b*b*b*b*b - 121.0 * b*b*b*b - 2.0) 
        + 5.5 * b*b*b*b*b*b*b*b + a /(2.0*b);
}

double rump2(double a, double b)
{
    double b2 = b*b;
    double b4 = b2*b2;
    double b6 = b4*b2;

    return 333.75 * b*b*b*b*b*b
        + a*a * (11.0 * a*a* b*b - b6 - 121.0 * b*b*b*b - 2.0)
        + 5.5  *b*b*b*b*b*b*b*b + a/(2.0*b);
}

double rump3(double a, double b)
{
    double b2 = b*b;
    double b4 = b2*b2;
    double b6 = b4*b2;
    double b8 = b4*b4;
    double a2 = a*a;

    return 333.75 * b6
        + a2 * (11.0 * a2* b2 - b6 - 121.0 * b4 - 2.0)
        + 5.5  *b8 + a/(2.0*b);
}

int main(int argc, char** argv)
{
    const double a = 77617.0;
    const double b = 33096.0;

    printf("%f\n", rump1(a, b)); /* ---->  1.172604 */
    printf("%f\n", rump2(a, b)); /* ----> -2361183241434822606848.000000 */
    printf("%f\n", rump3(a, b)); /* ----> -1180591620717411303424.000000 */

    return 0;
}

と、計算の仕方によってはさっきの Excel のようになったり、 激しく大きな負数になったりします。

困りましたね。

端折りますが、単精度・倍精度・4倍精度でも計算は合いません。

---

まさにカタストロフィック！

この元ネタは 「Algorithms for verified inclusion」という1988年の論文です。 興味がある方は是非御覧ください。

この記事の魅力としては

Excel を使ってもC言語を使っても 1.17260394 という計算結果を得られたけど、本当の答えは-0.827396...という圧倒的に変な結果が得られた というわけで

世の中の人はExcelに対してとても不安になったに違いない！

浮動小数点数って難しい…って思ってもらうための話でした。

第二弾もあるよ。いつか書きます！

第二弾も書いたよ！

math-koshimizu.hatenablog.jp