wolfram-mathematica - NDS 求解波动方程时的不稳定性

Question

我正在尝试使用 NDSolve求解波动方程以检查使用它是否更容易和/或更快，而不是我的旧特性方程。方法实现。

我得到了很多使用特征方法没有得到的不稳定性，而且由于这些是简单的方程，我想知道有什么问题......(希望不是问题的物理方面......)

ans = Flatten@NDSolve[{
u[t, x]*D[d[t, x], x] + d[t, x]*D[u[t, x], x] + D[d[t, x], t] == 0,
D[d[t, x], x] + u[t, x]/9.8*D[u[t, x], x] + 
 1/9.8*D[u[t, x], t] + 0.0001 u[t, x]*Abs[u[t, x]] == 0,
u[0, x] == 0,
d[0, x] == 3 + x/1000*1,
u[t, 0] == 0,
u[t, 1000] == 0
},
d, {t, 0, 1000}, {x, 0, 1000}, DependentVariables -> {u, d}
]

Animate[Plot[(d /. ans)[t, x], {x, 0, 1000}, 
        PlotRange -> {{0, 1000}, {0, 6}}], {t, 0, 1000}
]

有人能帮我吗？

编辑:

我已经放置了 NDSolve解决方案(遵循 JxB 的编辑)与我的特性解决方案，一起放在同一个动画上。除了最初的快速振荡外，它们匹配得足够近。随着时间的推移，它们往往会开始不同步，但我相信这可能是由于我们在推断特征时必须承认的一个小简化。



红色: NDsolve ;蓝色:“手动”特性法；

按 F5(刷新浏览器)，从 t=0 重新开始动画.

(xx 比例是我在“手动”方法中使用的点数，其中每个点代表 NDSolve/物理比例的 20 个单位)

一起玩 NDSolve网格采样，呈现完全不同的振荡效果。有没有人拥有或知道确保正确集成的技术？

Answer 1

通过将系数更改为无限精度(例如，1/9.8->10/98)，并设置 WorkingPrecision->5 (值 6 太高)，我不再收到错误消息:

ans = Flatten@
  NDSolve[{D[u[t, x] d[t, x], x] + D[d[t, x], t] == 0, 
    D[d[t, x], x] + u[t, x] 10/98*D[u[t, x], x] + 
      10/98*D[u[t, x], t] + 1/10000 u[t, x]*Abs[u[t, x]] == 0, 
    u[0, x] == 0, d[0, x] == 3 + x/1000, u[t, 0] == 0, 
    u[t, 1000] == 0}, d, {t, 0, 1000}, {x, 0, 1000}, 
   DependentVariables -> {u, d}, WorkingPrecision -> 5]

Animate[
 Plot[(d /. ans)[t, x], {x, 0, 1000}, 
  PlotRange -> {{0, 1000}, {0, 6}}], {t, 0, 1000}]

我不知道这个方程，所以我不相信解决方案:小规模振荡最初增长，然后减弱。