数値微分のステップ幅hの影響評価

プログラムの動作

diff_h.py は、数値微分におけるステップ幅 \(h\) の影響を評価するためのPythonスクリプトです。特定の関数 \(f(x)\) に対して、異なるステップ幅 \(h\) を用いた2点前方差分法による数値微分の精度を検証します。

主な機能は以下の通りです。

  • 関数の定義: 微分対象の関数 \(f(x) = x^3\) と、その厳密な導関数 \(f'(x) = 3x^2\) をプログラム内で定義します。

  • 数値微分の計算: さまざまなステップ幅 \(h\) (例: \(1.0, 0.1, \dots, 1.0 \times 10^{-6}\)) を用いて、2点前方差分法により数値微分値を計算します。

  • 結果の出力: 計算された \(x\) の値、関数値 \(f(x)\)、厳密な導関数 \(f'(x)\)、および各ステップ幅における数値微分値をCSVファイルに保存し、同時にコンソールにも表形式で表示します。

このプログラムは、数値微分の精度がステップ幅 \(h\) によってどのように変化するか、特に \(h\) が大きすぎると近似が粗くなり、小さすぎると浮動小数点演算の丸め誤差の影響を受けやすくなる現象を視覚的に理解するのに役立ちます。

原理

このプログラムでは、関数 \(f(x) = x^3\) の微分を評価します。この関数の厳密な導関数は \(f'(x) = 3x^2\) です。 数値微分には「2点前方差分法 (Two-point forward difference method)」が用いられます。これは、点 \(x\) における関数の導関数を、ステップ幅 \(h\) を用いて以下の近似式で計算する方法です。

\[f'(x) \approx \frac{f(x+h) - f(x)}{h}\]

この式は、点 \((x, f(x))\) と点 \((x+h, f(x+h))\) を結ぶ直線の傾きが、点 \(x\) における接線の傾き(導関数)の近似となることを利用しています。ステップ幅 \(h\) が十分に小さい場合、この近似は高精度になりますが、無限に小さくすると浮動小数点演算の精度限界により誤差が蓄積する可能性があります。

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

このプログラムの実行には、以下の非標準ライブラリが必要です。

  • numpy: 数値計算を効率的に行うためのライブラリ。

インストールは、以下の pip コマンドで実行できます。

pip install numpy

必要な入力ファイル

このプログラムは、微分対象の関数とその導関数をコード内で直接定義しており、外部からの入力ファイルを必要としません。

生成される出力ファイル

プログラムを実行すると、以下のファイルが生成されます。

  • diff_h.csv:

    • 内容: 計算された \(x\) の値、関数値 \(f(x)\)、厳密な導関数 \(df/dx\)、そして各ステップ幅 \(h\) で計算された数値微分値がカンマ区切り形式で保存されます。

    • ヘッダー行: 以下の形式で出力されます。

      x,f(x),df/dx,Df/Dx(h=1.0),Df/Dx(h=0.1),Df/Dx(h=0.01),Df/Dx(h=0.001),Df/Dx(h=1.0e-6)
      
    • データ行: 各行は、対応する \(x\) の点における上記ヘッダーの各項目に対応する数値を含みます。

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

diff_h.py は、引数を指定せずに実行します。

python diff_h.py

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

以下のコマンドでプログラムを実行します。

python diff_h.py

実行すると、以下のようなメッセージがコンソールに出力され、計算結果の表が表示されます。

Numerical differentiation using two-point forward difference method
Write to [diff_h.csv]

  x  :   f(x)    df/dx         h=1.0         h=0.1        h=0.01       h=0.001     h=1.0e-06
0.000:  0.000000    0.000000    1.00000000    0.01000000    0.00010000    0.00000100    0.00000000
0.100:  0.001000    0.030000    0.33100000    0.03310000    0.03030000    0.03003000    0.03000000
0.200:  0.008000    0.120000    0.78100000    0.12610000    0.12060000    0.12006000    0.12000000
0.300:  0.027000    0.270000    1.33100000    0.28210000    0.27330000    0.27033000    0.27000000
... (以下略)

この出力は、x の値、対応する \(f(x)\)、厳密な導関数 \(df/dx\)、そして異なるステップ幅 \(h\) で計算された数値微分値を並べています。 h=1.0 のようにステップ幅が大きい場合、数値微分値は厳密な \(df/dx\) から大きくずれています。 h が小さくなるにつれて (h=0.1, h=0.01, h=0.001)、数値微分値は厳密な \(df/dx\) に近づいていることが確認できます。 しかし、h=1.0e-6 のように非常に小さくなると、x=0.000 の行に見られるように、丸め誤差の影響で厳密解とわずかに異なる場合(この場合は厳密解が0なので誤差が目立つ)もあります。これは、浮動小数点数の精度限界による一般的な現象です。

また、同じディレクトリ内に diff_h.csv という名前のCSVファイルが生成され、このコンソール出力と同じデータが保存されます。