異なる浮動小数点精度における累積誤差の比較

プログラムの動作

sum_error.py は、異なる浮動小数点精度（float16, float32, float64）における累積誤差を比較するためのPythonプログラムです。指定された浮動小数点数 h を n 回繰り返し加算し、各精度で計算された合計値と真値との差（累積誤差）を評価します。

プログラムの主な機能は以下の通りです。

精度別の累積加算: NumPyライブラリの np.float16, np.float32, np.float64 データ型を使用して、同じ加算処理を異なる精度で実行します。
誤差の計算: 各ステップにおいて、正確な値（数学的に計算された真値）と、各精度で計算された累積和との差を誤差として算出します。誤差の計算自体は、より高い精度である np.float64 で行われます。
結果の出力:
- 計算過程で指定されたステップ間隔 iPrintStep ごとに、現在のステップの正確な値、各精度での累積和、および対応する誤差を標準出力に整形して表示します。
- すべての計算ステップにおける正確な値、各精度での累積和、および誤差を sum_error.csv というCSVファイルに保存します。

このプログラムは、浮動小数点演算における丸め誤差が、データ型の精度によってどのように蓄積されるかを具体的に示し、数値計算の信頼性を理解する上で役立ちます。

原理

浮動小数点数は、コンピュータ内部で有限のビット数で表現されます。このため、多くの実数は正確に表現できず、近似値として扱われます。この近似によって生じるわずかな誤差を丸め誤差と呼びます。加算などの演算を繰り返すと、各ステップで発生する小さな丸め誤差が累積され、最終的な結果に大きな影響を与えることがあります。

本プログラムでは、以下の原理に基づいて累積誤差を評価します。

真値の計算: プログラムでは、h を i+1 回加算した「真値」（exact value）を、より高い精度（通常はシステムがサポートする最高の精度、ここではfloat64）で計算します。これは、数学的には $E_i = (i+1) \cdot h$ と表されます。
精度別累積和の計算:
- 各浮動小数点精度 p （float16, float32, float64）について、初期値 $S_{p,0} = 0$ とし、ステップごとに $S_{p,i} = S_{p,i-1} + h$ の計算を繰り返します。ただし、この加算は各精度 p のデータ型で行われます。
- つまり、各ステップでの加算 $S_{p,i-1} + h$ の結果は、その精度 p の範囲内で最も近い値に丸められます。
誤差の評価: 各ステップ i において、真値 $E_i$ と、各精度で計算された累積和 $S_{p,i}$ との差を誤差として定義します。 $$Err_{p,i} = E_i - S_{p,i}$$ この誤差 $Err_{p,i}$ は、より高い精度（float64）で計算されます。

このプロセスを通じて、精度の低い float16 は float32 や float64 に比べてより大きな丸め誤差を生じやすく、その誤差が繰り返し加算によって顕著に蓄積されることを観察できます。逆に、float64 のように高精度なデータ型では、丸め誤差が小さく、累積誤差も抑制されることが期待されます。

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

このプログラムは、数値計算を効率的に行うために NumPy ライブラリを使用します。

ライブラリ名: numpy

インストールは、Pythonのパッケージマネージャ pip を使用してコマンドラインから実行できます。

pip install numpy

必要な入力ファイル

このプログラムは、外部の入力ファイルを必要としません。必要なすべてのパラメータ（加算する値 h、繰り返し回数 n、表示ステップ間隔 iPrintStep）は、コマンドライン引数として直接指定されるか、プログラム内部のデフォルト値が使用されます。

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

プログラムは、sum_error.py と同じディレクトリに sum_error.csv というCSVファイルを出力します。

ファイル名: sum_error.csv
内容: このファイルには、各計算ステップにおける以下のデータがカンマ区切りで保存されます。
1. exact: 現在のステップにおける真値（正確な累積和）。
2. float16: np.float16 型で計算された累積和。
3. float32: np.float32 型で計算された累積和。
4. float64: np.float64 型で計算された累積和。
5. error(float16): exact と float16 での累積和との差。
6. error(float32): exact と float32 での累積和との差。
7. error(float64): exact と float64 での累積和との差。

ファイルは以下のヘッダー行から始まります。 exact,float16,float32,float64,error(float16),error(float32),error(float64)

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

プログラムは以下の形式でコマンドラインから実行します。

python sum_error.py h n iPrintStep

h: 加算する浮動小数点数。
- 例: 0.01
- 省略した場合のデフォルト値: 0.01
n: 加算を繰り返す回数（整数）。
- 例: 1000
- 省略した場合のデフォルト値: 101
iPrintStep: 結果を標準出力に表示するステップの間隔（整数）。
- 例: 100
- 省略した場合のデフォルト値: 10

引数を省略した場合、プログラムはデフォルト値を使用します。引数を3つ未満で実行した場合、Usage メッセージが表示されます。

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

h を 0.01、繰り返し回数 n を 100、表示ステップ間隔 iPrintStep を 10 としてプログラムを実行します。

python sum_error.py 0.01 100 10

このコマンドを実行すると、以下のような出力が標準出力に表示されます（実際の数値は環境によって多少異なる場合があります）。

Summing up 0.01 for 100 times with different precision floating point types
Write to [sum_error.csv]

exact:	sum16 (error)               	sum32 (error)               	sum64 (error)               

0000: 0.000000000000000000 (+0.00e+00)  0.000000000000000000 (+0.00e+00)  0.000000000000000000 (+0.00e+00)  

1000: 0.100000001490116120 (-1.49e-08)  0.100000000000000010 (+1.39e-17)  0.100000000000000010 (+1.39e-17)  

2000: 0.200000002980232240 (-2.98e-08)  0.199999999999999980 (-2.78e-17)  0.200000000000000010 (+2.78e-17)  

3000: 0.300000011920928960 (-1.19e-08)  0.300000000000000040 (+4.16e-17)  0.300000000000000040 (+4.16e-17)  

4000: 0.400000005960464500 (-5.96e-09)  0.400000000000000020 (+2.78e-17)  0.400000000000000020 (+2.78e-17)  

5000: 0.500000000000000000 (+0.00e+00)  0.500000000000000000 (+0.00e+00)  0.500000000000000000 (+0.00e+00)  

6000: 0.600000023841857900 (-2.38e-08)  0.599999999999999980 (-2.78e-17)  0.599999999999999980 (-2.78e-17)  

7000: 0.699999988079071000 (+1.19e-08)  0.699999999999999960 (-5.55e-17)  0.699999999999999960 (-5.55e-17)  

8000: 0.800000011920929000 (-1.19e-08)  0.800000000000000040 (+5.55e-17)  0.800000000000000040 (+5.55e-17)  

9000: 0.900000035762786900 (-3.58e-08)  0.900000000000000020 (+2.78e-17)  0.900000000000000020 (+2.78e-17)  

0000: 1.000000000000000000 (+0.00e+00)  1.000000000000000000 (+0.00e+00)  1.000000000000000000 (+0.00e+00)

実行結果の説明:

Summing up 0.01 for 100 times with different precision floating point types: プログラムの実行設定を表示します。
Write to [sum_error.csv]: 結果が sum_error.csv に保存されることを示します。
exact: 左端の列は、数学的に正確な累積和 $(i+1) \cdot 0.01$ を示します。
sum16 (error): float16 型で計算された累積和と、それに対応する誤差（真値との差）が括弧内に指数表記で表示されます。float16 は精度が最も低いため、累積誤差が比較的大きくなり、e-08 オーダーで見られます。
sum32 (error): float32 型で計算された累積和と誤差。float16 よりは精度が高く、誤差は e-17 から e-16 オーダーと小さいですが、依然として誤差が存在します。
sum64 (error): float64 型で計算された累積和と誤差。最も精度が高いため、誤差は非常に小さく、e-17 から e-16 オーダーと、float32 と同等か、わずかに小さい場合もあります。

この出力例では、特に float16 の誤差が他の精度に比べて大きいことが明確に示され、浮動小数点精度が累積誤差に与える影響を視覚的に確認できます。