scipy 로 상미방의 수치 적분 수행하기

1. 개요

 상미분방정식(이하 상미방)의 수치 적분을 수행하기 위해서 scipy.integrate.odeint 함수의 사용법에 대해서 알아보도록 하겠다. scipy.integrate 클래스는 적분을 수행하는 다양한 함수들이 모여 있는데 그 중 연립 상미방의 수치적분을 구하는 ode, odeint, complex_ode 함수들이 마련되어 있다.

 

• odeint             : 상미방의 수치적분

• ode                : 상미방의 수치적분 (generic interface class)

• complex_ode  : 복소 상미방의 수치적분

 

이중 odeint 는 내부적으로 FORTRAN 라이브러리인 odepack 의 'lsoda'를 사용하여 상미방을 풀어내는 함수이다. 이 함수는 (연립) 상미분 방정식

• y' = f ( y, t )

의 수치해를 구해주며 stiff ode 와 non-stiff ode 모두에 적용된다.

2. 기본 문법

 일단 레퍼런스의 함수 문법은 다음과 같다.

y, infodic = scipy.integrate.odeint (    func,# callable(y, t, ...) : 시간 t 에서 y(t)의 미분값을 구할 수 있는 함수    y0,  # array_like : y 의 초기 조건 (벡터도 된다.)    t,   # array_like : y(t)값을 구할 t 점들. 초기값의 시간이 이 배열의 첫 요소여야 함.    args=()  # 튜플 : func() 에 넘길 추가적인 인수들    Dfun=None, col_deriv=0, full_output=0, ml=None, mu=None, rtol=None, atol=None,    tcrit=None, h0=0.0, hmax=0.0, hmin=0.0, ixpr=0, mxstep=0, mxhnil=0, mxordn=12,    mxords=5, printmessg=0 )

 

여기서 필수적인 인수는 func, y0, t 세 개이며 각각 ODE함수, 초기값, 시간(벡터)이다. 나머지 인수들은 모두 선택적으로 지정해 줄 수 있으며 지정하지 않으면 설정값이 사용된다. 그리고 y0와 t는 리스트, 튜플, ndarray (이것들을 묶어서 array_like 라고 한다.) 등이 될 수 있다.

 이 함수를 사용하기 전에 먼저 상미방을 기술하는 함수를 먼저 정의해야 한다. 이 함수는 array_like 를 반환해야 하며 ndarray y와 시간 t 를 받아서 그 미분을 구하는 함수이다.

 odeint() 함수의 출력 y는 배열(ndarray) 이고 shape 이 ( len(t), len(y0) ) 이며 2차 배열이다. 입력으로 넘어온 t 배열의 각 점에서의 y값들을 갖는다. 두 번째 출력인 infodict 는 full_output== True 로 지정되었을 경우 추가적인 정보가 넘어오게 된다.

3. 간단한 1차 ode 예제

 다음과 같이 가장 간단한 1차 미방을 시간 구간 [0, 5] 에서 수치해를 구하는 예제를 해 보도록 하겠다.

• y' =  -2y, y(0)=1

위 예제에서는 pylab 모듈에 numpy가 포함되어 있으므로 따로 numpy는 import하지 않았다.

4. 2차 ode : 조화진동자

  조화 진동자 방정식은 2차 상미방으로서 x''(t) = -x(t) 이다. 2차 방정식은 두 개의 1차 방정식으로 분해할 수 있다. 즉, y(t) = [ x(t) ; x'(t) ] =: [y1(t); y2(t)] 라고 정의하면

              y'(t) = [x'(t); -x(t)] = [y2(t); -y1(t)]

와 같이 기술할 수 있다. 초기값이 [0; 1] 일 경우 시간 [0, 5]에 대해서 이 상미방을 풀면 다음과 같다. 이 예제에서는 ode1()함수의 반환값도 튜플이고 odeint()함수의 초기값도 튜플로 주었다.

파이썬 강좌 전체 목록 (TOP) >>>

c{ns},n{ns001}