在 ode 求解器中使用的步长 - python
Step size to use in ode solver - python
我正在使用 scipy.integrate
中的 ode
求解器来求解我的微分方程。我想看看最终结果是否受到积分中步长选择的影响dt
,这是我得到的:
我希望 any 小积分步长的结果是正确的,但我看到的几乎相反...有人知道发生了什么吗?
--
代码:
from scipy.integrate import ode
import matplotlib.pyplot as plt
import numpy as np
y0, t0 = 0., 0
def f2(t, y, arg1):
t0 = 5
shape = np.piecewise(t, [t<t0, t>=t0], [arg1*t, arg1*(2*t0-t)])
return shape
t1 = 10.
dtlist = np.array([t1/j for j in range(1, 50)])
temp = []
for i in dtlist:
dt = i
r = ode(f2).set_integrator('zvode', method='bdf')
r.set_initial_value(y0, t0).set_f_params(2.0)
while r.successful() and r.t < t1:
r.integrate(r.t+dt)
temp.append(r.y)
fig, axs = plt.subplots(1, 1, figsize=(15, 6), facecolor='w', edgecolor='k')
xaxis = dtlist
yaxis = temp
axs.plot(xaxis, yaxis, 'bo', label = 'obtained')
axs.plot(xaxis, xaxis - xaxis + 50, 'r--', label = 'expected')
axs.set_xlabel('time increment $\delta t$')
axs.set_ylabel('Final $y$ value')
axs.legend(loc = 4)
您需要确保最后的积分步骤总是在 t1
结束。由于浮点错误,可能会发生所需的最后一个积分步骤 r.t
略小于 t1
,因此执行了大约 t1+dt
的额外不需要的步骤。你可以使用
r.integrate(min(r.t+dt, t1))
适时切断积分。或者做一些更复杂的事情来捕捉 r.t+dt
已经接近 t1
.
的情况
我正在使用 scipy.integrate
中的 ode
求解器来求解我的微分方程。我想看看最终结果是否受到积分中步长选择的影响dt
,这是我得到的:
我希望 any 小积分步长的结果是正确的,但我看到的几乎相反...有人知道发生了什么吗?
--
代码:
from scipy.integrate import ode
import matplotlib.pyplot as plt
import numpy as np
y0, t0 = 0., 0
def f2(t, y, arg1):
t0 = 5
shape = np.piecewise(t, [t<t0, t>=t0], [arg1*t, arg1*(2*t0-t)])
return shape
t1 = 10.
dtlist = np.array([t1/j for j in range(1, 50)])
temp = []
for i in dtlist:
dt = i
r = ode(f2).set_integrator('zvode', method='bdf')
r.set_initial_value(y0, t0).set_f_params(2.0)
while r.successful() and r.t < t1:
r.integrate(r.t+dt)
temp.append(r.y)
fig, axs = plt.subplots(1, 1, figsize=(15, 6), facecolor='w', edgecolor='k')
xaxis = dtlist
yaxis = temp
axs.plot(xaxis, yaxis, 'bo', label = 'obtained')
axs.plot(xaxis, xaxis - xaxis + 50, 'r--', label = 'expected')
axs.set_xlabel('time increment $\delta t$')
axs.set_ylabel('Final $y$ value')
axs.legend(loc = 4)
您需要确保最后的积分步骤总是在 t1
结束。由于浮点错误,可能会发生所需的最后一个积分步骤 r.t
略小于 t1
,因此执行了大约 t1+dt
的额外不需要的步骤。你可以使用
r.integrate(min(r.t+dt, t1))
适时切断积分。或者做一些更复杂的事情来捕捉 r.t+dt
已经接近 t1
.