📂数値解析

シンプソンの公式

定義

$f : [a,b] \to \mathbb{R}$が$[a,b]$で積分可能で、$[a,b]$を$\displaystyle h:= {{b-a} \over {n}}$の間隔で等しく分割した$a = x_{0} < \cdots < x_{n} = b$みたいなノードポイントで割ったとする。以下で定義された数値積分オペレータ$I_{n}^{2}$をシンプソンのルールという。 $$ I_{n}^{2} (f) := \sum_{k=1}^{n/2} {{h} \over {3}} \left[ f(x_{2k-2}) + 4 f( x_{2k-1} ) + f(x_{2k} ) \right] $$

定理

$f \in C^4 [a,b]$とする。シンプソンのルールのエラー$E_{1}^{2}$と漸近エラー$\tilde{E}_{n}^{2}$は以下の通り。

• [1]: $$E_{1}^{2} (f) = - {{h^5} \over {90}} f^{(4)} ( \xi )$$
• [2]: $$\tilde{E}_{n}^{2} (f) = - {{ h^4 } \over {180}} [ f^{(3)} (b) - f^{(3)} (a) ]$$

説明

$I_{n}^{2} (f)$を展開すると次のようになる。 $$ \begin{align*} I_{n}^{2} (f) =& {{h} \over {3}} [ f(x_{0}) + 4 f ( x_{1} ) + 2 f( x_{2} ) + 4 f ( x_{3} ) + 2 f ( x_{4} ) + \cdots \\ & + 2 f (x_{n-2} ) + 4 f ( x_{n-1} ) + f(x_{n} ) ] \end{align*} $$ 確定積分$\displaystyle I (f) = \int_{a}^{b} f(x) dx$の数値積分を求めるために、線形補間を使った台形則とは異なり、二次補間をした方法だ。

特に注目すべき点は、エラーが$\displaystyle E_{1}^{2} (f) = - {{h^5} \over {90}} f^{(4)} ( \xi )$として計算されることだが、これは積分のために$2$次元で補間したにも関わらず$f$が$3$次未満の多項式なら、エラーが正確に$0$になるという意味だ。

証明 ¹

[1]

戦略: 二次関数が与えられた関数の二次補間なので、多項式補間の性質を使える。

---

便宜上$\displaystyle c:= \left( {{a+b} \over {2}} \right)$と置こう。 $$ I_{1}^{2} (f) := \left( {{ b - a } \over { 6 }} \right) \left[ f(a) + 4 f \left( c \right) + f(b) \right] $$ これは区間$[a,b]$で$f$を二次補間してその関数の積分値として$I(f)$を近似したものだ。

ニュートンの差分公式: $$p_{n} (x) =\sum_{i=0}^{n} f [ x_{0} , \cdots , x_{i} ] \prod_{j=0}^{i-1} (x - x_{j} )$$

3点$\displaystyle a , c , b$とすべての$x \in [a,b]$に対して$f(x) = p_{2+1} (x)$とする。 $$ \begin{align*} \displaystyle E_{1}^{2} (f) :=& I(f) - I_{1}^{2} (f) \\ =& \int_{a}^{b} \left[ f(x) - p_{2} (x) \right] dx \\ =& \int_{a}^{b} \left[ p_{2+1} (x) - p_{2} (x) \right] dx \\ =& \int_{a}^{b} (x-a)(x-c)(x-b) f [a,c,b,x] dx \end{align*} $$ ここで$w$を次のように定義する。 $$ w(x) := \int_{a}^{b} (t-a)(t-c)(t-b) dt $$ $(t-a)(t-c)(t-b)$は$t = c$を中心とし奇関数であるため、$w(b)=0$であり、上端と下端が等しい場合、積分値が$0$になるので、$w(a) = 0$だ。すると部分積分法と差分に関する微分式に従って $$ \begin{align*} \displaystyle E_{1}^{2} (f) =& \int_{a}^{b} w’(x) f [a,c,b,c] dx \\ =& \left[ w(x) f [ a,c,b, x] \right]_{a}^{b} - \int_{a}^{b} w(x) {{d} \over {dx}} f [a,c,b,x] dx \\ =& - \int_{a}^{b} w(x) f [a,c,b,x,x] dx \end{align*} $$ 微積分学の基本定理により$w(x) \ge 0$であるので、積分の平均値定理を使える。

積分の平均値定理: 閉区間$[a,b]$で関数$f$が連続で$w(x) \ge 0$が積分可能ならば、$\displaystyle \int_{a}^{b} f(x) w(x) dx = f( \eta ) \int_{a}^{b} w(x) dx$を満たす$\eta$が$[a,b]$に少なくとも一つ存在する。

そして、差分の性質により、 $$ \begin{align*} \displaystyle E_{1}^{2} (f) =& - f [a,c,b,\eta,\eta] \int_{a}^{b} w(x) dx \\ =& - {{ f^{(4)} ( \xi ) } \over {24}} \left[ {{4} \over {15}} h^5 \right] \\ =& - {{h^5} \over {90}} f^{(4)} ( \xi ) \end{align*} $$ 何らかの$ \eta , \xi \in [a,b]$が存在することになる。

■

証明[2]

戦略: リーマン和を導き出せば、後は微積分学の基本定理によって自然と導かれる。ただし、シンプソンのルールは$\displaystyle \sum_{k=1}^{n/2}$のように加算を行うので、リーマン和を導く過程で$\displaystyle {{2} \over {2}}$を掛けるトリックが使われる。

---

定理[1]により、実際の$I(f)$と$I_{n}^{2} (f)$の誤差はいくつかの$\xi_{k} \in [x_{2(k-1)}, x_{2k} ]$に対して次のように計算される。 $$ \begin{align*} \displaystyle E_{n}^{2} (f) =& I (f) - I_{n}^{2} (f) \\ =& \sum_{k=1}^{n/2} \left( - {{ h^5 } \over { 90 }} f^{(4)} ( \xi_{k} ) \right) \end{align*} $$ これについて、 $$ \begin{align*} \lim_{n \to \infty} {{ E_{n}^{2} (f) } \over { h^4 }} =& \lim_{n \to \infty} {{1} \over {h^4}} \sum_{k=1}^{n/2} \left( - {{ h^5 } \over { 90 }} f^{(4)} ( \xi_{k} ) \right) \\ =& \lim_{n \to \infty} {{1} \over {h^4}} {{2} \over {2}} \sum_{k=1}^{n/2} \left( - {{ h^5 } \over { 90 }} f^{(4)} ( \xi_{k} ) \right) \\ =& - {{ 1 } \over { 180 }} \lim_{ n \to \infty} \sum_{k=1}^{n/2} 2 h f^{(4)} ( \xi_{k} ) \\ =& - {{ 1 } \over { 180 }} \lim_{ n \to \infty} \sum_{k=1}^{n/2} {{b-a} \over {n/2}} f^{(4)} ( \xi_{k} ) \\ =& - {{ 1 } \over { 180 }} \int_{a}^{b} f^{(4)} (x) dx \\ =& - {{ 1 } \over { 180 }} [ f^{(3)} (b) - f^{(3)} (a) ] \end{align*} $$ したがって、 $$ \lim_{n \to \infty} {{\tilde{E}_{n} (f) } \over { E_{n} (f) }} = 1 $$

$$ E_{n}^{2} (f) \approx \tilde{E}_{n}^{2} (f) = - {{ h^4 } \over { 180 }} [ f^{(3)}(b) - f^{(3)}(a) ] $$

■