オイラー法の誤差について - きゅーぶのブログ

期末試験でオイラー法に関する問題が出たんだけど、試験時間内に解けなかったので代わりにここで解くことにします。

オイラー法っていうのは、常微分方程式の近似解を数値的に求める方法で、ルンゲ・クッタ法のもうちょっと簡単バージョンみたいなやつです。具体的には、次のようなものです。

${\lbrack0,1\rbrack\times \mathbb{R}}$ 上の ${C^1}$ 級関数 ${f=f(x,y)}$ に対して、
微分方程式の初期値問題
\begin{align*}
\begin{cases}
y' = f(x,y) &(0\leq x\leq 1),\\
y(0) = y_0
\end{cases}
\end{align*}
が ${\lbrack0,1\rbrack}$ 上 ${C^2}$ 級の解 ${y=y(x)}$ を持つという状況を考えます。
いま、解 ${y(x)}$ の具体的な形はわかってないけど、分点 ${x_0,x_1,\dots , x_n}$ における近似値 ${Y_0,Y_1,\dots ,Y_n}$ を数値的に求めたいなあと思うことにします。
オイラー法と呼ばれる方法では、近似列 ${\{Y_i\}_{i=0}^{n}}$ を次で定めます。

まず、 ${\lbrack 0,1\rbrack}$ を ${n}$ 等分して、 ${x_i = i/n}$ とします。
また、 ${h=1/n}$ として、
\begin{align*}
\begin{cases}
Y_{i+1}= Y_i + hf(x_i , Y_i), \\
Y_0 = y_0
\end{cases}
\end{align*}
で近似列 ${\{Y_i\}_{i=0}^{n}}$ を定めます。

気持ちとしては、 ${y(0)}$ の値と ${y(x)}$ の傾きだけを頼りに、 ${x_{i}}$ での ${y(x)}$ の値を見積もっていこうという算段です。具体的な関数でグラフを書いたらわかったような気になると思います。
あ、Wikipediaに画像があったので貼っておきます。
Euler method - Wikipedia, the free encyclopedia
f:id:hinabit:20150213011059p:plain

それで、このオイラー法に関して期末試験で次のような問題が出題されました。

定数 ${K>0}$ が存在して、任意の ${x\in\lbrack0,1\rbrack}$ および ${y,z\in \mathbb{R}}$ に対し、
\begin{align*}
\lvert f(x,y)-f(x,z)\rvert \leq K\lvert y-z\rvert
\end{align*}
が成り立つとします。（以下、この条件をリプシッツ条件と呼ぶことにします。）
このとき、ある定数 ${B>0}$ （ ${h}$ に無関係）が存在し、
任意の ${h\geq 0}$ に対し、
\begin{align*}
\max_{i=1,\dots,n}\lvert y(x_i)-Y_i \rvert \leq Bh
\end{align*}
が成り立つことを示しなさい。

実際の試験ではもうちょっと誘導があったのですが、「そんなのがなくても解けるぜ！」という硬派な方のため、ここでは省略しました。
以下で証明を与えます。

証明

まず、 ${y(x)}$ に関するテイラーの定理より、各 ${x_i}$ について、
\begin{align*}
y(x_{i+1})=y(x_{i}) + hf(x_i,y(x_i)) + \frac{h^2}{2} y''(\xi)
\end{align*}
を満たす ${\xi \in ( x_{i} , x_{i+1} )}$ が存在することがわかります。

この ${\xi}$ を用いると、
\begin{align*}
y(x_{i+1})-Y_{i+1}=y(x_i) - Y_i + h \{f(x_i,y(x_i)) - f(x_i,Y_i) \} + \frac{h^2}{2}y''(\xi_i)
\end{align*}
となります。つまり、 ${y(x_{i+1})=y(x_{i}) + hf(x_i,y(x_i)) + \frac{h^2}{2} y''(\xi)}$ から ${Y_{i+1}= Y_i + hf(x_i , Y_i)}$ を引いて少し整理しただけです。

グラフで書くとこんな感じです。
f:id:hinabit:20150213014104j:plain

ここで、三角不等式を用いれば、
\begin{align*}
\lvert y(x_{i+1})-Y_{i+1}\rvert\leq\lvert y(x_i) - Y\rvert + h \lvert f(x_i,y(x_i)) - f(x_i,Y_i) \rvert + \frac{h^2}{2}\lvert y''(\xi_i)\rvert
\end{align*}
が成り立ちます。

いま、 ${\lvert y''(x) \rvert \leq M}$ ${(0\leq x \leq 1)}$ となる定数 ${M > 0}$ を取ります。
すると、リプシッツ条件（問題文の2行目に書いてあるやつ）と、ついさっき取った ${M}$ を用いて、
\begin{align*}
\lvert y(x_{i+1})-Y_{i+1}\rvert\leq(1+hK)\lvert y(x_i) - Y\rvert + \frac{h^2 M}{2}
\end{align*}
と評価できます。

これにより、 ${x_{i+1}}$ での誤差評価を ${x_i}$ での誤差評価に落としこむことが出来ました。
これを繰り返し用いることで、
\begin{align*}
\lvert y(x_{i+1})-Y_{i+1}\rvert&\leq(1+hK)\lvert y(x_i) - Y_i\rvert + \frac{h^2 M}{2} \\
&\leq(1+hK)\{(1+hL)\lvert y(x_{i-1}) - Y_{i-1}\rvert + \frac{h^2 M}{2}\}+\frac{h^2 M}{2} \\
&=(1+hK)^{2}\lvert y(x_{i-1}) - Y_{i-1}\rvert +(1+(1+hK)) \frac{h^2 M}{2}\} \\
&\leq(1+hK)^{3}\lvert y(x_{i-2}) - Y_{i-2}\rvert +(1+(1+hK)+(1+hK)^{2}) \frac{h^2 M}{2}\} \\
&\vdots \\
&\leq (1+hK)^{i+1}\lvert y(x_{0}) - Y_{0}\rvert +(1+(1+hK)+\dots (1+hK)^{i}) \frac{h^2 M}{2}\} \\
&=(1+(1+hK)+\dots (1+hK)^{i}) \frac{h^2 M}{2}\} \\
&=\frac{1}{hK}((1+hK)^{i+1} -1) \frac{h^2 M}{2}\} \\
&=\frac{hM}{2K}((1+hK)^{i+1} -1)
\end{align*}
と誤差評価ができます。長いですが、途中に ${\lvert y(x_{0}) - Y_{0}\rvert =0}$ であることや、等比数列の和の公式などを使って計算しています。

ところで、
\begin{align*}
(1+hK) \lt e^{hK}
\end{align*}
ですので、上の式はさらに次のように評価できます。
\begin{align*}
\frac{hM}{2K}((1+hK)^{i+1} -1) &\leq \frac{hM}{2K}(e^{hK(i+1)}-1) \\
&= \frac{hM}{2K}(e^{\frac{i+1}{n}K}-1) \\
&\leq \frac{hM}{2K}(e^{K}-1).
\end{align*}

まとめると、
\begin{align*}
\max_{i=1,\dots,n}\lvert y(x_i)-Y_i \rvert \leq Bh \leq \frac{hM}{2K}(e^{K}-1)
\end{align*}
がわかり、 ${\frac{M}{2K}(e^{K}-1)}$ は ${h}$ に無関係な定数なので、示したかったことが言えました。□