Mathematical Analysis and Numerical Solution of a Boundary Value
Problem for the Covid-19 SIR Model
SERDAR SALDIROĞLU1, SERDAL PAMUK2
1Department of Basic Sciences,
Naval Academy National Defense University,
Istanbul,
TURKEY
2Department of Mathematics, Faculty of Science and Arts,
Kocaeli University, Kocaeli,
TURKEY
Abstract: - This paper extends the work presented at IX. International Istanbul Scientific Research Congress
held on May, 14-15, 2022, Istanbul/Türkiye. In that Congress the Authors narrowly focused on the numerical
solutions of a boundary value problem for the Covid-19 SIR model appearing in literature. In this study this
boundary value problem is solved numerically for all cases and also the stability analysis of the equilibrium
points of the model is presented. The basic reproduction number is obtained and the importance of this
number for the stability and the instability of the equilibrium points is emphasized. Numerical solutions are
obtained using bvp4c, a boundary value problem solver in MATLAB and the results are presented in figures.
Key-Words: - Covid-19, Stability Analysis, Mathematical Modelling, MATLAB, Boundary Value Problem,
Mathematical Analysis.
Received: April 9, 2023. Revised: December 17, 2023. Accepted: February 12, 2024. Published: April 3, 2024.
1 Introduction
The Covid-19 disease, which first appeared in China
in December 2019 and spread all over the world in a
short time, has taken its place among the epidemics,
which has significantly affected life globally.
Explaining the effects of this epidemic with
mathematical models, as in other epidemics that
have radically changed the life of humanity, has an
important place in the literature.
Understanding past outbreaks can help us better
prepare for future ones. Communicable diseases
such as plague, malaria, smallpox, cholera, measles,
tuberculosis, AIDS and flu, which are transmitted
from animals, soil, water or human to human, have
affected social life throughout history, causing
demographic, social and economic problems. The
answer to the question of how we will deal with
future outbreaks can be obtained by examining past
outbreaks, [1].
Mathematical models are very important in
analyzing the spread and control of infectious
diseases, [2], [3], [4], [5], [6], [7], [8], [9], [10].
Various mathematical models were used to be able
to comment on these diseases and examine
infectious diseases. In the SIR model, the society is
divided into three groups, [11].
The well-known compartment model, consisting
of susceptible, infected and recovered
compartments, abbreviated as the SIR model, has
been commonly used in infectious disease spread
simulations for more than half a century although
the mathematical model is very simple, [12].
Susceptible individuals have never been
infected and therefore can catch the disease.
Individuals who are infected can spread the disease,
also individuals in the recovered state are assumed
to be immune for life, [13].
Epidemiologically SIR model is used to
determine the causes of diseases and health
problems, to eliminate the natural development of
diseases, to determine the health levels of the
population, to investigate the change in time and
when compared with other societies, to evaluate the
results of clinical research, to evaluate the
effectiveness of health services, and to determine
the risks of encountering certain health problems of
individuals in the population, [14].
PROOF
DOI: 10.37394/232020.2024.4.2
Serdar Saldiroğlu, Serdal Pamuk
E-ISSN: 2732-9941
11
Volume 4, 2024
2 The Model
The SIR model, [10], uses the system of ordinary
differential equations
(1)
(2)
(3)
where S, I and R mean the susceptible,
infectious, and removed populations, respectively,
and the parameter is the transmission rate and is
the recovery rate (or in other words, the duration of
infection ). In the SIR model, for
example, a person could change his or her condition
from susceptible to infected with a ratio , then to
removed with a ratio . Removed persons will never
become susceptible. From eqs. (1)-(3) we have:
(4)
which means that .
This is the total population size, and we denote it by
N. At the time t=0 we assume ,
and . From eqs. (1) and (2) we
get:
,
(5)
which yields that:
(6)
Integrating both sides we get:
,
(7)
where C is an arbitrary constant. Using initial
conditions we have:
(8)
If eq. (8) is substituted in eq. (7) we get:
.
(9)
Similarly, from eqs.(1) and (3) we obtain:
,
(10)
whose solution is:
.
(11)
3 Stability Analysis of the Equilibrium
Points of the Model
In this section [15], we find the equilibrium points
of the model by considering that the model has two
different set of equilibrium points, namely the
disease-free equilibrium points and the disease-
present equilibrium points.
Since R represents the number of the removed
populations, it is enough to consider only the eqs.
(1) and (2) to find the disease-free equilibrium
points of the model. Therefore, it is easy to see that
the disease-free equilibrium points of the system
(1)-(2) are
for any real number .
We must consider the whole model consisting
of eqs. (1)-(3) to find the disease-present
equilibrium points which we call it
.
If eq. (2) is set equal to zero we find
since . From eq. (9) one gets:
(12)
Therefore, we obtain:
(13)
(14)
Let
be the basic reproduction
number (replicate number) that measures the
average number of the new infected individuals
generated by a single infected individual in a
population of susceptible individuals. The value of
will indicate whether the epidemic could occur
or not.
As a result, the disease-present equilibrium
points become
,
,
.
(15)
Since
must be one obtains
.
PROOF
DOI: 10.37394/232020.2024.4.2
Serdar Saldiroğlu, Serdal Pamuk
E-ISSN: 2732-9941
12
Volume 4, 2024
3.1 Local Stability Analysis of the Equilibria
Now we proceed to study the stability behavior of
equilibria
and
.
3.1.1 Local Stability Analysis of the Disease-Free
Equilibrium
In this section, we analyze the local stability of the
COVID-19 disease-free equilibrium. Computing the
Jacobian matrix:
,
(16)
at the disease-free equilibrium points
, we get:
(17)
The characteristic equation of this matrix becomes:
,
(18)
whose solutions are and .
These are the eigenvalues of J. If we get
and , so that we have an attractive
equilibrium line. Therefore, the equilibrium points
are stable but not asymptotically.
If we obtain and , so
that we have a repulsive line of the equilibria, which
means that the equilibrium points are unstable, [16].
3.1.2 Local Stability Analysis of the Disease-
Present Equilibrium
In this section, we analyze the local stability of the
disease-present equilibrium. If we now use the
disease-present equilibrium point in eq. (16) we get:
(19)
In this case the characteristic equation becomes:
(20)
According to the Routh Hurwitz stability
criterion, if the equilibrium point is
asymptotically stable and it is unstable if
[16], [17], [18], [19], [20].
4 Solutions to the Boundary Value
Problem for the SIR Model
In the following computations we take
,,
, , ,
10, ,
.
4.1 Susceptible Individuals
If we differentiate eq. (1) we get:
(21)
Using eq. (2) one obtains:
(22)
From eq. (1) it is easy to see that:
(23)
We now solve this ordinary differential equation
(ODE) with the boundary conditions provided above
to study the number of healthy individuals who may
contract the disease in 30 days.
Fig. 1: Susceptible Individuals and the Change of
Them in Time.
As seen in Figure 1, the dashed curve is the
graph of susceptible individuals while the solid
curve is the graph of the derivative of susceptible
individuals over 30 days. While the number of the
population that can get sick is small at the
beginning, it increases over time and progresses to
the entire population. The solid curve shows that the
population that may be sick continues to increase
over time and progress toward the entire population.
4.2 Infected Individuals
Similarly, if we differentiate eq. (2) we get:
(24)
Using eq. (1) one obtains:
(25)
PROOF
DOI: 10.37394/232020.2024.4.2
Serdar Saldiroğlu, Serdal Pamuk
E-ISSN: 2732-9941
13
Volume 4, 2024
and
(26)
We now solve this ODE with the boundary
conditions given above to study the number of
infected individuals in 30 days.
Fig. 2: Infected Individuals and the Change of Them
in Time.
The dashed curve is the graph of infected
individuals while the solid curve is the graph of the
derivative of infected individuals over 30 days.
Although the number of infected populations is
small at the beginning, it increases over time. Also,
the solid curve appears as evidence that the
population that may become ill continues to increase
over time. Figure 2 also shows us that the change in
the number of infected individuals in society is
proportional to the change in the number of infected
individuals over time.
4.3 Removed Individuals
Recalling and differentiating
eq. (3) we get:
(27)
From eq. (1) one obtains:
(28)
After some tedious work we come up with the
differential equation
(29)
We now obtain the solution of this ODE with
the boundary conditions provided above to study the
number of removed individuals in 30 days.
Fig. 3: Removed Individuals and The Change of
Them in Time
The dashed curve is the graph of deceased
individuals while the solid curve is the graph of the
derivative of deceased individuals over 30 days.
Although the number of individuals who died due to
the epidemic is very low at the beginning, it
continues to increase over time. Both curves
continue linearly as they represent the number of
dead individuals over time. According to Figure 3, a
total of 2500 individuals died in the first 30 days.
Therefore, an average of 83-84 people died daily.
Indeed, Figure 3 reveals to us the whole reality
regarding the number of changes in society. The
daily number of dead individuals clearly expresses
the number of changes in total deaths over time.
Therefore, Figure 3 gives us descriptive information
about all possible situations.
5 Discussion
In this paper we have shown numerically that as the
virus infects individuals, the number of susceptible
individuals decrease, while the number of
individuals exposed to the virus increase. In this
model, transitions from one group to another are
unidirectional and there are no returns. In other
words, an individual who has recovered will not be
re-infected with the disease. Those who catch the
virus interact with people who cannot be isolated
and who do not have protection, and they transmit
the disease to them, and the epidemic continues in
this way. As the infected people regain their health
or die, the number of individuals in the non-
infectious group naturally increases.
PROOF
DOI: 10.37394/232020.2024.4.2
Serdar Saldiroğlu, Serdal Pamuk
E-ISSN: 2732-9941
14
Volume 4, 2024
In summary, among epidemiological models;
the susceptible-infected-recovered (SIR) model,
which says that infection provides permanent
immunity, they have been used to describe diseases
that spread where there is healing.
However, the SIR model does not have a latent
stage (no exposed individual) and in this case is not
suitable as a model for infectious exposure-
progressed diseases such as COVID-19, [2].
6 Conclusion
In this study we have solved numerically a boundary
value problem related to the Covid-19 SIR model
over a period of time. By doing this we have
observed the changes in the number of susceptible,
infected and removed individuals and shown these
changes in figures. Also, we have provided the
stability analysis of the equilibrium points of the
model and emphasized the importance of basic
reproduction number (replicate number) for the
stability and for the instability of the equilibrium
points. We have seen that if , there is an
increase in the epidemic growth rate, if
there is a decrease in the epidemic growth rate and if
the epidemic growth rate is traveling at a
constant speed.
References:
[1] Parildar, Hülya and Dikici, Mustafa Fevzi,
History of pandemics, Clinical Medicine
Family Medicine, 2020, 12, pp.1-8.
[2] Carcione, José M.; Santos, Juan E.; BA, Jing,
A simulation of a COVID-19 epidemic based
on a deterministic SEIR model, Frontiers in
Public Health, 8, 2020, pp. 230.
DOI:10.3389/fpubh.2020.00230
[3] Yang, Hyun Mo; Lombardi Junior, Luis
Pedro; YANG, Ariana Campos, Are the SIR
and SEIR models suitable to estimate the
basic reproduction number for the Covid-19
epidemic?, Medrxiv, 2020,
DOI:10.1101/2020.10.11.20210831
[4] Ahmed, Idris, Modu, Goni Umar, Yusuf,
Abdullahi, Kumam, Poom, Yusuf, Ibrahim, A
mathematical model of coronavirus disease
(COVID-19) containing asymptomatic and
symptomatic classes, Results in Physics, 21,
2021, DOI: 10.1016/j.rinp.2020.103776
[5] Derman, Orhan, Humanity's struggle with
epidemic diseases throughout history, Journal
of Child Health and Diseases, 63(1-4), 2020,
pp. 26.
[6] Javeed, Shumaila, Anjum, Subtain, Alimgeer,
Khurram Saleem, Atif, M., Khan, Mansoor
Shaukat, Farooq, W. Aslam, Hanif, Atıf,
Ahmad, Hijaz, YAO, Shao-Wen. A novel
mathematical model for COVID-19 with
remedial strategies, Results in Physics, 27,
2021. DOI: 10.1016/ j.rinp. 2021.104248
[7] Kadah, Nezir, Stability and passivity analysis
of nonlinear RLC circuits, Master's Thesis,
Van Yuzuncu Yıl University, Institute of
Science and Technology, 2019, 10257087.
[8] Labzai, Abderrahim, Kouidere, Abdelfatah,
Balatif, Omar, Rachik, Mostafa, Stability
analysis of mathematical model new corona
virus (COVID-19) disease spread in
population, Commun. Math. Biol. Neurosci,
41, 2020. DOI: 10.28919/cmbn/4551
[9] Çay, İrem, On the local and global stability of
an sirs epidemic model with logistic growth
and information intervention, Turkish Journal
of Mathematics, 45(4), 2021, pp. 1668-1677.
[10] Pamuk, Serdal, Saldiroğlu, Serdar, Boundary
Value Problems and Numerical Solutions for
the COVID-19 Secret Model, 9th
International Istanbul Scientific Research
Congress, Istanbul, 14-15 May 2022.
[11] Güzey, Nurbanu, Parameter estimation with
the RSS method in SIR-based modeling of the
COVID-19 spread in Turkey, Journal of
Gumushane University Institute of Science
and Technology, 11(3), 2021, pp. 956-963.
DOI: 10.17714/ gumusfenbil.757291
[12] Hirose, Hideo, A relationship between the SIR
model and the generalized logistic distribution
with applications to SARS and COVID-19,
10th International Congress on Advanced
Applied Informatics (IIAI-AAI). IEEE, 2021,
pp. 837-842. arXiv Preprint
arXiv:2009.09653.
[13] Zaman, Gul; Kang, Yong Han; Jung, Il Hyo,
Stability analysis and optimal vaccination of
an SIR epidemic model, Bio Systems, 93(3),
2008, pp. 240-249.
[14] Işik, Nurettin and Kaya, Abdullah,
Mathematical models and herd immunization
in the spread and control of infectious
diseases, Ataturk University Journal of
Veterinary Sciences, 15 (3), 2020, pp.301-
307. DOI: 10.17094/ataunivbd.715371
[15] Faniran, Taye, Bakare, Emmanuel Afolabi,
Potucek, Radovan, Ayoola, Olusola Ezekiel,
Global and Sensitivity Analyses of
Unconcerned COVID-19 Cases in Nigeria: A
Mathematical Modeling Approach, WSEAS
Transactions on Mathematics, 2021, 20, pp.
PROOF
DOI: 10.37394/232020.2024.4.2
Serdar Saldiroğlu, Serdal Pamuk
E-ISSN: 2732-9941
15
Volume 4, 2024
218-234,
https://doi.org/10.37394/23206.2021.20.23.
[16] Assimakis, Nicholas, Adam, Maria, Closed
form solutions of lyapunov equations using
the vech and veck Operators, WSEAS
Transactions on Mathematics, 2021, 20, pp.
276-282,
https://doi.org/10.37394/23206.2021.20.28.
[17] Zhu, Qun, Lin, Shijia, Wu, Runxin, Chen,
Fengde, Dynamic Behaviors of a
Commensalism Model Incorporating
Nonselective Harvesting in a Partial Closure,
WSEAS Transactions on Mathematics, 2023,
22, pp. 798-806,
https://doi.org/10.37394/23206.2023.22.88.
[18] Chen, Fengde, Lin, Sijia, Chen, Shangming,
CHON, Yanbo, Global stability of leslie-
gower predator-prey model with density
dependent birth rate on prey species and prey
refuge, WSEAS Transactions on Systems and
Control, 2023, 22, pp. 41-48,
https://doi.org/10.37394/23202.2023.22.5.
[19] Krokavec, Dušan, Ostensible metzler linear
uncertain systems: goals, LMI synthesis,
constraints and quadratic stability, WSEAS
Transactions on Systems and Control, 2023,
18, pp. 255-262,
https://doi.org/10.37394/23203.2023.18.25.
[20] Sánchez, Fabián Toledo, Alzate, Pedro Pablo
Cardenas, Salcedo, Carlos Arturo Escudero, A
note on the stability of a modified Lotka-
Volterra model using Hurwitz polynomials,
WSEAS Transactions on Mathematics, 2021,
20, pp. 431-441,
https://doi.org/10.37394/23206.2021.20.44.
APPENDIX
MATLAB CODES
MATLAB Code to obtain the number of Susceptible
Individuals
N = 60480000; % total population
r0 = 10; % breeding value
R0 = 0; % first population to recover
i_period = 5; % infectious period duration
beta = 1/i_period ;% recycling rate
alpha = r0*beta/N; % infection rate
ya(1) = 6048000;
yb(1) = 60446500;
f=@(x,y)[y(2);y(2)^2/y(1)+(alpha*y(1)-beta)*y(2)];
bc=@(ya,yb) [ya(1)-6048000;yb(1)-60446500];
%boundary conditions
xmesh= linspace (0,30,100); %create a network
solinit = bvpinit(xmesh, [1 0]); % first guess of the
solution
sol = bvp4c(f, bc, solinit); % run solver
figure;
plot(sol.x, sol.y(1,:), '--')
hold on;
plot(sol.x, sol.y(2,:), '-')
xlabel('Day');
ylabel('Susceptible Individuals S');
legend('S', 'dS/dt');
hold off;
MATLAB Code to obtain the number of Infected
Individuals
N = 60480000; % total population
r0 = 10; % breeding value
R0 = 0; % first population to recover
i_period = 5; % infectious period duration
beta = 1/i_period ;% recycling rate
alpha = r0*beta/N; % infection rate
ya(1) = 1;
yb(1) = 31000;
f=@(x,y)[y(2);-alpha*y(2)*y(1)+y(2)^2/y(1)-
alpha*beta*y(1)^2];
bc=@(ya,yb) [ya(1)-1;yb(1)-31000]; % boundary
conditions
xmesh= linspace (0,30,200); % create a network
solinit = bvpinit(xmesh, [100 0]); % first guess of
the solution
sol = bvp4c(f, bc, solinit); % run solver
figure;
plot(sol.x, sol.y(1,:), '--')
hold on;
plot(sol.x, sol.y(2,:), '-')
xlabel('Day');
ylabel('Infected Individuals');
PROOF
DOI: 10.37394/232020.2024.4.2
Serdar Saldiroğlu, Serdal Pamuk
E-ISSN: 2732-9941
16
Volume 4, 2024
legend('I', 'dI/dt');
hold off;
MATLAB Code to obtain the number of Removed
Individuals
N = 60480000; % total population
H0 = 6048000;
r0 = 10; % breeding value
R0 = 0; % first population to recover
i_period = 5; % infectious period duration
beta = 1/i_period ;% recycling rate
alpha = r0*beta/N; % infection rate
ya(1) = 0;
yb(1) = 2500;
f=@(x,y)[y(2);(alpha*H0*exp(-(alpha/beta)*y(1))-
beta)*y(2)];
bc=@(ya,yb) [ya(1);yb(1)-2500]; %boundary
conditions
xmesh = linspace (0,30,31); % create a network
solinit = bvpinit(xmesh, [1 0]); %% first guess of
the solution
sol = bvp4c(f, bc, solinit); % run solver
figure;
plot(sol.x, sol.y(1,:), '--')
hold on;
plot(sol.x, sol.y(2,:), '-')
xlabel('Day');
ylabel('Removed Individuals');
legend('R', 'dR/dt'); hold off;
Contribution of Individual Authors to the
Creation of a Scientific Article (Ghostwriting
Policy)
The authors equally contributed in the present
research, at all stages from the formulation of the
problem to the final findings and solution.
Sources of Funding for Research Presented in a
Scientific Article or Scientific Article Itself
No funding was received for conducting this study.
Conflict of Interest
The authors have no conflicts of interest to declare.
Creative Commons Attribution License 4.0
(Attribution 4.0 International, CC BY 4.0)
This article is published under the terms of the
Creative Commons Attribution License 4.0
https://creativecommons.org/licenses/by/4.0/deed.en
_US
PROOF
DOI: 10.37394/232020.2024.4.2
Serdar Saldiroğlu, Serdal Pamuk
E-ISSN: 2732-9941
17
Volume 4, 2024
