1 常微分方程式

1.1 常微分方程式のイメージ

微分方程式は，物理や工学の分野で問題を解く強力なツールばかりか，生物や経済などで も広く応用されている．自然科学や工学の問題を定量的に考察する場合，微分方程式ほど 強力な道具はない．この微分方程式を使うためには，方程式をる作ことと解くことが必要 である．ここでは，微分方程式を解くこと，特に数値計算により非常に精度の良い近似値 を求める方法を学習する．微分方程式には解析解が無いのが普通であるが，理工学上の諸 問題では精度良く解の近似値を求めたい状況にしばしば遭遇する．このような時，数値計 算の出番となる．数学に無い面白さがありますので，楽しんでください．

すでに学習したように，独立変数が二つ以上の多変数の関数の微分(偏微分)を含む微分方 程式を偏微分方程式(partial differential equation)という．それに対して，一変数の 関数の微分を含む方程式を常微分方程式(ordinary differential equation)という．ここ では，常微分方程式，特に1 階の場合の解の近似値を求める方法を学習する．学習する方 程式の形は

$$\if 11 \frac{\mathrm{d}y}{\mathrm{d}x} \else \frac{\mathrm{d}^{1} y}{\mathrm{d}x^{1}}\fi =f(x,y)$$ (1)

である．1階だといってバカにはできない．後で述べることになるが，これが数値計算 できると，どんな高階の常微分方程式も同じ方法で計算ができるのである．数学だと1 解 が解ければ高階の微分方程式が解けるわけではないが，数値計算では可能なのである．

ここでの主題は，この微分方程式を満たす$y(x)$を求めることになる．計算を進める前に， この方程式が何を表すか考えることにする．式(1) の左辺は，解 $y(x)$ の導関数となっている．即ち，解の曲線の接線を表す．導関数の値が座標$(x,y)$ の関数になっているので，座標が決まれば，その場の曲線の傾きが決きまることになる．

それでは，この常微分微分方程式のイメージをつかむことにする．それには，実際 の微分方程式を考えるのが良いであろう．例えば，

$$\frac{\mathrm{d}y}{\mathrm{d}x}=\sin x \cos x -y\cos x$$ (2)

のような常微分方程式を考えることにする．いかにも難しげな微分方程式であるが，これ には解析解がある．解析解はとりあえずおいておくことにして，この式の右辺を考える． 先ほど述べたように，これは接線の傾きを表す．場所ごとに接線の傾きが決まっているの で，それを $\mathrm{xy}$平面に図示することができる．式(2)の右辺の 値である各座標の傾きを線の傾きで表すと，図1のようになる． この傾きを方向場と言う．方向場から，大体の解の様子がわかる．

この微分方程式の解析解は，

$$y=\sin x -1+c_1e^{-\sin x}$$ (3)

である．1階の微分方程式なので，1個の未知数を含む．この未知数の値が異なる5本の曲 線と，先ほどの方向場を重ねて書きすると，図2のようになる．微分方 程式の解である曲線$y(x)$が方向場に沿うことが理解できるであろう．元の微分方程式が 傾きを表すので，あたりまえのことである．

式(2)の微分方程式から，関数$y(x)$の値を得るにはもう一つ条件が必 要である．通常この条件は， $y(x_0)=y_0$のように与えられる．これを初期値といい，初 期値が与えられるものを初期値問題という．一方，2 点以上のxで定めるyの値が決まって いるような問題を境界値問題という²．ここでは，もっぱら初期値問題 を解くことにする．

図: 微分方程式 $\frac{\mathrm{d}y}{\mathrm{d}x}=\sin x \cos x -y\cos x$の方向場

図 2: 方向場と解曲線

1.2 数値計算のイメージ

初期値問題を計算するルーチンの基礎的な考え方はどれも似通っており，次に述べるとお りである．まず，(1)式の微分方程式を極限の $\mathrm{d}$の代わりに有 限な$\Delta$ に置き換える．$\Delta$が小さければ，元の微分方程式の良い近似になる はずである．すると，式(1)の微分方程式は，

$$\Delta y=f(x,y) \Delta x$$ (4)

のように近似できる．これを用いて，$x_i$から$\Delta x$離れた$y$の値$y_{i+1}$を計 算する．

$$y_{i+1} = y(x_i+\Delta x)$$

$$= y_i+\Delta y$$

$$= y_i+f(x_i,y_i)\Delta x$$ (5)

この式と初期値${x_0,y_0}$を用いると，次々に $(x_1,y_1), (x_2,y_2), (x_3,y_3), \dots$が計算できる．

式(5)は，

• 次の$y_{i+1}$値は，もとの$y_i$に，そこでの傾き $f(x_i,y_i)$に$x$の増分$\Delta x$を乗じたものを加えることにより求められる．

と言っているのである．イメージにすると，図3のようになる． この図からも分かるようにこの方法をそのまま適用した場合(オイラー法)，精度がよくな い．出発点のみの導関数を用いているため，終点付近では傾きが異なるからである．刻み 巾$\Delta x$を小さくすることにより解決できるが，その分，計算時間が必要になる．そ のため，$x_i$と$x_{i+1}$の間で，出来るだけ精度よく，この間の導関数の平均を計算する工 夫がいろいろ考えられている．これから，以降その方法を示すことになる．

図 3: 方向場と微分方程式の解 $(x_i, y_i)$と $(x_{i+1}, y_{i+1})$

ホームページ: Yamamoto's laboratory
著者: 山本昌志