使用solve_ivp代替odeint求解初始问题值

Question

目前，我使用 odeint 求解以下 ODE 方程组

dx/dt = (-x + u)/2.0

dy/dt = (-y + x)/5.0

初始条件：x = 0，y = 0

但是，我想使用 solve_ivp 这似乎是此类问题的推荐选项，但老实说我不知道如何调整代码......

这是我与 odeint 一起使用的代码：

import numpy as np

from scipy.integrate import odeint, solve_ivp

import matplotlib.pyplot as plt

def model(z, t, u):

    x = z[0]

    y = z[1]

    dxdt = (-x + u)/2.0

    dydt = (-y + x)/5.0

    dzdt = [dxdt, dydt]

    return dzdt

def main():

    # initial condition

    z0 = [0, 0]

    # number of time points

    n = 401

    # time points

    t = np.linspace(0, 40, n)

    # step input

    u = np.zeros(n)

    # change to 2.0 at time = 5.0

    u[51:] = 2.0

    # store solution

    x = np.empty_like(t)

    y = np.empty_like(t)

    # record initial conditions

    x[0] = z0[0]

    y[0] = z0[1]

    # solve ODE

    for i in range(1, n):

        # span for next time step

        tspan = [t[i-1], t[i]]

        # solve for next step

        z = odeint(model, z0, tspan, args=(u[i],))

        # store solution for plotting

        x[i] = z[1][0]

        y[i] = z[1][1]

        # next initial condition

        z0 = z[1]

    # plot results

    plt.plot(t,u,'g:',label='u(t)')

    plt.plot(t,x,'b-',label='x(t)')

    plt.plot(t,y,'r--',label='y(t)')

    plt.ylabel('values')

    plt.xlabel('time')

    plt.legend(loc='best')

    plt.show()

main()