tig*_*gro 5 python-3.x differential-equations numerical-stability jitcode-jitcdde-jitcsde
我正在尝试使用 Python 研究以下延迟微分方程的行为:
\ny\'\'(t) = -y(t)/\xcf\x84^2 - 2y\'(t)/\xcf\x84 - Nd*f(y(t-T))/\xcf\x84^2,\nRun Code Online (Sandbox Code Playgroud)\n其中f是截止函数,当其参数的绝对值在 1 到 10 之间时,它本质上等于恒等式,否则等于 0(见图 1),并且 、 和Nd是\xcf\x84常数T。

为此,我使用 JiTCDDE 包。这为上式提供了合理的解。尽管如此,当我尝试在方程右侧添加噪声时,我得到了一个在几次振荡后稳定为非零常数的解。这不是方程的数学解(唯一可能的常数解等于零)。我不明白为什么会出现这个问题以及是否可以解决它。
\n我在下面重现我的代码。这里,为了简单起见,我用高频余弦代替噪声,将其引入方程组作为虚拟变量的初始条件(余弦可以直接引入系统,但对于一般噪音这似乎不可能)。为了进一步简化问题,我还删除了涉及该函数的术语f,因为没有它也会出现问题。图 2 显示了代码给出的函数图。

from jitcdde import jitcdde, y, t\nimport numpy as np\nfrom matplotlib import pyplot as plt\nimport math\nfrom chspy import CubicHermiteSpline\n\n\n# Definition of function f:\ndef functionf(x):\n return x/4*(1+symengine.erf(x**2-Bmin**2))*(1-symengine.erf(x**2-Bmax**2))\n\n#parameters:\n\xcf\x84 = 42.9\nT = 35.33\nNd = 8.32\n\n# Definition of the initial conditions:\ndt = .01 # Time step.\ntotT = 10000. # Total time.\nNmax = int(totT / dt) # Number of time steps.\nVt = np.linspace(0., totT, Nmax) # Vector of times.\n\n# Definition of the "noise"\nX = np.zeros(Nmax)\nfor i in range(Nmax):\n X[i]=math.cos(Vt[i])\npast=CubicHermiteSpline(n=3)\nfor time, datum in zip(Vt,X):\n regular_past = [10.,0.]\n past.append((\n time-totT,\n np.hstack((regular_past,datum)),\n np.zeros(3)\n ))\nnoise= lambda t: y(2,t-totT)\n\n\n# Integration of the DDE\ng = [\n y(1),\n -y(0)/\xcf\x84**2-2*y(1)/\xcf\x84+0.008*noise(t)\n ]\ng.append(0)\n\n\nDDE = jitcdde(g) \nDDE.add_past_points(past) \nDDE.adjust_diff()\n\ndata = []\nfor time in np.arange(DDE.t, DDE.t+totT, 1):\n data.append( DDE.integrate(time)[0] )\n\nplt.plot(data)\nplt.show()\nRun Code Online (Sandbox Code Playgroud)\n顺便说一句,我注意到即使没有噪声,解在零点似乎也是不连续的(负数时 y 设置为零),我不明白为什么。
\n随着评论的揭晓,你的问题最终归结为:
\nstep_on_discontinuities假设相对于积分时间的延迟很小,并执行放置在延迟组件指向积分开始的时间上的步骤(0在您的情况下)。通过这种方式处理初始不连续性。
totT但是,在您的情况下,使用延迟虚拟变量实现输入会给系统带来很大的延迟。\n相应的步骤step_on_discontinuities将在totT其本身,即在所需的积分时间之后。\n因此,当您进入for time in np.arange(DDE.t, DDE.t+totT, 1):代码时,DDE.t是totT\n因此,在真正开始集成和观察之前,您已经迈出了一大步,这可能看起来像是不连续性并导致奇怪的结果,特别是您看不到输入的效果,因为它已经 \xe2\x80\x9cending \xe2\x80\x9d 此时。\n要避免这种情况,请使用adjust_diff或integrate_blindly代替step_on_discontinuities。
| 归档时间: |
|
| 查看次数: |
127 次 |
| 最近记录: |