3.2 Euler法, Runge-Kutta法はベクトル値関数でもOK

Euler法の場合は

(14) $\displaystyle \vec{x}_{i+1}=\vec{x}_i+\Delta t \vec{f}(\vec{x}_i,t_i).$

もしベクトルを扱えるプログラミング言語ならば、 この形のままでプログラムが書ける。

そうでない場合もこの式を成分表示した

\begin{subequations}% 2021-05-6 01:14の式群
\begin{align}&x_{i+1}=x_{j}+\Delt...
..._i,t_i) &y_{i+1}=y_{j}+\Delta t f_y(x_i,y_i,t_i) \end{align}\end{subequations}

を使えば良い。(ここで $ f_x$, $ f_y$$ \vec f$ の第1成分($ x$成分)、 第2成分($ y$成分) を表す。偏微分ではない。)


実際には、$ x$, $ y$ は配列 x[], y[] ではなく、 スカラー変数 x, y に記憶することが多いであろうから、 次のどちらかになるだろうか?5
  dx = dt * fx(x,y,t);
  dy = dt * fy(x,y,t);
  x += dx;
  y += dy;
  newx = x + dt * fx(x,y,t);
  newy = y + dt * fy(x,y,t);
  x = newx;
  y = newy;
(要点は x, y の更新は、 $ f$ の値の計算が終わってから、ということである。)


Runge-Kutta法の場合は

\begin{subequations}\begin{align}&\vec{k}_1=\Delta t\vec{f}(\vec{x}_i,t_i), &\...
...4=\Delta t\vec{f}(\vec{x}_i+\vec{k}_3,t_i+\Delta t)\end{align}\end{subequations}

$ \vec{k}_1$, $ \vec{k}_2$, $ \vec{k}_3$, $ \vec{k}_4$ を計算して 6

(17) $\displaystyle \vec{x}_{i+1}=\vec{x}_i+\frac{1}{6} \left( \vec{k}_1+2\vec{k}_2+2\vec{k}_3+\vec{k}_4 \right)$

とすれば良い。 やはりベクトルを扱えるプログラミング言語ならば、 この形のままでプログラムが書ける。

ベクトルが扱えないプログラミング言語ならば、 (16a)-(16d) を成分表示した

\begin{subequations}\begin{align}&k_{1,x}=\Delta t\vec{f}_{x}(x_i,y_i,t_i), &k...
...}_{y}(x_i+k_{3,y}/2,y_i+k_{3,y}/2,t_i+\Delta t/2), \end{align}\end{subequations}

$ k_{1,x}$, $ k_{1,y}$, $ k_{2,x}$, $ k_{2,y}$, $ k_{3,x}$, $ k_{3,y}$, $ k_{4,x}$, $ k_{4,y}$ を計算して
\begin{subequations}\begin{align}&x_{i+1}=x_i+\frac{1}{6} \left( k_{1,x}+2k_{2,x...
...} \left( k_{1,y}+2k_{2,y}+2k_{3,y}+k_{4,y} \right) \end{align}\end{subequations}

とすれば良い。 -- これはやっていることはベクトルを成分表示しただけなのだが、 私の経験によると間違える人がとても多い。 間違えることを見込んで、対応するための時間を授業で用意したこともあるが、 今は「そういうプログラミング言語を使うのは拙い、 プログラミング言語を換えるべき」と考えている。


この例では $ \vec{f}$$ 2$次元であるから、 (18a)-(18h) という8つの式で済んでいるが、 ベクトル $ \vec{x}$ の次元 $ n$ が高くなると面倒になる。 その場合は、 $ \vec{x}$, $ \vec{k}_1$, $ \vec{k}_2$, $ \vec{k}_3$, $ \vec{k}_4$ などを配列で取り扱うプログラムが多い (この文書ではそれは説明しない)。

桂田 祐史