scipy.integrate.ode.set_solout工作吗？

Question

scipy.integrate.ode 接口提供了一种方法，用于在任何步骤上违反约束条件时停止积分，称之为 set_solout。然而，即使在最简单的示例中，我也无法让此方法正常工作。以下是一个尝试：

import numpy as np
from scipy.integrate import ode

def f(t, y):
    """Exponential decay."""
    return -y

def solout(t, y):
    if y[0] < 0.5:
        return -1
    else:
        return 0

y_initial = 1
t_initial = 0

r = ode(f).set_integrator('dopri5') # Integrator that supports solout
r.set_initial_value(y_initial, t_initial)
r.set_solout(solout)

# Integrate until t = 5, but stop when solout constraint violated
r.integrate(5)

# The time when solout should have terminated integration:
intersection_time = np.log(2)

当 t = log(2) = 0.693... 时，积分应该被solout停止，但实际上直到 t = 5，y = 0.007 才停止。这是 scipy 的一个bug吗？还是我没有正确使用 set_solout？

- Ted Pudlik

Answer 1

事实证明，在调用set_initial_value之前需要先调用set_solout。（我是通过研究scipy测试套件中set_solout测试的代码找到的这个问题。）因此，将我提出的代码中这两个调用的顺序颠倒一下就可以得到正确的结果。

即使这个行为是正确的，文档中也应该提到set_solout的这个用法。我已经在GitHub上向SciPy提交了一个问题。

更新：在SciPy 0.17.0中修复了这个问题；即使在调用set_initial_value之后调用set_solout，set_solout仍然可以正常工作，而且问题代码将产生正确的结果。

- Ted Pudlik