前進差分
この章の目標
前進差分演算子 $\Delta$ の定義と性質、前進差分表の構成、ニュートン前進差分補間公式を理解する。
目次
はじめに:なぜ前進差分か
導関数 $f'(x)$ を知りたいのに、式で微分できない場面は多い。測定データのように飛び飛びの値しか無い場合や、関数がプログラムの中でしか計算できない場合である。そのようなときに使える最も単純な道具が前進差分である。隣り合う 2 点の値の差を取るだけで、変化の大きさを見積もれる。
差を取る操作は、数列の規則を見つけるときにも役立つ。$1, 4, 9, 16, 25$ の隣り合う差は $3, 5, 7, 9$、その差は $2, 2, 2$ と一定になる。2 回目で一定になることから、元の数列が 2 次式($n^2$)で表せることが分かる。この章では、この「差を取る」操作を体系的に扱う。
1. 定義
関数 $f(x)$ に対して、前進差分(forward difference)は次のように定義される。
前進差分の定義
$$\Delta f(x) = f(x + h) - f(x)$$
ここで $h > 0$ は差分幅(step size)である。
直感的には、$\Delta f(x)$ は区間 $[x,\; x+h]$ における $f$ の「変化量」を表す。これは微分における $f'(x)\,h$ の離散的な対応物である。
2. 高階差分
差分演算を繰り返し適用することで高階差分が得られる。
$n$ 階前進差分
$$\begin{aligned} \Delta^n f(x) &= \sum_{k=0}^{n} (-1)^{n-k} \\ &\qquad \times \binom{n}{k} f(x + kh) \end{aligned}$$
最初のいくつかの高階差分を具体的に書き下すと以下のようになる。
- $\Delta^0 f(x) = f(x)$(恒等作用素)
- $\Delta^1 f(x) = f(x+h) - f(x)$
- $\Delta^2 f(x) = f(x+2h) - 2f(x+h) + f(x)$
- $\Delta^3 f(x) = f(x+3h) - 3f(x+2h)$ $+\, 3f(x+h) - f(x)$
二項係数が交互符号で現れる構造に注目されたい。$n$ 階差分は $n+1$ 個の関数値を使って計算される。
差分表(Difference Table)
等間隔点 $x_i = x_0 + ih$ での関数値から差分表を構成できる。
| $i$ | $x_i$ | $f_i$ | $\Delta f_i$ | $\Delta^2 f_i$ | $\Delta^3 f_i$ |
|---|---|---|---|---|---|
| 0 | $x_0$ | $f_0$ | |||
| $\Delta f_0$ | |||||
| 1 | $x_1$ | $f_1$ | $\Delta^2 f_0$ | ||
| $\Delta f_1$ | $\Delta^3 f_0$ | ||||
| 2 | $x_2$ | $f_2$ | $\Delta^2 f_1$ | ||
| $\Delta f_2$ | |||||
| 3 | $x_3$ | $f_3$ |
各列は左隣の列の隣接する2値の差として計算される。$\Delta f_i = f_{i+1} - f_i$、$\Delta^2 f_i = \Delta f_{i+1} - \Delta f_i$ である。
3. 差分作用素
前進差分を形式的な作用素として扱うと、演算間の関係が明確になる。
主要な作用素
- シフト作用素 $E$: $Ef(x) = f(x+h)$
- 前進差分作用素 $\Delta$: $\Delta f(x) = f(x+h) - f(x)$
- 恒等作用素 $I$: $If(x) = f(x)$
- 微分作用素 $D$: $Df(x) = f'(x)$
作用素間の基本関係
$$\Delta = E - I$$
すなわち $E = I + \Delta$ である。また形式的べき級数展開により
$$E = e^{hD}$$
が成り立つ。これらを組み合わせると
$$\Delta = e^{hD} - I$$
が得られ、差分と微分の橋渡しが可能になる。
$n$ 階差分は作用素の $n$ 乗として $\Delta^n = (E - I)^n$ と表され、二項展開すると
$$\Delta^n = \displaystyle\sum_{k=0}^{n} (-1)^{n-k} \binom{n}{k} E^k$$
となり、前節の二項公式が再導出される。
4. Newton の前進差分補間公式
等間隔データ点 $x_i = x_0 + ih$($i = 0, 1, \ldots, n$)における関数値 $f_i$ が与えられたとき、Newton の前進差分補間公式は次のように表される。
Newton の前進差分補間公式
$s = (x - x_0)/h$ とおくと
$$P_n(x) = \displaystyle\sum_{k=0}^{n} \binom{s}{k} \Delta^k f_0$$
ここで $\displaystyle\binom{s}{k} = \dfrac{s(s-1)(s-2)\cdots(s-k+1)}{k!}$ は一般化二項係数である。
展開すると
$$\begin{aligned} P_n(x) = {} & f_0 + s\,\Delta f_0 \\ & + \dfrac{s(s-1)}{2!}\,\Delta^2 f_0 \\ & + \dfrac{s(s-1)(s-2)}{3!}\,\Delta^3 f_0 + \cdots \end{aligned}$$
となる。差分表の最上行($\Delta^k f_0$)だけで補間多項式が構成される点が特徴である。
計算例
$x_0 = 0,\; h = 1$ で $f(0)=1,\; f(1)=2,\; f(2)=5,\; f(3)=10$ が与えられたとする。差分表は
$\Delta f_0 = 1,\; \Delta f_1 = 3,\; \Delta f_2 = 5$
$\Delta^2 f_0 = 2,\; \Delta^2 f_1 = 2$
$\Delta^3 f_0 = 0$
よって
$$\begin{aligned} P_3(x) &= 1 + x \cdot 1 \\ &\quad + \dfrac{x(x-1)}{2}\cdot 2 \\ &\quad + \dfrac{x(x-1)(x-2)}{6}\cdot 0 \\ &= 1 + x + x(x-1) \\ &= x^2 + 1 \end{aligned}$$
3 階差分が 0 であるため、データは 2 次多項式 $x^2 + 1$ で完全に表現される。
5. 導関数との関係
前進差分と導関数の間には、$h \to 0$ の極限で次の基本的な関係がある。
前進差分商と導関数
$f$ が微分可能であるとき
$$\begin{aligned} \lim_{h \to +0} \dfrac{\Delta f(x)}{h} &= \lim_{h \to +0} \dfrac{f(x+h) - f(x)}{h} \\ &= f'(x) \end{aligned}$$
($h > 0$ と定義したので右側からの極限をとる。$f$ が微分可能なら、その値は $f'(x)$ に一致する。)
これは微分の定義そのものであり、前進差分商 $\Delta f(x)/h$ が導関数 $f'(x)$ の1次近似であることを示している。
高階についても同様に
$$\dfrac{\Delta^n f(x)}{h^n} \approx f^{(n)}(x)$$
が成り立つ。作用素の言葉では $\Delta/h \to D$($h \to 0$)に対応する。
6. 打ち切り誤差
$f$ が十分に滑らかであるとき、Taylor 展開により前進差分商の誤差を評価できる。
$$\begin{aligned} f(x+h) ={} & f(x) + hf'(x) \\ & + \dfrac{h^2}{2}f''(x) + O(h^3) \end{aligned}$$
より
$$\begin{aligned} \dfrac{\Delta f(x)}{h} &= \dfrac{f(x+h) - f(x)}{h} \\ &= f'(x) + \dfrac{h}{2}f''(x) + O(h^2) \end{aligned}$$
打ち切り誤差
$$\begin{gathered} \dfrac{\Delta f(x)}{h} - f'(x) = \dfrac{h}{2}f''(\xi), \\ \xi \in (x,\; x+h) \end{gathered}$$
打ち切り誤差は $O(h)$ である。これは1次精度(first-order accurate)と呼ばれる。
比較のため、中心差分の打ち切り誤差は $O(h^2)$(2次精度)であり、前進差分より一段高い精度を持つ。しかし前進差分は片側の点だけで計算できるので、領域の左端のように片側にしか点が無い場所でも使える。また時間発展の問題では、現在の値から次の時刻の値を直接計算する陽的解法(前進 Euler 法など)の土台になる。
丸め誤差との兼ね合い
$h$ を小さくすると打ち切り誤差は減少するが、$f(x+h)$ と $f(x)$ の差が小さくなることで桁落ち(cancellation error)が増大する。浮動小数点数の丸め単位を $\varepsilon_{\mathrm{mach}}$ とすると、丸め誤差は $O(\varepsilon_{\mathrm{mach}}/h)$ 程度であり、総合誤差を最小化する最適な $h$ は
$$h_{\mathrm{opt}} \sim \sqrt{\varepsilon_{\mathrm{mach}}}$$
のオーダーとなる。倍精度($\varepsilon_{\mathrm{mach}} \approx 10^{-16}$)では $h \approx 10^{-8}$ が目安である。
注意:$h$ は小さいほど正確、ではない
打ち切り誤差が $O(h)$ なので「$h$ を小さくするほど正確になる」と思いたくなるが、実は違う。$f(x) = e^x$ の $x = 1$ での導関数 $e$ を前進差分商で倍精度計算したときの誤差は次のとおりである。
| $h$ | $10^{-1}$ | $10^{-4}$ | $10^{-8}$ | $10^{-12}$ | $10^{-14}$ |
|---|---|---|---|---|---|
| 誤差 | $1.4 \times 10^{-1}$ | $1.4 \times 10^{-4}$ | $6.6 \times 10^{-9}$ | $4.3 \times 10^{-4}$ | $9.3 \times 10^{-3}$ |
$h = 10^{-8}$ 付近で誤差が最小になり、それより小さくすると桁落ちで誤差が増える。$h = 10^{-12}$ では $h = 10^{-4}$ のときより悪い。
7. 応用
7.1 数値微分
解析的に導関数を求めることが困難な場合、前進差分商 $\Delta f/h$ は最も単純な数値微分の手段である。精度が不足する場合は中心差分やRichardson の外挿法と組み合わせることが多い。
7.2 有限差分法(PDE の数値解法)
偏微分方程式の数値解法において、前進差分は時間方向の離散化に広く用いられる。
FTCS スキーム(熱方程式)
1次元熱方程式 $\dfrac{\partial u}{\partial t} = \alpha \dfrac{\partial^2 u}{\partial x^2}$ に対して、時間微分を前進差分、空間2階微分を中心差分で近似すると
$$\dfrac{u_i^{n+1} - u_i^n}{\Delta t} = \alpha \dfrac{u_{i+1}^n - 2u_i^n + u_{i-1}^n}{(\Delta x)^2}$$
が得られる(FTCS: Forward Time, Central Space)。この陽的スキームは $\alpha\,\Delta t / (\Delta x)^2 \le 1/2$ のとき安定である。
7.3 数列と差分方程式
離散系での差分方程式 $\Delta y_n = g(n, y_n)$ は、連続系の微分方程式 $y' = g(x, y)$ の離散的対応物である。再帰式 $y_{n+1} = y_n + g(n, y_n)$ として逐次計算でき、Euler 法との対応関係は明確である。
7.4 組合せ論・多項式の次数判定
$f$ が $n$ 次多項式であるとき、$\Delta^n f$ は定数となり $\Delta^{n+1} f = 0$ である。したがって差分表の列が 0 になる段階から元の関数の多項式次数を判定できる。この性質は組合せ論における数列の一般項の推定にも応用される。
8. 練習問題
問1(基礎):$f(x) = x^3$、$h = 1$ として $x = 0, 1, 2, 3, 4$ の差分表を作り、$\Delta^3 f$ が定数になることを確かめよ。その定数の値を求めよ。
ヒント
関数値は $0, 1, 8, 27, 64$ である。差を 3 回取る。7.4 節より、$n$ 次多項式の $n$ 階差分は定数になる。
問2(標準):$f(x) = e^x$ の $x = 0$ での導関数を、$h = 0.1$ の前進差分商で近似せよ。誤差が 6 節の式 $\dfrac{h}{2}f''(\xi)$($0 < \xi < 0.1$)の範囲に入ることを確かめよ。
ヒント
$\dfrac{e^{0.1} - 1}{0.1} \approx 1.0517$ である。真の値は $f'(0) = 1$。$f''(\xi) = e^{\xi}$ は $1$ と $e^{0.1}$ の間にある。
問3(応用):Taylor 展開を使って $\dfrac{\Delta^2 f(x)}{h^2}$ と $f''(x)$ の差を求め、誤差の次数を答えよ。
ヒント
$f(x+2h)$ と $f(x+h)$ をそれぞれ $h^3$ の項まで展開し、$\Delta^2 f(x) = f(x+2h) - 2f(x+h) + f(x)$ に代入する。$h^2$ で割ったあとに残る最初の項に注目する。
9. よくある質問
Q1. 前進差分と中心差分の違いは何か
前進差分 $\Delta f(x) = f(x+h) - f(x)$ は片側の情報のみを使い、打ち切り誤差が $O(h)$ である。一方、中心差分 $\delta f(x) = f(x+h/2) - f(x-h/2)$ は両側の情報を使い、打ち切り誤差が $O(h^2)$ と高精度である。ただし前進差分は境界点でも適用でき、実装が単純であるという利点がある。
Q2. Newton の前進差分補間公式はどのような場面で使うか
等間隔に配置されたデータ点が与えられたとき、差分表を用いて多項式補間を行う手法である。特にデータ列の先頭付近(始端側)での補間精度が高く、数表の補間や数値積分の基礎として広く使われている。
Q3. 前進差分は偏微分方程式の数値解法でどのように使われるか
偏微分方程式の時間方向の離散化に前進差分を使う手法が陽的有限差分法である。例えば熱方程式の時間微分を前進差分で、空間微分を中心差分で近似する FTCS スキームが代表例である。実装は単純だが、安定のためには $\alpha\,\Delta t / (\Delta x)^2 \le 1/2$ を満たす必要がある。
10. 参考資料
- Wikipedia「差分」(日本語版)
- Wikipedia「Finite difference」(英語版)
- Wikipedia「Newton's forward difference formula」(英語版)
- R. L. Burden & J. D. Faires, Numerical Analysis, 10th ed., Cengage, 2016.
- K. E. Atkinson, An Introduction to Numerical Analysis, 2nd ed., Wiley, 1989.
- C. Jordan, Calculus of Finite Differences, 3rd ed., Chelsea, 1965.