[미분방정식] Euler's Method

한경민·2023년 5월 9일

미분방정식

목록 보기
1/3
post-thumbnail

미분방정식이 무엇인지 이해하고, 파이썬으로 해를 근사시키는 방법을 알아봅시다.

대한민국에서 의무교육을 받은 분이라면 방정식이라는 용어가 생소하진 않을 것입니다. 방정식은 미지수의 값에 따라 참이 되기도, 거짓이 되기도 하는 등식을 의미하죠. 우린 일차방정식을 풀 줄 알고, 이차방정식은 근의 공식(이에이분에)을 달달달 외운 적도 있습니다. 이러한 방정식은 수학에서 굉장히 중요한 의미를 갖습니다. 세상에서 일어나는 여러 현상들을 방정식을 통해 수학적 모델링을 할 수 있기 때문입니다.

방정식에 대한 리마인드는 이 정도면 된 거 같네요. 본격적으로 이번 포스트의 주제로 들어가보도록 합시다.

미분방정식?

식에 미분의 꼴이 나타나는 방정식을 미분방정식이라고 합니다.
예를 들어,

y=7x+2yy'=7x+2y

꼴의 방정식은 미분방정식입니다.

특히 미분방정식에선 차수(Order)라는 개념이 있는데, 이는 다항식에서의 차수하고는 의미가 조금 다릅니다.
미분방정식은 미분의 꼴이 나타나는 방정식이라고 했죠?? 바로 변수 중 가장 많이 미분이 된 그 횟수를 미분방정식의 차수라고 합니다. 위 예시의 경우는 차수가 1이라고 할 수 있겠네요.

미분방정식은 여러 사회현상이나 자연현상을 수학적으로 모델링할 때 굉장히 유용하게 사용됩니다. 세상의 모든 것은 조금씩 변화하기 때문입니다. 미분의 의미는 어느 한 시점에서의 순간변화율을 구하는 것이기 때문에 변화를 수학적으로 서술하기에 미분만한 것이 없는 것이죠!

Euler's Method?

그래서 난데없이 등장한 Euler's Method는 무엇일까요?
모든 미분방정식의 해를 간단히 구할 수 있다면 얼마나 좋을까요...? 안타깝게도 그렇지 않습니다. 만약 그랬다면 나비에-스토크스 방정식이 밀레니엄 문제로 채택되는 일은 없었을 겁니다.

이런 경우 근사(Approximation)란 방법을 사용할 수 있습니다!(이러한 아이디어는 수치해석에서 굉장히 유용합니다.) 수학적 방법을 이용해서 실제값에 가깝게 값을 맞추어 가는 과정입니다. 컴퓨터를 조금 배우신 분이면 컴퓨터가 소수점을 나타내는 방법 중 부동소수점 구조를 들어보셨을텐데, 이 또한 근사입니다.

한편, 굉장히 구하기 힘든 미분방정식이 있다고 해봅시다. 예를 들어,

y=f(x, y), y(x0)=y0y'=f(x,\space y),\space y(x_0)=y_0

로 미분방정식과 초기 조건이 주어졌다고 합시다. 이 경우, Euler's Method란 다음과 같습니다.

yn=yn1+hf(xn1, yn1), xn=x0+nhy_n=y_{n-1}+hf(x_{n-1},\space y_{n-1}),\space x_n=x_0+nh

이 때 hhxn1x_{n-1}xnx_n의 간격입니다. xnx_n에서 접선의 기울기를 xn1x_{n-1}으로부터 hh만큼 떨어진 거리에 있는 점의 기울기로 설정하여 그래프를 그려나간다면, 그 그래프가 실제 방정식의 해의 그래프와 비슷할 것이라는 아이디어입니다.

식을 그림으로 표현하면 다음과 같습니다.

코드

그림으로 살펴보았으니 직접 코드를 작성해서 살펴보도록 합시다! 코드는 Python으로 작성하였습니다.

import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt

def model(y, t):
    k = -5
    dydt = t + y
    return dydt
    
approximation = []
y_n = y0
for i, val in enumerate(t):
    approximation.append(y_n)
    print(f"{i + 1}번째 시도")
    print(f"x_n = {val:.1f}")
    print(f"y_n = {y_n:.3f}") # 근사
    print(f"y(x_n) = {y[i][0]:.3f}") # 해
    print(f"Error : {y[i][0] - y_n:.3f}")
    print()
    
    y_n += h * model(y_n, val)
    
plt.plot(t, y, label='solution curve')
plt.plot(t, approximation, label='approximation curve')
plt.legend()
plt.show()

scipy를 이용하여 미분방정식을 직접 풀고, 이를 Euler's Method를 푼 결과와 비교해서 그래프를 그립니다. 결과는 아래와 같습니다.

profile
📁예비 데이터 엔지니어

0개의 댓글