diff_order.py 技術ドキュメント
プログラムの動作
diff_order.py は、指数関数 \(f(x) = \exp(x)\) の導関数を、様々なステップサイズ \(h\) と異なる数値微分近似法を用いて計算し、その精度を評価するPythonスクリプトです。
主な機能は以下の通りです。
関数の定義: 微分対象の関数 \(f(x) = \exp(x)\) と、その解析的な導関数 \(f'(x) = \exp(x)\) を定義します。
数値微分法の実装: 以下の5種類の数値微分近似法を関数として実装しています。
2点前方差分
2点後方差分
3点中心差分
5点中心差分
7点中心差分
誤差評価: 初期ステップサイズ
hから、kh倍ずつ減少させながらnh回にわたってステップサイズを変化させます。各ステップサイズにおいて、定義された各数値微分法で導関数を計算し、その結果と解析解との絶対誤差を求めます。結果の出力:
計算された数値微分値と誤差をコンソールに表示します。
すべての計算結果を
diff_order.csvというCSVファイルに保存します。matplotlibを用いて、ステップサイズhに対する各数値微分値と、各数値微分法の誤差を視覚的に表示するグラフを生成します。
このプログラムは、数値微分におけるステップサイズと近似法の選択が、結果の精度にどのように影響するかを理解し、可視化することを目的としています。
原理
このプログラムは、関数のテイラー展開に基づいた数値微分近似法を利用しています。関数 \(f(x)\) の点 \(x\) における導関数 \(f'(x)\) は、以下のテイラー展開から導き出されます。
これを変形することで、様々な差分近似式が得られます。
2点前方差分 (2-point forward difference): \(f(x+h) = f(x) + hf'(x) + O(h^2)\) から \(f'(x)\) について解きます。
\[ f'(x) \approx \frac{f(x+h) - f(x)}{h} \]この近似は \(O(h)\) のオーダーで、ステップサイズ \(h\) が減少すると線形に誤差が減少します。
2点後方差分 (2-point backward difference): \(f(x-h) = f(x) - hf'(x) + O(h^2)\) から \(f'(x)\) について解きます。
\[ f'(x) \approx \frac{f(x) - f(x-h)}{h} \]この近似も \(O(h)\) のオーダーです。
3点中心差分 (3-point central difference): \(f(x+h) = f(x) + hf'(x) + \frac{h^2}{2!}f''(x) + \frac{h^3}{3!}f'''(x) + O(h^4)\) と \(f(x-h) = f(x) - hf'(x) + \frac{h^2}{2!}f''(x) - \frac{h^3}{3!}f'''(x) + O(h^4)\) を用います。 これらを引くと、偶数次の項が消去され、高次の精度が得られます。
\[ f'(x) \approx \frac{f(x+h) - f(x-h)}{2h} \]この近似は \(O(h^2)\) のオーダーで、2点差分法よりも高精度です。
5点中心差分 (5-point central difference): さらに多くの点を使い、テイラー展開のより高次項を考慮して導出される、より高精度な近似です。
\[ f'(x) \approx \frac{-f(x+2h) + 8f(x+h) - 8f(x-h) + f(x-2h)}{12h} \]この近似は \(O(h^4)\) のオーダーです。
7点中心差分 (7-point central difference): 5点中心差分よりもさらに多くの点を使用し、非常に高精度な近似を提供します。
\[ f'(x) \approx \frac{f(x+3h) - 9f(x+2h) + 45f(x+h) - 45f(x-h) + 9f(x-2h) - f(x-3h)}{60h} \]この近似は \(O(h^6)\) のオーダーです。
プログラムでは、関数 func(x) として \(e^x\) を、解析解 dfdx(x) としても \(e^x\) を用いてこれらの数値微分法を評価し、それぞれの誤差を解析解との絶対差として計算します。
必要な非標準ライブラリとインストール方法
diff_order.py の実行には、以下の非標準ライブラリが必要です。
numpy: 数値計算を効率的に行うためのライブラリです。
matplotlib: グラフ描画のためのライブラリです。
これらのライブラリは、Pythonのパッケージインストーラ pip を使ってインストールできます。コマンドプロンプトやターミナルで以下のコマンドを実行してください。
pip install numpy matplotlib
必要な入力ファイル
diff_order.py は、実行に必要な外部入力ファイルはありません。
すべての計算パラメータ(初期ステップサイズ h、ステップサイズ減少率 kh、計算回数 nh、微分中心点 x0)はスクリプトのソースコード内部で定義されています。
もしこれらのパラメータを変更したい場合は、diff_order.py のソースコードを直接編集する必要があります。
生成される出力ファイル
diff_order.py は、以下のファイルと視覚出力を生成します。
diff_order.csv:このファイルは、各ステップサイズ
hにおける数値微分結果をカンマ区切り形式で保存します。ファイルはプログラム実行ディレクトリに生成されます。
内容は以下の列で構成されます。
Ndiv: ステップサイズの計算回数(1からnhまで)。h: その回のステップサイズ。Df/Dx(2-point): 2点前方差分法による数値微分値。Df/Dx(3-point): 3点中心差分法による数値微分値。Df/Dx(5-point): 5点中心差分法による数値微分値。Df/Dx(7-point): 7点中心差分法による数値微分値。
Matplotlibグラフ:
プログラムの実行が完了すると、
matplotlibによって生成されたグラフウィンドウが表示されます。このウィンドウには2つのサブプロットがあります。
上段: ステップサイズ
hに対する各数値微分法の計算結果と、解析解df/dx (exact)がプロットされます。下段: ステップサイズ
hに対する各数値微分法の誤差(解析解との絶対差)が両対数スケールでプロットされます。これにより、誤差がhに対してどのように減少するかが視覚的に確認できます。
このウィンドウは、ユーザーがターミナルで
Enterキーを押すまで開いたままになります。
コマンドラインでの使用例 (Usage)
diff_order.py はコマンドライン引数を受け付けません。そのため、特別なオプションなしで直接実行します。
基本的な実行コマンドは以下の通りです。
python diff_order.py
プログラムは、ソースコード内で定義されているデフォルトのパラメータ(x0 = 1.0, h = 0.5, kh = 0.8, nh = 101)を使用して計算を実行します。
コマンドラインでの具体的な使用例
diff_order.py を実行すると、まずコンソールに出力が行われ、その後グラフウィンドウが表示されます。
実行コマンド:
python diff_order.py実行結果の説明:
コンソール出力: プログラムはまず、解析値と、各ステップサイズ
hに対する各数値微分法の計算結果を標準出力に表示します。 出力例(一部抜粋):Numerical differentiation using differnet approximations Write to [diff_order.csv] Analytical values: f(1.0)=2.718281828459045 df/dx(1.0)=2.718281828459045 Ndiv: h Df/Dx(2-point) Df/Dx(3-point) Df/Dx(5-point) Df/Dx(7-point) 1 : 5.0000e-01 3.376839353998 2.839849503463 2.718362635445 2.718282710323 2 : 4.0000e-01 3.136890356616 2.771986422204 2.718296316270 2.718281830493 ... 101 : 1.6385e-08 2.718281804797 2.718281828461 2.718281828459 2.718281828459
diff_order.csvファイルの生成: プログラムが実行されたディレクトリにdiff_order.csvという名前のファイルが生成されます。このファイルには、上記のコンソール出力と対応するデータがCSV形式で格納されます。 内容例(一部抜粋):Ndiv,h,Df/Dx(2-point),Df/Dx(3-point),Df/Dx(5-point),Df/Dx(7-point) 1,0.5,3.376839353998135,2.839849503463375,2.718362635444983,2.718282710322521 2,0.4,3.136890356616147,2.771986422203649,2.718296316270388,2.718281830492801 ... 101,1.6384597395029463e-08,2.7182818047968594,2.718281828461159,2.718281828459045,2.718281828459045
Matplotlibグラフの表示: コンソール出力に続いて、
matplotlibによって生成されたグラフウィンドウが表示されます。 このグラフは、上段に各差分法の数値微分値(\(y\)軸)と解析解がステップサイズh(\(x\)軸)に対してプロットされ、下段には各差分法の誤差(\(y\)軸)がステップサイズh(\(x\)軸)に対して両対数プロットされます。 グラフの例(テキストでの説明):上段プロット: \(h\) が大きいときは各差分法の結果が大きく異なり、解析解から離れています。\(h\) が小さくなるにつれて、より高次の差分法(5点、7点)の結果が解析解に急速に近づくことが示されます。
下段プロット: 誤差が両対数スケールでプロットされており、各差分法のオーダー(\(O(h^1)\), \(O(h^2)\), \(O(h^4)\), \(O(h^6)\))に対応する傾きで誤差が減少していく様子が確認できます。ただし、非常に小さい \(h\) では浮動小数点演算の丸め誤差により、誤差が再び増加する傾向が見られることがあります。
グラフウィンドウは、コンソールで
Press ENTER to exit>>と表示された後、ユーザーがEnterキーを押すまで表示され続けます。Enterを押すと、プログラムは終了します。