5

我正在尝试使用 Python 研究以下延迟微分方程的行为：

y''(t) = -y(t)/τ^2 - 2y'(t)/τ - Nd*f(y(t-T))/τ^2,

其中f是一个截止函数，当其参数的绝对值在 1 和 10 之间时，它基本上等于恒等式，否则等于 0（见图 1），并且Nd和τ是T常数。

图 1：函数 f 的图

为此，我使用包 JiTCDDE。这为上述方程提供了一个合理的解。然而，当我尝试在等式右侧添加噪声时，我得到的解在几次振荡后稳定到非零常数。这不是方程的数学解（唯一可能的常数解等于零）。我不明白为什么会出现这个问题以及是否有可能解决它。

我在下面重现我的代码。在这里，为了简单起见，我用高频余弦代替了噪声，它被引入方程系统作为虚拟变量的初始条件（余弦可以直接引入系统，但对于一般噪音，这似乎是不可能的）。为了进一步简化问题，我还删除了涉及f函数的术语，因为没有它也会出现问题。图 2 显示了代码给出的函数图。

图 2：代码给出的解决方案图

from jitcdde import jitcdde, y, t
import numpy as np
from matplotlib import pyplot as plt
import math
from chspy import CubicHermiteSpline


# Definition of function f:
def functionf(x):
    return x/4*(1+symengine.erf(x**2-Bmin**2))*(1-symengine.erf(x**2-Bmax**2))

#parameters:
τ = 42.9
T = 35.33
Nd = 8.32

# Definition of the initial conditions:
dt = .01  # Time step.
totT = 10000.  # Total time.
Nmax = int(totT / dt)  # Number of time steps.
Vt = np.linspace(0., totT, Nmax)  # Vector of times.

# Definition of the "noise"
X = np.zeros(Nmax)
for i in range(Nmax):
    X[i]=math.cos(Vt[i])
past=CubicHermiteSpline(n=3)
for time, datum in zip(Vt,X):
    regular_past = [10.,0.]
    past.append((
        time-totT,
        np.hstack((regular_past,datum)),
        np.zeros(3)
        ))
noise= lambda t: y(2,t-totT)


# Integration of the DDE
g = [
     y(1),
     -y(0)/τ**2-2*y(1)/τ+0.008*noise(t)
     ]
g.append(0)


DDE = jitcdde(g) 
DDE.add_past_points(past)   
DDE.adjust_diff()

data = []
for time in np.arange(DDE.t, DDE.t+totT, 1):
    data.append( DDE.integrate(time)[0] )

plt.plot(data)
plt.show()

顺便说一句，我注意到即使没有噪音，解决方案在零点似乎也是不连续的（对于负时间，y 设置为零），我不明白为什么。

4

1 回答 1

2

正如评论所揭示的那样，您的问题最终归结为：

step_on_discontinuities假设延迟相对于集成时间很小，并执行放置在延迟组件指向集成开始的时间（0在您的情况下）的步骤。以这种方式处理初始不连续性。

totT但是，在您的情况下，使用延迟的虚拟变量实现输入会给系统带来很大的延迟。相应的步骤step_on_discontinuities将在totT其本身，即在所需的积分时间之后。因此，当您进入for time in np.arange(DDE.t, DDE.t+totT, 1):代码时，DDE.t是totT. 因此，在您真正开始整合和观察之前，您已经迈出了一大步，这可能看起来像是不连续性并导致奇怪的结果，特别是您看不到输入的效果，因为此时它已经“结束”了。为避免这种情况，请使用adjust_difforintegrate_blindly代替step_on_discontinuities.

于 2020-06-24T09:03:48.770 回答