[수치해석] Numerical solution of ODE (5) Predictor-Corrector Method

글에서 수치해석 방법을 줄여서 씁니다.

EE : Explicit Euler

IE : Implicit Euler

TR : Trapezoidal method

 

Explicit method가 implicit method 대신 쓰이는 이유

 

지금까지 Explicit Euler, Implicit Euler, Trapezoidal method를 살펴보았는데

trapezoidal method는 accuracy에서나, stability에서나 좋은 점이 많습니다. EE, IE 방법보다 정확도도 1차 높고 stability 영역도 exact solution과 동일합니다.

 

하지만 그럼에도 불구하고 Trapezoidal method를 쓰지 않고 Predictor-Corrector method(PC), 또는 Runge-Kutta method(RK)를 쓰는 이유가 있습니다. PC와 RK는 Explicit method로 분류됩니다. 즉, stability를 따질 때 conditionally stable(C-stable)합니다. 이 말을 또 다르게 표현하면 $\Delta t$ 크기가 기준에 만족하지 않으면 solution이 발산할 수 있다는 것입니다.

 

Predictor-Corrector 방법을 소개하기 전에 왜 conditionally stable한 explicit method를 쓰는지 이야기하려고 합니다.

 

예를 들어

$y^{\prime}=\lambda y^{2}$

다음과 같은 nonlinear system을 수치해석적인 방법으로 풀어보도록 하겠습니다.

Trapezoidal method로 문제를 풀면

 

$\frac{y^{n+1}-y^{n}}{\Delta t}=\frac{1}{2} \lambda\left( (y^{n})^{2}+(y^{n+1})^{2}\right)$

 

이 식을 좌변은 (n+1)항에 대해서, 우변은 (n)항에 대해서 정리하면 다음과 같다.

 

$y^{n+1}-\frac{1}{2} \lambda \Delta t (y^{n+1})^{2}=y^{n} + \frac{1}{2}\lambda \Delta t (y^{n})^{2}$

 

여기서 문제는 우리가 알고 있는 값이 $y^{n}$이고 $y^{n+1}$을 구해야 한다는 것이다.

위 식 대로라면 이차방정식을 풀어야 $y^{n+1}$을 구할 수 있다.

 

이차방정식 정도는 풀 수 있다고 생각할 수 있지만 더 복잡한 nonlinear DE를 풀어야 할 경우에는 어려워집니다.

하지만 Explicit method를 사용하면 훨씬 쉽게 풀 수 있습니다.

예를 들어 위 문제를 Explicit euler 방법을 쓰면

 

$\frac{y^{n+1}-y^{n}}{\Delta t}=\lambda (y^{n})^{2}$
$y^{n+1}=\lambda \Delta t (y^{n})^{2}+y^{n}$

 

이 경우에는 우변에 $y^{n}$를 대입함으로써 쉽게 $y^{n+1}$을 구할 수 있습니다.

이런 explicit method의 장점으로 더 정확한 explicit method를 사용하는 것이 nonlinear일 때 좋습니다.

---

Predictor-Corrector method (줄여서 PC method로 표현)

 

이름이 Predictor-Corrector인데 Explicit Euler 방법으로 해를 Predict하고 Trapezoidal 방법으로 Correct(수정)하기 때문입니다.

 

Predictor : $y^{n+1*}=y^{n}+\Delta t f\left(y^{n}, t^{n}\right)\leftarrow$ 처음에 EE로 임시 $y^{n+1*}$값을 가집니다.

Corrector : $y^{n+1}=y^{n}+\frac{1}{2} \Delta t\left[f\left(y^{n}, t^{n}\right)+f\left(y^{n+1*}, t^{n+1}\right)\right]$

Predictor에서 구한 임시값으로 TR을 적용합니다.

 

이 방법을 model problem으로 accuracy, stability, phase error를 살펴보겠습니다.

 

 

$\text{model problem : }y^{\prime}=\lambda y$
$\text{Predictor : }y^{n+1*}=y^{n}+\lambda \Delta t y^{n}$

$\text{Corrector : }y^{n+1}=y^{n}+\frac{\lambda \Delta t}{2}\left[y^{n}+y^{n+1*}\right]$

 

Predictor에서 구한 식을 Corrector에 대입하면 다음과 같습니다.

 

$\begin{aligned} y^{n+1} &=y^{n}+\frac{\lambda \Delta t}{2}\left[y^{n}+y^{n}+\lambda \Delta t y^{n}\right] \\ &=\left(y^{n}+\lambda \Delta t+\frac{(\lambda \Delta t)^{2}}{2}\right) y^{n} \end{aligned}$

 

1. Accuracy

 

이 식을 amplification factor 관점에서 보면,

$y^{n}=\sigma_{pc}^{n} y^{0} \quad\rightarrow \sigma_{pc}=1+\lambda \Delta t+\frac{(\lambda \Delta t)^{2}}{2}$

 

exact solution과 비교하겠습니다.

$\text{exact solution : }y^{n}=e^{\lambda \Delta t} y^{0} \quad\rightarrow e^{\lambda \Delta t}=1+\lambda \Delta t+\frac{(\lambda \Delta t)^{2}}{2}+\frac{\lambda \Delta t}{6}+\cdot$

 

exact solution은 $\Delta t^{2}$항까지 동일하므로 accuracy는 locally 3rd order accuracy라고 결론내릴 수 있습니다.

globally 2nd order accuracy 

 

2. Stability

 

$\lambda=\lambda_{R+i} \lambda_{I}$
$\left|\sigma_{pc}\right|=\left|1+\left(\lambda_{R} \Delta t+i \lambda_{I} \Delta t\right)+\frac{\left(\lambda_{R} \Delta t+i \lambda_{I} \Delta t\right)}{2}\right| \leq 1$

 

다음을 만족하는 영역을 찾으면,

EE방법과 비교했을 때 빨간 타원과 같은 stability 영역이 나옵니다.

이 방법은 EE와 마찬가지로, imaginary 성분만 있는 시스템에 대해서 unstable한 결과가 나옵니다.

 

대신 EE보다는 더 넓은 영역에서 stable하고 EE보다 정확도가 높다는 장점을 가집니다.

 

매트랩으로 나타낸 stability 영역입니다.

 

3. Phase error

 

$\lambda=i \omega$일 때 unstable하므로 따로 phase error를 다루지 않겠습니다.

 

---

다음에는 가장 많이 쓰이는 수치해석 방법인 Runge-Kutta method에 대해 다루겠습니다.