我正在使用 scipy 的solve_ivp来求解一个刚性微分方程。我正在使用方法BDF来解决这个问题。我已经使用了 MATLAB 的ode23s,并且在 MATLAB 中得到了正确的结果。但是,solve_ivp 返回错误的结果。问题可能在于计算X。我的代码是:
X0=np.array([0,0,0,0,0,0])
swd=0
start =timer()
for k in range(0,100):
tspan=[k*(1e-4), (k+.7)*(1e-4)]
tspan1=[(k+.7)*(1e-4), (k+1)*(1e-4)]
swm=1
sol=solve_ivp(conv,tspan,X0,method='BDF')
TT=sol.t
YY=sol.y
swm=0
sol1=solve_ivp(conv,tspan1,YY[:,-1],method='BDF')
TT1=sol1.t
YY1=sol1.y
X0=YY1[:,-1]
plt.plot(TT,YY[0,:],'r')
plt.plot(TT,YY[1,:],'b')
end=timer()
print(end-start)
对于定义为的函数
def conv(t,X):
global swd
if (swd==0) or (swd==1 and X[0]>0):
Xdot = A11.dot(X).reshape(6,1) + B11.dot(U)
swd=1
else:
Xdot = A10.dot(X).reshape(6,1) + B10.dot(U)
swd=0
return np.squeeze(np.asarray(Xdot))
矩阵A是6 6* ,B分别是6 2* 。X是6 1*。X值之一计算不正确,然后整个图返回不正确。仅正确返回了第一个周期的图。 详细信息:给定的代码是升压转换器。swm和swd分别是 MOSFET 和二极管的开关状态。据我调试,二极管在第一个周期正确切换。但是对于第二个循环和其他循环,二极管不会关闭,而是保持打开状态。