scipy.integrate.ode 放弃集成
Posted
技术标签:
【中文标题】scipy.integrate.ode 放弃集成【英文标题】:scipy.integrate.ode gives up on integration 【发布时间】:2021-03-13 09:58:40 【问题描述】:我遇到了scipy.integrate.ode 的一个奇怪问题。这是一个最小的工作示例:
import sys
import time
import numpy as np
import matplotlib.pyplot as plt
import matplotlib
from scipy.integrate import ode, complex_ode
def fun(t):
return np.exp( -t**2 / 2. )
def ode_fun(t, y):
a, b = y
f = fun(t)
c = np.conjugate(b)
dt_a = -2j*f*c + 2j*f*b
dt_b = 1j*f*a
return [dt_a, dt_b]
t_range = np.linspace(-10., 10., 10000)
init_cond = [-1, 0]
trajectory = np.empty((len(t_range), len(init_cond)), dtype=np.complex128)
### setup ###
r = ode(ode_fun).set_integrator('zvode', method='adams', with_jacobian=False)
r.set_initial_value(init_cond, t_range[0])
dt = t_range[1] - t_range[0]
### integration ###
for i, t_i in enumerate(t_range):
trajectory[i,:] = r.integrate(r.t+dt)
a_traj = trajectory[:,0]
b_traj = trajectory[:,1]
fun_traj = fun(t_range)
### plot ###
plt.figure(figsize=(10,5))
plt.subplot(121, title='ODE solution')
plt.plot(t_range, np.real(a_traj))
plt.subplot(122, title='Input')
plt.plot(t_range, fun_traj)
plt.show()
此代码工作正常,输出图为(ODE 显式依赖于输入变量,右侧面板显示输入,左侧面板显示第一个变量的 ode 解决方案)。
所以原则上我的代码是有效的。奇怪的是,如果我只是简单地替换积分范围
t_range = np.linspace(-10., 10., 10000)
通过
t_range = np.linspace(-20., 20., 10000)
我得到了输出
因此,集成商不知何故放弃了集成,并将我的解决方案保持不变。 为什么会这样?我该如何解决?
我测试过的一些东西:这显然不是分辨率问题,集成步骤已经很小了。相反,似乎集成商在几步之后甚至不再费心调用 ode 函数了。我已经通过在 ode_fun() 中包含一个 print 语句进行了测试。
我目前的怀疑是,在最初的几个集成步骤中,我的函数没有发生显着变化,因此集成商认为我的函数是恒定的。我可能需要在某个地方设置一些容忍度吗?
任何帮助表示赞赏!
【问题讨论】:
【参考方案1】:“我目前的怀疑是,在最初的几个集成步骤中,我的函数没有发生显着变化后,集成商认为我的函数是恒定的。”你的怀疑是正确的。 ODE 求解器通常具有内部步长,该步长根据求解器计算的误差估计进行自适应调整。这些步长与请求输出的时间无关;使用在内部步骤计算的点处的解插值计算请求时间的输出。
当您在 t = -20
启动求解器时,显然输入变化如此缓慢,以至于求解器的内部步长变得足够大,以至于当求解器接近 t = 0
时,求解器会直接跳过输入脉冲。
您可以使用set_integrator
方法的选项max_step
限制内部步长。如果我将max_step
设置为2.0(例如),
r = ode(ode_fun).set_integrator('zvode', method='adams', with_jacobian=False,
max_step=2.0)
我得到了你期望的输出。
【讨论】:
啊,我明白了!感谢您的回答,这解决了问题。也许是一个幼稚的后续问题:是否有办法强制积分器“不假设我的函数是恒定的”? (例如,用于自动扫描不同的参数)也许不同的积分器会做?还是我必须想办法确保选择 max_step 以免发生这种情况? 如果您的输入包含脉冲或尖峰(例如您的fun
),我认为您需要想办法将max_step
设置得足够小,以使求解器能够“看到”脉搏。
非常感谢!实际上,我发现这种违反直觉的做法是,我最初假设“看到脉冲的求解器”将由积分网格决定。显然情况并非如此!感谢您为我解决这个问题。以上是关于scipy.integrate.ode 放弃集成的主要内容,如果未能解决你的问题,请参考以下文章
scipy.integrate.ode 与两个耦合的 ODE?
scipy.integrate.ode.set_solout 有效吗?
通过 scipy.integrate.ode 使用自适应步长