1

有谁知道将跟踪下图的 Mathematica 代码？

这是该图的方程，具有常数系数的二阶线性微分方程：

在此处输入图像描述

下面是这个等式所描绘的图表：

在此处输入图像描述

引用《时间序列分析和示例预测》一书：

... 其中 δ(t ) 是一个脉冲 (delta) 函数，它像豌豆弹一样，在时间 t = 0 迫使钟摆远离其平衡，a 是豌豆冲击的大小。很容易想象，这个二阶微分方程所描绘的曲线是时间的阻尼正弦函数，尽管如果摩擦或粘度足够大，（过阻尼）钟摆可能会沿着指数曲线逐渐停止而不会交叉中心线。

4

1 回答 1

6 
eq = m z''[t] + c z'[t] + k z[t] == a DiracDelta[t];
parms = {m -> 1, c -> .1, k -> 1, a -> 1};
sol = First@DSolve[{eq /. parms, z[0] == 1, z'[0] == 0}, z[t], t];
Plot[z[t] /. sol, {t, 0, 70}, PlotRange -> All, Frame -> True, 
 FrameLabel -> {{z[t], None}, {Row[{t, " (sec)"}], eq}}, 
 GridLines -> Automatic]

数学图形

请注意，对于零初始条件，另一种选择是使用 Mathematica 中的控制系统函数，如下所示

parms = {m -> 10, c -> 1.2, k -> 4.3, a -> 1};
tf = TransferFunctionModel[a/(m s^2 + c s + k) /. parms, s]
sol = OutputResponse[tf, DiracDelta[t], t];

Plot[sol, {t, 0, 60}, PlotRange -> All, Frame -> True, 
 FrameLabel -> {{z[t], None}, {Row[{t, " (sec)"}], eq}}, 
 GridLines -> Automatic]

数学图形

更新

严格来说，DSolve上面的结果并不是这个问题的手工推导可以找到的。正确的解决方案应该如下

（也请参阅此内容以供参考）

正确的解析解由下式给出

数学图形

我在这里（第一章）针对这个问题和类似情况得出了这个结论。

使用上述解决方案，正确的响应将如下所示：

parms = {m -> 1, c -> .1, k -> 1, a -> 1};
w = Sqrt[k/m];
z = c/(2 m w);
wd = w Sqrt[1 - z^2];
analytical = 
  Exp[-z w t] (u0 Cos[wd t] + (v0 + (u0 z w))/wd Sin[wd t] + 
     a/(m wd) Sin[wd t]);
analytical /. parms /. {u0 -> 1, v0 -> 0}

 (* E^(-0.05 t) (Cos[0.998749 t] + 1.05131 Sin[0.998749 t]) *)

绘制它：

Plot[analytical /. parms /. {u0 -> 1, v0 -> 0}, {t, 0, 70}, 
 PlotRange -> All, Frame -> True, 
 FrameLabel -> {{y[t], None}, {Row[{t, " (sec)"}], 
    "analytical solution"}}, GridLines -> Automatic, ImageSize -> 300]

数学图形

如果您将上面的图与上面显示的第一个图进行比较，DSolve您可以看到 附近的差异t=0。

于 2013-01-06T02:13:52.560 回答