使用二阶 ODE 在 Python 中实现 3/8 Runge Kutta

Question

我想求解方程 y'' + 5y' + 6y = cos(t)。由于这是二阶的，我首先创建了一个一阶 ODE 系统，其中 dy/dt = z = f(t,y,z) 和 dz/dt = cos(t) - 5z - 6y = g( t、y、z）。我是 python 的新手，不确定如何使用 ODE 系统准确地实现程序，但对于我的输入，我写了 func = [f,g] 和 ic = [y0,dy0].

其次，由于 Runge Kutta 评估其值仅在运行时可用的变量表达式，因此我正在使用 python 的 eval 函数，尽管我收到错误。下面是我的代码和错误信息：

def f(t,y,z): 
    return z
def g(t,y,z): 
    return np.cos(t)-5*z-6*y
t = np.linspace(0,20,100)
y = np.zeros((1,len(t)))
z = np.zeros((1,len(t)))

fun = [f(t,y,z),g(t,y,z)]
ic = [1,0] # y0,dy0

def ruku(fun,h):
    t0=0
    tf=20
    
    t=t0
    y=ic
    fc=0

    while t < tf:
        if t+h > tf:
            exit()
        if h == tf-t:
            exit()
        k1=eval('fun',t,y)
        k2=eval('fun',t+h/3,y+h*k1/3)
        k3=eval('fun',t+2*h/3,y+h*(k2-k1/3))
        k4=eval('fun',t+h,y+h*(k3-k2+k1))
        
        y=y+h*(k1+3*(k2+k3)+k4)/8
        fc=fc+4
        t=t+h
    return y

ruku(fun,0.01)

-->TypeError: globals must be a dict

在此先感谢您的任何建议；我真的很努力让这个工作。

Answer 1

原来我低估了Python的简单性！将我的循环替换为：

    while t +h <= tf:
            
        k1=fun(t,y)
        k2=fun(t+h/3,y+h*k1/3)
        k3=fun(t+2*h/3,y+h*(k2-k1/3))
        k4=fun(t+h,y+h*(k3-k2+k1))
        
        y=y+h*(k1+3*(k2+k3)+k4)/8
        fc=fc+4
        t=t+h

并将我的输入转换为 np.arrays

使用二阶 ODE 在 Python 中实现 3/8 Runge Kutta

Implementing the 3/8 Runge Kutta in Python with a 2nd order ODE

python

global

typeerror

runge-kutta