目前,我使用 odeint 求解以下 ODE 方程组dx/dt = (-x + u)/2.0dy/dt = (-y + x)/5.0初始条件:x = 0,y = 0但是,我想使用 solve_ivp 这似乎是此类问题的推荐选项,但老实说我不知道如何调整代码......这是我与 odeint 一起使用的代码:import numpy as npfrom scipy.integrate import odeint, solve_ivpimport matplotlib.pyplot as pltdef model(z, t, u): x = z[0] y = z[1] dxdt = (-x + u)/2.0 dydt = (-y + x)/5.0 dzdt = [dxdt, dydt] return dzdtdef 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()
添加回答
举报
0/150
提交
取消