Odeint 在多个时期不一致，建模驱动摆

Question

似乎有很多关于摆建模和 odeint 的问题。我相信这个问题足够具体，可以独立存在。它涉及将时间数组传递给 odeint。

我正在建模一个受驱动的阻尼摆。我预计瞬态行为会在一段时间后消失，并且一直在绘制我的角速度与时间的关系来观察这一点。 我的问题是，改变周期数似乎并不能提供一致的结果。我看不出代码或我的假设在哪里失败。

from numpy import *
from scipy.integrate import odeint
import matplotlib.pyplot as plt

#pendulum diff eq
def pendulum(y,t,b,gamma,drivefreq):
    phi,omega=y
    dydt = [omega,-b*omega - sin(phi) + g*cos(drivefreq*t)]
    return dydt

#pendulum parameters: dampening, force amplitude, drivefreq
b=0.05;     g=0.4;      drivefreq=0.7 
args=(b,g,drivefreq) 

#num pts per period, num periods, time array
N=256;  nT=200;
t=linspace(0,nT*2*pi/drivefreq,nT*N)

上面的行有问题吗？在这里使用非整数值是不好的形式吗？ Linspace 仍然应该给出一个恒定间隔的数组。我已经看到其他示例成功地做到了这一点......我的想法是根据驾驶时间段的时间并设置一些数字，即每时间段 256 分。这有问题吗？

#initial conditions
y0= [0,0] #[phi0,omega0]

#run odeint
out=odeint(pendulum, y0,t,args)
omega  = out[:,1]

#plot ang velocity vs time
fig=plt.figure('ang velocity vs time')
plt.plot(t,omega)

下面是周期数 (nT) 等于 140,180 和 200 的图。我希望​​看到相同行为的继续，但是 180 周期结果不会失去其瞬态，并且 200 结果达到稳定最快的状态行为！我的逻辑哪里错了？

Answer 1

您的 Lipschitz 常数约为 L=1，它给出的时间差 dT 的误差放大系数为 exp(L*dT)=exp(dT)。仅考虑大约1e-16 的正常数值噪声，它只需要dT=37 将这个初始误差放大到大约1 的贡献，如exp(37)*1e-16 = 1.17。

如您所见，在从 0 到 1200 或更大的时间跨度内，即使算法执行中的最轻微变化也会导致轨迹看似随机变化。在从 0 到大约 30 的时间跨度下，您只能保证在这些程序变化下至少具有图形相似性。