UP
Scilab/Xcosでシミュレーション
【練習】ニュートンの運動方程式
【COVID19】SIRモデル
【練習】ニュートンの運動方程式
調和振動子の運動方程式：
\[
\begin{array}{rl}
\frac{dx}{dt}&=v\\
m\frac{dv}{dt}&=-m\omega^2x\\
\end{array}
\]
積分形に直すと
\[
\begin{array}{ll}
x(t)&=x_0+\int_0v(t)dt\\
v(t)&=v_0-\int_0\omega^2x(t)dt\\
\end{array}
\]
となるのでXcos上で
などと組めばよい。ただし\(x_0=3.0,\ v_0=1.5,\ \omega^2=0.4\)としている。
実行結果：
EM.zcos
【COVID19】SIRモデル
\(S:\ \)未感染者の割合
\(I:\ \)感染者の割合（その時点で感染している人）
\(R:\ \)回復者の割合
普通これらは割合ではなく絶対数に取ることが多いようだがここでは
割合とする。
\(S\rightarrow 感染\rightarrow I\rightarrow 回復\rightarrow R\)
\(S\)と\(I\)が出会うと感染がおこる。
\[
\begin{array}{rl}
\frac{dS}{dt}&=-\beta IS\\
\frac{dI}{dt}&=(\beta S-\gamma)I\\
\frac{dR}{dt}&=\gamma I\\
\end{array}
\]
\(\beta:\ \)感染速度
\(\gamma:\ \)回復速度
\(\beta S>\gamma\)のときに\(I\)が増大する（流行）。
積分型に書き直す：
\[
\begin{array}{rl}
S&=S_0-\int_0\beta IS\,dt\\
I&=I_0+\int_0(\beta S-\gamma)I\,dt\\
R&=R_0+\int_0\gamma I\,dt\\
\end{array}
\]
\(\beta=\epsilon\,r_0/\tau_0\)
\(\gamma=1/\tau_0\)
\(r_0:\ \)（\(S=1\)のときに）一人の感染者が何人に感染させるか
\(\tau_0:\ \)感染期間
外出制限等で接触機会を変化させるパラメータとして接触係数\(\epsilon\)を導入した。
\(\epsilon=1\)で通常通り、\(\epsilon=0.5\)で接触機会を半減。
パラメータ、初期条件等：\(\tau_0=7,\ r_0=2,\ S_0=0.999,\ I_0=0.001,\ R=0.0\)
SIRm.zcos
\(S:\,黒,\ I:\,緑,\ R:\,赤,\ \epsilon:\,青\)
感染制御なし
感染制御あり
一旦接触機会を減らしてゆっくり平常に戻すと収束後の未感染者数\(S\)が多くなる
（累計感染者\(R\)が少なくて済む）。
非常にゆっくり戻すと\(S\)は\(1/r_0\)に漸近するかというと、傾きを減らして
いったとき戻す過程で流行が複数回起こったりして必ずしもよい結果となるとは
限らない（極限自体は\(1/r_0\)である可能性もあるが）。
流行が2つに分かれる例：\(\tau_0=3.3,\ r_0=3,\ S_0=0.99,\ I_0=0.01,\ R=0.0\)
2020.04
2020/04/13 21:18:50
更新