228
8 Solving Ordinary Differential Equations
Since there is no loss in the R category (people are either recovered and immune,
or dead), we are done with the modeling of this category. In fact, we do not strictly
need Eq. (8.14) for R, but extensions of the model later will need an equation for R.
Dividing by Δt in (8.13) and (8.14) and letting Δt → 0, results in the
corresponding differential equations
I
= βSI − γ I,
(8.15)
and
R
= γ I .
(8.16)
To summarize, we have derived three difference equations and three differential
equations, which we list here for easy reference. The difference equations are:
S
n+1
= S
n
− βΔtS
n I
n ,
(8.17)
I
n+1
= I
n
+ βΔtS
n I
n
− γ ΔtI
n ,
(8.18)
R
n+1
= R
n
+ γ ΔtI
n .
(8.19)
Note that we have isolated the new unknown quantities S n+1 , I n+1 , and R n+1 on
the left-hand side, such that these can readily be computed if S n , I n , and R n are
known. To get such a procedure started, we need to know S 0 , I 0 , R 0 . Obviously, we
also need to have values for the parameters β and γ .
The three differential equations are:
S
= −βSI,
(8.20)
I
= βSI − γ I,
(8.21)
R
= γ I .
(8.22)
This differential equation model (and also its discrete counterpart above) is known
as an SIR model. The input data to the differential equation model consist of the
parameter values for β and γ , as well as the initial conditions S(0) = S 0 , I (0) = I 0 ,
and R(0) = R 0 .
8.3.2 A FE Method for the System of ODEs
Let us apply the same principles as we did in Sect. 8.2.2 to discretize the differential
equation system by the Forward Euler method. We already have a time mesh and
time-discrete quantities S n , I n , R n , n = 0, . . . , N t . The three differential equations
are assumed to be valid at the mesh points. At the point t n we then have
S
(t n ) = −βS(t n )I (t n ),
(8.23)
I
(t n ) = βS(t n )I (t n ) − γ I (t n ),
(8.24)
R
(t n ) = γ I (t n ),
(8.25)