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
著者: 山本昌志