使用 SciPy 以数值方式求解 ODE答案

【问题标题】：Solving an ODE Numerically with SciPy使用 SciPy 以数值方式求解 ODE
【发布时间】：2019-04-13 12:03:22
【问题描述】：

我正在尝试找到一个数值解并最终绘制 Gyllenberg-Webb 模型（癌细胞生长模型）。这个模型看起来像：

其中β是增殖细胞的繁殖率，µp是增殖细胞的死亡率，µq是静止细胞的死亡率，r0和ri是函数（转换率) 的N(t)。还有N(t) = P(t)+Q(t)。

出于我的目的，我在这里定义了 r_0(N) = bN 和 r_i(N) = aN 以使事情更简单。

我的问题是当我尝试用 pyplot 绘制我的解决方案时，我得到了

ValueError: x 和 y 必须具有相同的第一维

我想这是不言自明的，但我不确定如何在不破坏其他所有内容的情况下修复它。

我的代码是：

import numpy as np
import matplotlib.pyplot as plt
import scipy.integrate

def fun(P,t, params):
    beta, mp,b,N,Q = params
    return(beta-mp-(b*N))*P+(a*N)*Q

params = (0.5,0.6,0.7,0.8,0.9)

tvec = np.arange(0,6,0.1)
s1 = scipy.integrate.odeint(
    fun,
    y0 = 1,
    t = tvec,
    args = (params,))

#print(s1)

plt.plot(fun,tvec)

【问题讨论】：

标签： python matplotlib scipy differential-equations

【解决方案1】：

您当前正在绘制fun 与tvec。你真正想要的是绘制tvec vs s1。您还必须在fun 中定义参数a；我在下面的代码中将其设置为1：

import numpy as np
import matplotlib.pyplot as plt
import scipy.integrate


def fun(P, t, params):
    beta, mp, b, N, Q = params
    return (beta-mp-(b*N))*P + (1.0 * N)*Q

params = (0.5, 0.6, 0.7, 0.8, 0.9)

tvec = np.arange(0, 6, 0.1)
s1 = scipy.integrate.odeint(
    fun,
    y0=1.,
    t=tvec,
    args=(params,))

plt.plot(tvec, s1)
plt.show()

这将绘制：

【讨论】：

【解决方案2】：

最后你会想要解决耦合系统。这并不复杂，只需将状态对象做成一个向量，并以正确的顺序返回导数即可。

def fun(state,t, params):
    P, Q = state
    beta, mp, mq, a, b = params
    N = P+Q
    r0N, riN = b*N, a*N
    return [ (beta-mp-r0N)*P + riN*Q, r0N*P - (riN+mq)*Q ]

params = (0.5,0.6,0.7,0.8,0.9)

tsol = np.arange(0,6,0.1)
sol = odeint( fun, y0 = [ 1, 0],  t = tsol,  args = (params,))

Psol, Qsol = sol.T; plt.plot(tsol, Psol, tsol, Qsol)

【讨论】：

这真的很有帮助，谢谢你好心的陌生人