Scipy solve_ivp 非常慢/卡顿

Question

我最近在解决一个问题，涉及到两个成对的常微分方程（ode）。有人建议我使用scipy库里的solve_ivp来解决这个问题，但我发现这个方法要么卡住，要么运行得非常慢（我试过让它运行半个小时，结果还是没有任何结果）。

我在数值解常微分方程方面的经验有限，而且之前从来没有用过solve_ivp，所以可能是我哪里出错了，但我一直找不到问题所在。

import numpy as np
from matplotlib import pyplot as plt
from scipy import integrate as sc

param_dict = {
    'r1': 0.36,
    'r2': 0.15,
    'k1': 4000000,
    'b12': 0.0000002,
    'b21': 0.69,
    'p1': 80,
    'p2': 275,
    'q1': 0.0002,
    'q2': 0.0004,
    'phi1': 0.00095,
    'phi2': 0.000075,
    'c1': 20000,
    'c2': 30000
}

# INITIAL VALUES
init_X = np.array([4000000,275000])
init_E = np.array([800,200])
y0 = np.concatenate((init_X, init_E))

def functions(t, X, r1, r2, k1, b12, b21, q1, q2, phi1, phi2, c1, c2, p1, p2):
    x1, x2, y1, y2 = X
    dx1dt = (r1 * x1) * (1 - (x1 / k1) - (b12 * x2)) - (q1 * y1 * x1)
    dx2dt = (r2 * x2) * (1 - (x2 / (b21 * x1))) - (q2 * y2 * x2)
    dy1dt = phi1 * y1 * (p1*q1 * x1 - c1)
    dy2dt = phi2 * y2 * (p2*q2 * x2 - c2)
    return [dx1dt, dx2dt, dy1dt, dy2dt]

t_span = [0,20]
t = np.linspace(0,100,1000)

res = sc.solve_ivp(functions, t_span, y0, args=tuple(param_dict.values()), dense_output=True)

zz = res.sol(t)
x1 = []
x2 = []
e1 = []
e2 = []
for ii in range(len(zz.T)):
    x1.append(zz.T[ii,0])
    x2.append(zz.T[ii,1])
    e1.append(zz.T[ii,2])
    e2.append(zz.T[ii,3])

plt.plot(x1)
plt.plot(x2)

如果我把init_E = [0,0]设置成这个值，它就能正常运行，但只要换成其他值，就会出现卡住或者运行极慢的问题。

这个模型是关于捕食者和猎物的关系，主要用于渔业管理，感兴趣的朋友可以了解一下。

Answer 1

我觉得问题可能出在你的字典里的键（也就是里面的内容）和你函数里的参数（也就是你传给函数的东西）顺序不一致。

如果我这样做，就能立刻解决这个问题。

param_names = ['r1', 'r2', 'k1', 'b12', 'b21', 'q1', 'q2', 'phi1', 
               'phi2', 'c1', 'c2', 'p1', 'p2']
args = tuple(param_dict[p] for p in param_names)
res = sc.solve_ivp(functions, t_span, y0, args=args, dense_output=True)

（我不确定这是不是正确的解决办法，但它不会卡住）。

给点建议：

• 一般来说，在把你的函数交给求解器之前，最好先调用一次这个函数，看看它的输出结果。比如在这样的初始条件下：

f0 = functions(0, y0, *args)
assert np.array_equal(f0, [-719200.0, 15139.945652173912, 33440.0, 3.75])