scipy.solve_ivp：在第二次发生事件时终止轨迹

Question

我想求解 return 到初始条件所需的时间。作为玩具示例，让我们以一个简单的谐波振荡器 d/dt (x, y) = (y, -x) 为例，初始条件 (x,y) = (1,0).

我希望代码为 return T = 6.28... (= 2π)，因为这是 return 到 y 为正的初始条件所花费的时间。当达到这一点时，它应该停止积分。

import numpy as np
from scipy.integrate import solve_ivp

def f(t,x):
    return [x[1], -x[0]]

def event(t,x):
    return x[0]

event.direction = 1
event.terminal = True # Try false, this gives the result but integrates beyond the event I need.

initial_condition = [0, 1]

sol = solve_ivp(f, [0, 10], initial_condition, method='RK45',dense_output=True, \
                events = event, rtol=1e-6, atol=3e-6)

print(sol.t_events)

问题是设置 event.terminal = True 会立即启动积分，因为事件在 t = 0 时得到满足。有没有什么办法可以让积分仅在事件满足两次后停止？

Answer 1

为了将来参考，这里是我一直在使用的解决方法：

将初始条件设置为事件前的一小部分，例如initial_condition = [1e-100, 1] 而不是 [0, 1]。它似乎不会影响轨迹，即使对于更高的 degree-of-freedom 混沌系统也是如此，但它不会在 t=0 时触发事件。

另一种可行的解决方法是在没有事件的情况下集成一段时间，然后使用结果作为进一步集成的初始条件，这次有事件。我怀疑这个答案对任何人都有用，但我发布这个只是为了完整性。

Answer 2

与您的解决方法类似，我将 1e-100 添加到事件而不是初始状态。这样做的好处是完全不影响轨迹，甚至在理论上也不影响。

当然，这确实意味着您未来的事件会提前一小部分发生，但如果这对您很重要，则取决于您的情况。

对于您编写的示例代码，这将使您的事件如下：

def event(t,x):
    return x[0] + 1e-100

我试过的一件事没有工作是在集成期间修改终止参数。我想如果我从 event.terminal = False 开始，然后在事件定义集 event.terminal = True 内，那么它将在第一个事件期间（在 t=0）继续，并在第二个事件处终止。但是，如果在集成过程中更改了参数，它似乎会忽略更改的参数。

（此问题与相关，但我没有标记任何内容的名誉。）

scipy.solve_ivp：在第二次发生事件时终止轨迹

scipy.solve_ivp: Terminate trajectory on second occurrence of event

python

solver

scipy

numerical-integration