用欧拉方法处理 ODE 不起作用

Question

所以我必须用不同的方法和 Δt 来处理 ODE，但突然间我的函数停止工作，Δt<1。 ODE 是：T(t)'=R*T(t) - K*T(t-δ)，R 和 K 常量，十进制和小常量。

第一部分是在 δ=0 的情况下处理的，所以它很顺利，没有问题，但现在我在 δ=5 的情况下处理它，突然 ODE 就不起作用了。当 Δt=1 时，程序可以使用 Euler 方法很好地处理 ODE，但当 Δt=.1 时，它突然停止工作。

它没有崩溃或其他任何事情，但它没有正确计算数字，而是似乎无法识别 δ，或者假定它为零，然后继续这样，直到图形等于δ=0.

我正在使用 Numpy 和 Matplotlib，这是代码：

r = 0.1
R = 0.1
CeCw = 0.27
γ1 = 1
γ2 = 0.164
α = 0.612
δ = 5
Te0 = 1
tf = 240
Δt2 = np.arange(0, tf+.1, .1)

def d(Te, n):
    return R*Te[n] - α*γ1*CeCw*Te[n-δ]

Te2 = np.array([Te0])
for n in range(len(Δt2)-1):
    if n*.1<=δ:
        Te2 = np.append(Te2, 1)
    else:
        Te2 = np.append(Te2, Te2[n] + .1*d(Te2, n))

如果我将 .1 更改为简单的 1，它会突然再次开始工作（尽管它生成了一个不正确的图形，因为它不是正确的方法）。我尝试用 Decimal() 包装数字并使用不同的函数，但不知何故 .1 回来删除 δ。

编辑：向前移动并尝试使用 δ=5 实施 Heun 的方法，并且在 Δt=.1 或更低的情况下，它的行为就像 δ=0

Answer 1

这种方程通常称为延迟微分方程并且不是 ODE。专门的求解器通常使用足够高的误差阶插值来获取过去时间的值。

更具体地说，您的代码没有正确实现延迟。 delay δ 是时间差，不是索引差。对于索引差异 m_δ，您想要 m_δ*h = δ 的实现类型。然后用m_δ确定初始段的长度和导数函数中的偏移量。

这看起来像

δ = 5
mδ = 50
h = δ/mδ
Te0 = 1
tf = 240
Δt2 = np.arange(0, tf+h, h)

def d(Te, n):
    return R*Te[n] - α*γ1*CeCw*Te[n-mδ]

Te2 = np.array(len(Δt2)*[Te0], dtype=float)

for n in range(mδ,len(Δt2)-1):
    Te2[n+1] = Te2[n] + h*d(Te2, n)

这确实会产生振荡

使用时间循环的 Heun 方法获得相同的视觉图像

for n in range(mδ,len(Δt2)-1):
    k1 = d(Te2, n)
    Te2[n+1] = Te2[n] + h*k1
    k2 = d(Te2, n+1)
    Te2[n+1] = Te2[n] + 0.5*h*(k1+k2)

数值与小数点后第二位不同，正如一阶和二阶方法所预期的那样 h=0.1。

用欧拉方法处理 ODE 不起作用

crunching a ODE with Euler's Method doesn't work

python-3.x

numpy

differential-equations