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

Euler法の場合は

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

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

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

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

実際には、$x$, $y$ は配列 x[], y[] ではなく、 スカラー変数 x, y に記憶することが多いであろうから、 次のどちらかになるだろうか？⁵

  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法の場合は

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

$$\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)$$ (17)

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

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

で $k_{1,x}$, $k_{1,y}$, $k_{2,x}$, $k_{2,y}$, $k_{3,x}$, $k_{3,y}$, $k_{4,x}$, $k_{4,y}$ を計算して

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

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