简单的ODE函数，但值爆炸

Question

伙计们，我正在尝试运行一些简单的SIR模拟，其想法是让人们适应他们的环境，这样，如果没有人感染，他们会与其他7个人互动，但是如果被感染的人群很高，他们就会减少互动。这个特定的变量是b，我用dbdt表示其演变，但是由于某种原因，当我对其进行绘制时，它是单调的。当我生成dbdt时，我得到正确的值7,7,7,6,6,5,4,4,7,7，这意味着随着流行病的传播，它们会减少接触，但随着它们增加得更好，但是由于某种原因，该函数对其进行计算并得到荒谬的数字，例如140。

import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt
from math import log, exp
def g(s,t):
    g = 0.2
    k = 1
    N = 100 
    theta = 0.1
    delta = 0.08
    S = s[0]
    I = s[1]
    R = s[2]
    D = s[3]
    C = s[4]
    b = s[5]
    dbdt = b* exp(-k*I/N)
    if dbdt >= 7:
        dbdt = 7
    dSdt = - b*S*I/N 
    dIdt = b*S*I/N - g*I
    dRdt = g*I - theta*R
    dDdt = delta*theta*R 
    dCdt = (1-delta)*theta*R
    return [dSdt, dIdt,dRdt,dDdt,dCdt,dbdt]

#Graph x axis
t = np.linspace(0,20)
#Vector for initial conditions of each equation
s0=[99,1,0,0,0,7]
#Ordinary differential equation
s = odeint(g,s0,t)

plt.subplot(2, 1, 1)
plt.plot(t,s[:,0],'r--', linewidth=2.0)
plt.plot(t,s[:,1],'b-', linewidth=2.0)
plt.plot(t,s[:,2],'g-', linewidth=2.0)
plt.plot(t,s[:,3],'y-', linewidth=2.0)
plt.plot(t,s[:,4],'o-', linewidth=2.0)
plt.xlabel("t")
plt.ylabel("Number")
plt.legend(["Susceptible","Infected","Resolving","Dead","Recovered"])

plt.subplot(2, 1, 2)
plt.plot(t,s[:,5],'m-', linewidth=2.0)
plt.legend(["Distancing"])


plt.show()
所以问题是：我必须修复什么才能使其正常工作？请注意，如果我保持b为外源代码，则该代码有效。 

伙计们，我正在尝试运行一些简单的SIR模拟，其目的是让人们适应他们的环境，这样，如果没有人感染，他们会与其他7个人互动，但是如果被感染的人群是......]

Answer 1

import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt
from math import log, exp

beta = 7
def g(s,t):
    g = 0.2
    k = 1
    N = 100
    theta = 0.1
    delta = 0.08
    S = s[0]
    I = s[1]
    R = s[2]
    D = s[3]
    C = s[4]
    #Current value of beta
    b = s[5]
    #New Value of beta
    newbeta = exp(log(beta)-k*I/N)
    #New difference
    dbdt = newbeta - b
    dSdt = - b*S*I/N 
    dIdt = b*S*I/N - g*I
    dRdt = g*I - theta*R
    dDdt = delta*theta*R 
    dCdt = (1-delta)*theta*R
    return [dSdt, dIdt,dRdt,dDdt,dCdt,dbdt]

#Graph x axis
t = np.linspace(0,20)
#Vector for initial conditions of each equation
s0=[99,1,0,0,0,beta]
#Ordinary differential equation
s = odeint(g,s0,t)


f = plt.figure(figsize=(20,5))

plt.subplot(1, 2, 1)
plt.plot(t,s[:,0],'r--', linewidth=2.0)
plt.plot(t,s[:,1],'b-', linewidth=2.0)
plt.plot(t,s[:,2],'g-', linewidth=2.0)
plt.plot(t,s[:,3],'y-', linewidth=2.0)
plt.plot(t,s[:,4],'o-', linewidth=2.0)
plt.xlabel("t")
plt.ylabel("Number")
plt.legend(["Susceptible","Infected","Resolving","Dead","Recovered"])

plt.subplot(1, 2, 2)
plt.plot(t,s[:,5],'m-', linewidth=2.0)
plt.legend(["Distancing"])

plt.savefig('simple.png', dpi=1000)
plt.show()