Python中的自适应ODE算法

Question

我想在2D空间中整合一个粒子的路径，使用的是integrate.ode模块。我的情况有点特别，我只想整合到一个特定的位置，这个位置是由粒子允许的最大x坐标x_max决定的。

我遇到的主要问题是，粒子可能一开始移动得很慢，后来又加速。因此，我不想在这个区域浪费时间步长太小的计算。这个算法应该能够调整，当粒子的速度变快时，使用更小的时间步长。

为了这个目的，我下面有一些粗略的伪代码：

backend = "dopri5"
x_max = 1
solver = ode(f)
solver.set_integrator(backend)  
solver.set_initial_value(y0, t0)
t, y = [t0], [y0]
k = 1.2
while solver.successful() and solver.y[0] < x_max
   solver.integrate(solver.t+dt)
   t.append(solver.t)
   y.append(solver.y)
   v_current = numpy.linalg.norm(y[-1])
   v_previous = numpy.linalg.norm(y[-2])
   if numpy.abs( v_current-v_previous ) > k * v_previous:
       dt = 0.8*dt
       del y[-1]
   else:
       dt = dt*1.2

问题是这个算法可能不够稳健，因为选择k, 1.2, 0.8这些值有点随意，可能会导致算法的不稳定。

有没有人能建议一个更好的方法呢？

Answer 1

如果你想把一些离散事件（比如交叉点）和常微分方程（ODE）结合起来，可能需要试试一些专门为这类问题设计的库。看起来有几个不错的选择：

• PyDSTool
• Assimulo2
• scikits.ode

如果你的模型结构表明你总是能检测到事件，那你可以尝试大步走，并用二分查找来找到事件发生的时间或位置。不过，这种方法通常在多维情况下不太适用。

这类问题在游戏编程、实时物理模拟和碰撞检测中很常见，你可以尝试借鉴这些领域的技术（比如想象快速移动的子弹穿过慢速移动的物体）。

如果你想对某个事件做出反应，比如一个弹跳的球，那你实际上是在处理混合系统的问题。