SHIN, Byeong-Chun's Board

Name	신병춘

Subject	RK2, 4 methods

# Runge-Kutta Method

import numpy as np
import matplotlib.pyplot as plt

def RK2(f, h=1./16, y0=0, a=0, b=1):
    t = np.arange(a,b+h,h); y = np.zeros_like(t)
    y[0] = y0; N = t.size
    for n in range(0,N-1):
        v1 = f( t[n],y[n] ); v2 = f( t[n]+h, y[n]+h*v1 )
        y[n+1] = y[n] + h/2*( v1 + v2 )

    EY = Y(t)
    err = h**(1/2)*np.linalg.norm(EY-y)
    plt.plot(t,EY, t,y,'r--'); plt.pause(.5)
    print(' RK2-err = ', '%12.6e' %err, '  with h = ', h)
    return err

def RK4(f, h=1/16, y0=0, a=0, b=1):
    t = np.arange(a,b+h,h); y = np.zeros_like(t)
    y[0] = y0; N = t.size
    for n in range(0,N-1):
        v1 = f( t[n],y[n] )
        v2 = f( t[n]+h/2, y[n]+h/2*v1 )
        v3 = f( t[n]+h/2, y[n]+h/2*v2 )
        v4 = f( t[n]+h, y[n]+h*v3 )
        y[n+1] = y[n] + h/6*( v1 + 2*v2 + 2*v3 + v4 )

    EY = Y(t)
    err = h**(1/2)*np.linalg.norm(EY-y)
    plt.plot(t,EY, t,y,'r--'); plt.pause(.5)
    print(' RK4-err = ', '%12.6e' %err, '  with h = ', h)
    return err

def f(t,y):
    return -y + 2*np.cos(t)
#    return y/4*( 1- y/20 )

def Y(t):
    return np.sin(t) + np.cos(t)
#    return 20/( 1+19*np.exp(-t/4) )

# Solve
y0 = 1; a= 0; b = np.pi;
#y0 = 1; a= 0; b = 30;
plt.clf()
e1 = []; e2 = []
for k in range(1,10):
    h = 2**(-k)
    e1.append( RK2(f,h,y0,a,b) )

e1 = np.array(e1)
print('RK2 감축 비율 = ', e1[:-1]/e1[1:])

for k in range(1,10):
    h = 2**(-k)
    e2.append( RK4(f,h,y0,a,b) )

e2 = np.array(e2)
print('RK4 감축 비율 = ', e2[:-1]/e2[1:])

DATE: 2018.12.05 - 12:51

class 예제

파이썬 기초문법(강의용)