我使用的是scipy.integrate.odeint
,当切换到bdf
方法时遇到了麻烦(因为系统看起来很僵硬)。不幸的是,bdf
需要计算雅可比行列式。我的系统有成千上万个方程,因此这会杀死Python(分段错误!)
我有充分的理由相信亚当斯在这个系统上会做的很好。所以我想强迫scipy使用它。我认为odeint
无法实现,因此我正在尝试使用ode
。我已经编写了自己的函数,该函数应该与odeint
相同。看起来我的代码不正确:
def my_odeint(dfunc, V0, times, args=()):
r = integrate.ode(dfunc)
r.set_integrator('vode', method='adams')
r.set_initial_value(V0,times[0]).set_f_params(*args)
V=[V0]
for time in times[1:]:
V.append(r.integrate(time))
V = scipy.array(V)
return V
当我测试它时,我得到一个错误:
def dx(X, t):
x=X[0]
y=X[1]
return (-x, x+y)
X0 = [1,1]
times = scipy.linspace(0,2,10)
my_odeint(dx, X0, times)
> capi_return is NULL
> Call-back cb_f_in_dvode__user__routines failed.
然后
TypeError Traceback (most recent call last)
<ipython-input-6-cc286a9956a2> in <module>()
----> 1 my_odeint(dx, X0, times)
in my_odeint(dfunc, V0, times, args)
1352 V=[V0]
1353 for time in times[1:]:
-> 1354 V.append(r.integrate(time))
in integrate(self, t, step, relax)
406 self._y, self.t = mth(self.f, self.jac or (lambda: None),
407 self._y, self.t, t,
--> 408 self.f_params, self.jac_params)
in run(self, f, jac, y0, t0, t1, f_params, jac_params)
863 args = ((f, jac, y0, t0, t1) + tuple(self.call_args) +
864 (f_params, jac_params))
--> 865 y1, t, istate = self.runner(*args)
<ipython-input-4-1770d3974ec9> in dx(X, t)
1 def dx(X, t):
----> 2 x=X[0]
3 y=X[1]
> TypeError: 'float' object has no attribute '__getitem__'
所以我的问题是-如何修复我的代码,以使我的函数接受与
odeint
相同的输入,但使用adams
并且从不切换到bdf
? 最佳答案
scipy.ode要求右侧为f(t, y, *f_args)
形式,odeint要求其右侧为f(y, t, *f_args)
形式。
所以这将工作:
def my_odeint(dfunc, V0, times, args=()):
r = integrate.ode(lambda t, X: scipy.array(dfunc(X, t, *args)))
r.set_integrator('vode', method='adams')
r.set_initial_value(V0,times[0])
V=[V0]
for time in times[1:]:
V.append(r.integrate(time))
V = scipy.array(V)
return V
请注意,根据文档,vode实现的是Adams-Moulton,而不是Adams-Bashforth。也就是说,您仍然有一个隐式方法。如果集成继续失败,请尝试在
integrate.ode
中也实现的Dormand / Prince方法。关于python - 如何使用Adams强制scipy进行整合?,我们在Stack Overflow上找到一个类似的问题:https://stackoverflow.com/questions/40317096/