scipy中具有时变系数的ODE答案

【问题标题】：ODE with time-varying coefficients in scipyscipy中具有时变系数的ODE
【发布时间】：2021-07-12 08:34:27
【问题描述】：

我正在评估一组具有时变系数的 ODE

def deriv(y, t, N, coefficients):
    S, I, R = y
    dSdt = coefficients['beta'](t) * S * I / N * -1
    dIdt = coefficients['beta'](t) * S * I / N - coefficients['gamma']* I
    dRdt = coefficients['gamma'] * I
    return dSdt, dIdt, dRdt

特别是，我在一个预先计算的数组中有“beta”值，其大小等于 int(max(t))。

coefficients = {'beta' : beta_f,'gamma':0.1}
def beta_f(t):
     return  mybetas.iloc[int(t)]

# Initial conditions vector
y0 = (S0, I0, R0)
# Integrate the SIR equations over the time grid, t.
ret = odeint(deriv,y0,t,args=(N,coefficients))

当我运行 odeint 时，它还会评估超出 max(t) 的值，从而在 beta_f 中引发索引越界错误。如何限制 odeint 的评估跨度？

【问题讨论】：

标签： python numpy scipy ode

【解决方案1】：

从len(mybetas) == int(max(t)) 开始，即使 t 的值不超过 max(t)，也会出现越界错误。例如，mybetas.iloc[int(max(t))] 会给出越界错误，即使 int(max(t)) <= max(t) 的正值是 t。

但就您而言，odeint 确实检查了集成域之外的一些值。就在几周前，我不得不处理一个与您类似的问题，以下关于 stackoverflow 的两次讨论非常有帮助：

integrate.ode sets t0 values outside of my data range

Solve ODEs with discontinuous input/forcing data

第二个链接解释了为什么使用odeint 在for 循环中一个接一个地在每个单独的整数时间步上求解 ODE，而不是让 odeint 处理您的由您的 beta 值跳跃引起的导数。

否则，如果这适合您的研究案例，您可以插值您的 beta，并让函数 beta_f 返回 beta 的插值。当然，您必须将插值域稍微扩展到积分域之外，因为odeint 可能想要评估一些大于max(t) 的t 的导数。

【讨论】：