next up previous
: グラフのプロット : 12/3 : 簡単な数値計算 : Newton 法

差分法

微分方程式を数値的に解く (近似解を求める) 方法として差分法がある. $f(x)$ が微分可能のとき,

\begin{displaymath}f(x+h) = f(x)+hf'(x)+o(h)\end{displaymath}

より, $h$ が小さいとき

\begin{displaymath}f'(x) \simeq {(f(x+h)-f(x)) \over h}\end{displaymath}

. 同様に $f(x)$ が二階微分可能のとき

\begin{displaymath}f(x+h) = f(x)+hf'(x)+h^2/2 f''(x)+o(h^2)\end{displaymath}


\begin{displaymath}f(x-h) = f(x)-hf'(x)+h^2/2 f''(x)+o(h^2)\end{displaymath}

より

\begin{displaymath}f(x+h)+f(x-h) = 2f(x)+h^2 f''(x)+o(h^2)\end{displaymath}

よって, $h$ が小さいとき

\begin{displaymath}f''(x) \simeq {{f(x+h)-2f(x)+f(x-h)} \over h^2}\end{displaymath}

微分方程式中の未知関数の微分をこれらの差分で置き換えて得られる差分方程式 (漸化式) を解くことで, もとの方程式の近似解を得る方法が差分法である.

例 13.2   なめらかな水平面上におかれた, ばね定数 $k$ のばねの一方を壁に固定し, 他方に質量 $m$ のおもりをつける. のびちぢみのないときのおもりの位置 を原点とし, ばねののびが $y$ のときおもりの位置が $y$ であるとする. このときおもりが受ける力は $-ky$ だから, おもりの運動方程式は

\begin{displaymath}my'' = -ky\end{displaymath}

となる. $y = y(t)$ で, 微分は時刻 $t$ に関する微分. これを単振動の方程式と呼ぶ. 簡単のために $m=k$ とすると

\begin{displaymath}y'' + y = 0\end{displaymath}

一般解は $A$, $B$ を定数として $y(t) = A\cos(t)+B\sin(t)$ だが, これを 差分化して解いてみよう. $h$ を時刻のきざみ幅とする. すなわち $t = kh$ ($k=0,1,\ldots$) で考える. そして $y_k = y(kh)$ とおき, 上で述べたような置き換えを行うと

\begin{displaymath}{{y_{k+1}-2y_k+y_{k-1}} \over h^2} + y_k = 0\end{displaymath}

すなわち

\begin{displaymath}y_{k+1}=(2-h^2)y_k - y_{k-1}\end{displaymath}

初期条件 $y(0) = a$, $y'(0) = b$$y_0 = a$, $(y_1-y_0)/h = b$ で 置き換えればよい. これで得られる解が真の解に近いかどうかは別の解析が必要.

http://www.math.kobe-u.ac.jp/noro/main/node48.html 参照.



Masayuki Noro 平成14年2月25日