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)

集成应该在solout时停止t = log(2) = 0.693...,但是很快就会继续,直到t = 5,何时y = 0.007.

这是一个错误scipy,还是我没有set_solout正确使用？

Answer 1

事实证明你需要在打电话set_solout 之前打电话set_initial_value.(我通过研究测试套件中的set_solout测试来解决这个问题scipy.)因此,在我的问题代码中反转两个调用的顺序会产生正确的结果.

即使这种行为是正确的,也应该在文档中提到set_solout.我在GitHub上发布了SciPy的问题.

更新:此问题已在SciPy 0.17.0中修复; set_solout即使被调用也会工作set_initial_value,问题代码将产生正确的结果.