next up previous
Next: 8 Stiff problem (硬い問題) Up: 7 数値的安定性 Previous: 7 数値的安定性

7.1 線形安定性解析

$ \lambda\in\C$ $ {\rm Re}\lambda<0$ なる定数として

$\displaystyle \left\{
\begin{array}{lcl}
x'(t) &=& \lambda x(t) \quad\hbox{($t\in I\DefEq(a,+\infty)$)} \\
x(a) &=& x_0
\end{array}\right.
$

について適用してみる。この解は $ x(t)=x_0 e^{\lambda t}$ であり、 $ \left\vert x(t)\right\vert=\left\vert x_0\right\vert e^{\mathrm{Re}\;\lambda t}$ であるから、 $ \dsp\lim_{t\to\infty}x(t)=0$ である。 ところが、数値解法ではしばしば
$ \exists h_1>0$ , $ \exists h_0\in (0,h_1)$     s.t.     $ \left\{
\begin{array}{l}
\forall h>h_1\quad\dsp\lim_{j\to\infty}\vert x_j\ve...
...\forall h<h_0 \quad \dsp\lim_{j\to\infty}\vert x_j\vert=0
\end{array} \right.
$
ということが起こる。

1)
中点則の場合 (leapfrog method とも言うらしい). これはいわゆる中心差分近似

$\displaystyle x'(t)=\frac{x(t+h)-x(t-h)}{2h}+O(h^2)
$

を用いて作った公式

のことを指す。 これは任意の $ h>0$ に対して $ \lim_{j\to\infty} \vert x_j\vert=+\infty$ を満たす。

中心差分近似に基づく中点公式 $ x_{n+1}=x_{n-1}+h f(t_n,x_n)$ $ x'=\lambda x$ , $ x(0)=1$ を解く( $ \lambda = -4.35$ )。左側は $ x_1=e^{\lambda h}$ としたもの、右側は $ x_1=x_0+h f(t_0,x_0)$ とした もの。

\includegraphics[width=7cm]{figure/test-4.35.eps}      \includegraphics[width=7cm]{figure/test-4.35-by-midpoint2.eps}
2)
R-K 型公式の場合. 公式の段数を $ s$ , 次数を $ m$ とすると、

$\displaystyle x_{j+1}=R(\lambda h)x_j
$

のように書ける。 ここで $ R(z)$

なる $ z$ の有理式である (特に公式が explicit の場合には、 $ R(z)$ $ z$ の多項式になる)。 例えば Euler 法の場合 $ R(z)=1+z$ , 古典的 Runge-Kutta 法の場合、 $ R(z) = 1+z +\Dfrac{z^2}{2} +\Dfrac{z^3}{3!}+\Dfrac{z^4}{4!}$ . この $ R(z)$ に対して

$\displaystyle {\cal R}\DefEq
\left\{
z\in\C; \vert R(z)\vert<1
\right\}
$

絶対安定領域といい、

$\displaystyle \lambda h\in{\cal R} \Then \lim_{j\to+\infty}x_j=0
$

がなりたつ。 絶対安定領域が左半平面 $ \{z\in\C; {\rm Re}\;z<0\}$ を 含むような公式を A-安定 (A-stable) という。 A-安定な公式では、 $ {\rm Re}\;\lambda<0$ なるとき、 任意の $ h>0$ に対して数値解が 0 に収束する。 次に示すのは Runge-Kutta 法の絶対安定領域である。 カラーでお見せできないのが残念ですが。
\includegraphics[width=6cm]{figure/myrunge.eps}      \includegraphics[width=8cm]{figure/Absolute-Stable-Region-for-Runge.eps}

なお、レベル $ 1$ の等高線が負軸と交わる点の座標は(巾根で書けはするの だけれど...)、

$\displaystyle -2.7852935634052816235297591900854630347647578967169\cdots.$

なお、後述の硬い方程式の項を参照せよ。

3)
LM 法の場合. 公式

$\displaystyle (\alpha_0-z\beta_0)x_j+
(\alpha_0-z\beta_1)x_{j+1}+\cdots
(\alpha_0-z\beta_k)x_{j+k}=0
$

に対して特性方程式を

$\displaystyle (\alpha_0-z\beta_0)+
(\alpha_0-z\beta_1)\xi+\cdots
(\alpha_0-z\beta_k)\xi^k=0
$

で定義し、その根を $ \xi_\ell(z)$ ( $ \ell=1,\cdots,k$ ) とする。 この場合

$\displaystyle {\cal R}\DefEq\left\{z\in\C; \vert\xi_\ell(z)\vert<1\right\}
$

とおくと

$\displaystyle \lambda h\in{\cal R}\Then \lim_{j\to\infty}x_j=0.
$


next up previous
Next: 8 Stiff problem (硬い問題) Up: 7 数値的安定性 Previous: 7 数値的安定性
桂田 祐史
2015-05-30