Реализация 3/8 Runge Kutta в Python с ODE 2-го порядка

Я хочу решить уравнение y'' + 5y' + 6y = cos(t). Поскольку это 2-й порядок, я сначала создал систему ОДУ первого порядка, где dy/dt = z = f(t,y,z) и dz/dt = cos(t) - 5z - 6y = g(t,y, я). Я новичок в python и не знаю, как именно реализовать программу с системой ODE, но для входных данных я написал func = [f,g] и ic = [y0,dy0].

Во-вторых, поскольку Рунге Кутта оценивает выражение переменных, значения которых доступны только во время выполнения, я использую функцию eval python, хотя и получаю сообщение об ошибке. Ниже мой код и сообщение об ошибке:

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

Спасибо заранее за любые предложения; Я действительно пытаюсь сделать эту работу.