information_buried.py 技術ドキュメント

プログラムの動作

プログラム information_buried.py は、浮動小数点演算における「情報の埋没 (information buried)」という数値誤差の現象をデモンストレーションするPythonスクリプトです。

具体的には、大きな負の値 \(x\)（デフォルトでは \(-40.0\)）に対する指数関数 \(e^x\) の値を計算する際に、以下の2つの異なるテイラー級数展開の方法を比較します。

直接 \(e^x = \sum_{n=0}^{\infty} \frac{x^n}{n!}\) を計算する方法。
\(e^x = \frac{1}{e^{-x}} = \frac{1}{\sum_{n=0}^{\infty} \frac{(-x)^n}{n!}}\) を計算する方法。

これらの計算結果を、math.exp(x) で得られる正確な値と比較しながら、級数の項数 \(n\) の増加とともに各計算方法の精度がどのように変化するかを、標準出力に表示します。

このプログラムは、数値計算において同じ数学的結果を得る複数の方法が存在する場合でも、数値的安定性（特に桁落ちや情報の埋没）の観点から最適な方法を選択することの重要性を示しています。

原理

本プログラムは、指数関数 \(e^x\) のテイラー級数展開に基づいています。

\[e^x = \sum_{n=0}^{\infty} \frac{x^n}{n!} = 1 + x + \frac{x^2}{2!} + \frac{x^3}{3!} + \dots\]

大きな負の数 \(x\)（例えば \(x = -40\)）の場合、直接この級数を用いると、項 \(x^n\) は正と負の非常に大きな値が交互に現れます。例えば、\(x = -40\) のとき、

\(x^2 = (-40)^2 = 1600\)
\(x^3 = (-40)^3 = -64000\)
\(x^4 = (-40)^4 = 2560000\) といったように、絶対値が非常に大きい正負の項が交互に加算されることになります。最終的に \(e^x\) は非常に小さな正の値（例: \(e^{-40} \approx 4.25 \times 10^{-18}\)）になるため、これらの大きな正負の項を加算する過程で、有効桁が失われる「桁落ち」または「情報の埋没」と呼ばれる現象が発生しやすくなります。

これを回避するため、プログラムでは以下の等価な式も評価しています。

\[e^x = \frac{1}{e^{-x}}\]

ここで \(-x\) は正の大きな数となるため、\(e^{-x} = \sum_{n=0}^{\infty} \frac{(-x)^n}{n!}\) の級数展開では、すべての項 \(\frac{(-x)^n}{n!}\) が正の値となります。例えば \(x = -40\) のとき、\(-x = 40\) なので、

\((-x)^2 = 40^2 = 1600\)
\((-x)^3 = 40^3 = 64000\)
\((-x)^4 = 40^4 = 2560000\) といったように、すべての項が正であり、大きな正の値を小さな正の値に加算していく形になります。この加算では桁落ちが発生しにくく、より正確な \(e^{-x}\) の値が得られます。

したがって、1.0 / s2 の計算（ここで s2 は \(e^{-x}\) の級数和）は、s1 の直接的な計算（ここで s1 は \(e^x\) の級数和）よりも、x が大きな負の値の場合に \(e^x\) のより正確な近似値を与えることを示しています。

必要な非標準ライブラリとインストール方法

このプログラムは、Pythonの標準ライブラリである math モジュール（exp 関数と factorial 関数）のみを使用しています。したがって、追加の非標準ライブラリをインストールする必要はありません。

必要な入力ファイル

このプログラムは、外部の入力ファイルを必要としません。計算に使用されるパラメータ x および nmax は、プログラムのソースコード内に直接ハードコードされています。

生成される出力ファイル

このプログラムは、ファイルを生成しません。すべての結果は標準出力（コンソール）に直接表示されます。

コマンドラインでの使用例 (Usage)

このプログラムは、コマンドライン引数を受け付けません。設定されたパラメータ（x と nmax）を使用して、直接実行されます。

基本使用法:

python information_buried.py

コマンドラインでの具体的な使用例

現在の x = -40.0 と nmax = 120 の設定でプログラムを実行する例を以下に示します。

python information_buried.py

実行結果の冒頭部分と抜粋（実行環境によって値はわずかに異なる可能性があります）：

Error due to information buried: exp(-x) for large x

exact: exp(-40.0) = 4.248354255018654e-18
{'n':^3}:	{'sum(x^n/n!)':^18}	{'1 / sum((-x)^n/n!)':^18}
                 1.0	                 1.0	
               -39.0	  0.024390243902439025	
               761.0	 0.00012987012987012986	
             -9806.0	 2.607411684496245e-07	
            89906.67	 4.316830589139268e-09	
           -640153.8	 5.753360447385566e-11	
...
 2.4533161e-15	 4.2483542e-18	
 2.4533161e-15	 4.2483542e-18	
 2.4533161e-15	 4.2483542e-18	
 2.4533161e-15	 4.2483542e-18	
 2.4533161e-15	 4.2483542e-18	
 2.4533161e-15	 4.2483542e-18	
 2.4533161e-15	 4.2483542e-18	
 2.4533161e-15	 4.2483542e-18	
 2.4533161e-15	 4.2483542e-18	
 2.4533161e-15	 4.2483542e-18	

実行結果の説明:

プログラムはまず、x = -40.0 に対する exp(x) の正確な値 (exact) を表示します。この値は約 \(4.248 \times 10^{-18}\) と非常に小さいです。
その後、ループカウンタ n、直接級数展開による部分和を示す {'sum(x^n/n!)':^18} 列、および逆数級数展開による部分和の逆数を示す {'1 / sum((-x)^n/n!)':^18} 列の3つのデータが表示されます。
{'sum(x^n/n!)':^18} の列では、n の値が小さいときには exact の値から大きく乖離し、正負に大きく変動していることが確認できます。これは、大きな正負の項が交互に加算され、その際に桁落ち（情報の埋没）が発生しているためです。n が増えても、この列の値は exact の値には収束せず、最終的に \(2.45 \times 10^{-15}\) 程度の値に落ち着きますが、これは正確な値と比較して約1000倍ものオーダーの誤差があります。
一方、{'1 / sum((-x)^n/n!)':^18} の列では、n が増えるにつれて徐々に正確な値 \(4.2483542 \times 10^{-18}\) に近づいていくことがわかります。最終的には、ほぼ正確な値に近い精度で収束していることが確認できます。
この結果は、大きな負の引数 \(x\) に対する \(e^x\) の計算において、直接的な級数展開 \(e^x = \sum \frac{x^n}{n!}\) は数値的に不安定であるのに対し、\(e^x = \frac{1}{\sum \frac{(-x)^n}{n!}}\) の形式を用いる方が、情報の埋没による桁落ちを避け、より高い精度で計算できることを明確に示しています。