명시적 룽게-쿠타 메소드
개요
상미분 방정식 솔버ODE solver인 룽게-쿠타 메소드를 소개한다. 흔히 사용되는 4차 룽게-쿠타 메소드 RK4에 대한 구체적인 설명은 글이 따로 발행되어있다.
빌드업1
다음과 같이 주어진 상미분 방정식을 생각하자. $y$는 $t$에 대한 함수이며, $^{\prime}$은 $t$에 대한 미분을 의미한다. $$ y^{\prime} = f(t,y),\quad t \ge t_{0},\quad y(t_{0}) = y_{0} $$ 이를 $t_{n}$부터 $t_{n+1} = t_{n} + h$까지 적분하면(헷갈리므로 적분 변수는 $t$ 대신 $\tau$로 두자), $$ \int_{t_{n}}^{t_{n+1}}y^{\prime}d\tau = \int_{t_{n}}^{t_{n+1}} f(\tau,y)d\tau $$ 좌변을 정리하고 $\tau \equiv t_{n} + h\tau$로 치환하면 다음의 식을 얻는다. $$ \big[y(\tau)\big]_{t_{n}}^{t_{n+1}} = \int_{t_{n}}^{t_{n+1}} f(\tau,y)d\tau $$ $$ \Downarrow \\ y(t_{n+1}) = y(t_{n}) + \int_{t_{n}}^{t_{n+1}} f(\tau,y)d\tau = y(t_{n}) + h\int_{0}^{1} f\big(t_{n} + h\tau,y(t_{n} + h\tau)\big)d\tau $$
이제 뒤의 적분을 합으로 근사하여 다음의 식을 얻는다.
$$ y(t_{n+1}) = y(t_{n}) + h\sum_{j=1}^{\nu} b_{j} f\big(t_{n} + c_{j}h,y(t_{n} + c_{j}h)\big),\quad n=0,1,2,\cdots \tag{1} $$
이는 $y_{n} = y(t_{n})$을 알 때, $y_{n+1} = y(t_{n+1})$을 알 수 있는 근사식이다. 물론 우변에 모르는 값인 $y(t_{n} + c_{j}h)$가 포함되어 있으니 이도 근사해야한다. 식이 복잡하므로 $j = 1,2,\dots,\nu$에 대해서 다음과 같이 표기하자.
$$ \xi_{j} = y(t_{n} + c_{j}h) $$
우선 $c_{1} = 0$이라 두자. 그러면 $\xi_{1} = y_{n}$이므로 아는 값이다. 이제 멀티스텝 메소드를 적용하여 $\xi_{j}$들을 다음과 같이 순차적으로 근사한다고 하자.
$$ \begin{align*} \xi_{2} &= y(t_{n} + c_{2}h) \\ &= y_{n} + h a_{2,1}f(t_{n} + c_{1}h, y_{n} + c_{1}h) \\ &= y_{n} + h a_{2,1}f(t_{n}, y_{n}) & (\because c_{1} = 0)\\ &= y_{n} + h a_{2,1}f(t_{n}, \xi_{1}) & (\because \xi_{1}=y_{n}) \\[2em] \xi_{3} &= y(t_{n} + c_{3}h) \\ &= y_{n} + h a_{3,1}f(t_{n} + c_{1}h, y_{n} + c_{1}h) + h a_{3,2}f(t_{n} + c_{2}h, y_{n} + c_{2}h) \\ &= y_{n} + h a_{3,1}f(t_{n}, \xi_{1}) + h a_{3,2}f(t_{n} + c_{2}h, \xi_{2}) \\ &\vdots \\ \xi_{\nu} &= y_{n} + h\sum_{j=1}^{n-1}a_{\nu,j} f(t_{n} + c_{j}h, \xi_{j}) \end{align*} $$
이제 $(1)$의 우변의 값 $y_{n}$, $f(t_{n}, \xi_{1})$, $f(t_{n} + c_{2}h, \xi_{2})$, $\dots$, $f(t_{n} + c_{j-1}h, \xi_{j-1})$들을 모두 구한 셈이다. 명시적 룽게-쿠타explicit Runge-Kutta, ERK 메소드는 $y_{n+1}$의 값을 $y_{n}$, $f(t_{n}, \xi_{1})$, $f(t_{n} + c_{2}h, \xi_{2})$, $\dots$, $f(t_{n} + c_{j-1}h, \xi_{j-1})$의 선형결합으로 근사하여 얻는 방법을 말한다.
정의
명시적 룽게-쿠타 메소드란, 주어진 $y_{n}$에 대해서 $y_{n+1}$을 다음과 같이 근사하는 방법이다.
$$ y_{n+1} = y_{n} + h\sum_{j=1}^{\nu} b_{j} f(t_{n} + c_{j}h, \xi_{j}),\quad n=0,1,2,\cdots \tag{2} $$
여기서 각각의 $\xi_{j}$는,
$$ \begin{align*} \xi_{1} &= y_{n} \\ \xi_{2} &= y_{n} + ha_{2,1}f(t_{n}, \xi_{1}) \\ \xi_{3} &= y_{n} + ha_{3,1}f(t_{n}, \xi_{1}) + ha_{3,2}f(t_{n} + c_{2}h, \xi_{2}) \\ \vdots & \\ \xi_{\nu} &= y_{n} + h\sum_{i=1}^{\nu-1}a_{\nu, i}f(t_{n} + c_{i}h, \xi_{i}) \end{align*} $$
이때 $\nu \times \nu$ 행렬 $A = [a_{j,i}]$를 룽게-쿠타 행렬RK matrix이라 한다. 수식에서 비어있는 성분은 $0$이다.
$$ A = [a_{j,i}] = \begin{bmatrix} \\ a_{2,1} & \\ a_{3,1} & a_{3,2} \\ \vdots & \\ a_{\nu,1} & a_{\nu,2} & \cdots & a_{\nu,\nu-1} & \end{bmatrix} $$
또한 아래의 $\mathbf{b}$와 $\mathbf{c}$를 각각 RK 웨이트RK weights, Rk 노드RK nodes라 한다.
$$ \mathbf{b} = \begin{bmatrix} b_{1} \\ b_{2} \\ \vdots \\ b_{\nu} \end{bmatrix} \quad \text{and} \quad \mathbf{c} = \begin{bmatrix} c_{1} \\ c_{2} \\ \vdots \\ c_{\nu} \end{bmatrix} $$
또한 $(2)$가 $\nu$ 스테이지$\nu$ stages를 가진다고 말하며, 이러한 방법을 $\nu$ 스테이지(혹은 $\nu$ 차) 룽게-쿠타 메소드라 한다.
설명
정의에 따라 자명하게도 RK 행렬은 하삼각행렬이다.
위의 방법에 따라서, $\xi_{j}$는 $\xi_{1}$부터 $\xi_{\nu}$까지 순차적으로 계산할 수 있다. 그러면 결과적으로 $(2)$를 계산할 수 있고, $y_{n+1}$을 얻는다. 이를 반복하여 다음 스텝에서의 $y$값 $y_{n+2}$, $y_{n+3}$, $\dots$ 들을 구할 수 있다.
RK 메소드에서 계수들 $a_{j,i}$, $b_{j}$, $c_{j}$들은 구하는 것이 아니라 선택하여 사용하는 것이다. 계수는 $\begin{array}{c|c} \mathbf{c} & A \\ \hline & \mathbf{b}^{t} \end{array}$와 같은 꼴로 표기하는데 이를 RK 테이블RK tableaux라 한다. 흔히 사용되는 RK4의 RK 테이블은 다음과 같다.
$$ \begin{array}{c|cccc} 0 & \\[0.5em] \frac{1}{2} & \frac{1}{2} \\[0.5em] \frac{1}{2} & 0 & \frac{1}{2} \\[0.5em] 1 & 0 & 0 & 1 \\[0.5em] \hline {\displaystyle {\color{white}\dfrac{0}{0}}} & \frac{1}{6} & \frac{1}{3} & \frac{1}{3} & \frac{1}{6} \end{array} $$
같이보기
Arieh Iserles, A First Course in the Numerical Analysis of Differential Equations (2nd Edition, 2009), p38-41 ↩︎