# Mathematical Model and Analysis on the Impacts of Vaccination and Treatment in the Control of the COVID-19 Pandemic with Optimal Control

Authors: Alemzewde Ayalew; Yezbalem Molla; Tenaw Tilahun; Tadele Tesfa
Journal: Journal of Applied Mathematics
 (2023)

---

## Abstract

In this study, a nonlinear deterministic mathematical model that evaluates two important therapeutic measures of the COVID-19 pandemic: vaccination of susceptible and treatment for infected people who are in quarantine, is formulated and rigorously analyzed. Some of the fundamental properties of the model system including existence and uniqueness, positivity, and invariant region of solutions are proved under a certain meaningful set. The model exhibits two equilibrium points: disease-free and endemic equilibrium points under certain conditions. The basic reproduction number,R0, is derived via the next-generation matrix approach, and the dynamical behavior of the model is explored in detail. The analytical analysis reveals that the disease-free equilibrium solution is locally as well as globally asymptotically stable when the associated basic reproduction number is less than unity which indicates that COVID-19 dies out in the population. Also, the endemic equilibrium point is globally asymptotically stable whenever the associated basic reproduction number exceeds a unity which implies that COVID-19 establishes itself in the population. The sensitivity analysis of the basic reproduction number is computed to identify the most dominant parameters for the spreading out as well as control of infection and should be targeted by intervention strategies. Furthermore, we extended the considered model to optimal control problem system by introducing two time-dependent variables that represent the educational campaign to susceptibles and continuous treatment for quarantined individuals. Finally, some numerical results are illustrated to supplement the analytical results of the model using MATLAB ode45.


---

## Body

## 1. Introduction

Nowadays, the world is faced with a number of issues, including climatic change, racism, terrorism, unemployment, poverty, and public health crises. Among them, the global coronavirus disease 2019 (COVID-19) outbreak is still a major public health concern. Severe acute respiratory syndrome-2 (SARS-CoV-2) is the virus that causes COVID-19 [1–4]. The first case of this new coronavirus, which frequently manifests as pneumonia, was found in Wuhan, China, in December 2019 [3, 5–7]. One of the things which make this disease so dangerous is that it spreads very quickly between people. Direct contact with contaminated surfaces and inhalation of respiratory droplets from infected individuals are the two ways that COVID-19 spreads from person to person [3, 6, 8, 9]. COVID-19 typically takes 2 to 14 days to incubate after being exposed to the virus [3, 10–12]. Fever, a dry cough, shortness of breath, and fatigue are some of the most prevalent clinical signs of COVID-19 [5, 13].The outbreak was classified as a Public Health Emergency of International Concern by the World Health Organization (WHO) in January 2020 and as a pandemic on March 11, 2020 [14, 15]. More than 2.58 million fatalities and over 116.17 million confirmed cases were recorded as of March 7th, 2021 [16]. Globally, there have been more than 197 million new cases, and of these, more than 4 million people have died until the end of July 2021 [17, 18]. The WHO report of COVID-19 indicated that on 10 August 2021, more than 202.14 million people were infected with the virus, and over 4.28 million have died.The development of applicable intervention strategies to mitigate the COVID-19 outbreak requires a deep understanding of its transmission process and prevention programs. Mathematical modeling in epidemiology has become an effective tool for understanding and describing the dynamics of infectious disease [19]. The basic and important research subjects for epidemiological models are the existence of the threshold value which distinguishes whether the infectious disease will die out or persist via the local and global stability of the disease-free and endemic equilibrium [20].Since the COVID-19 outbreak, different scholars have studied the transmission dynamics of the disease by considering different scenarios using mathematical modeling [1–3, 7, 10, 18, 21–40]. To the best of our knowledge, mathematical models that take the possibility of using additional preventative measures in addition to vaccination to stop the spread of COVID-19 are extremely uncommon. As a result, in this study, we are motivated to develop a model for COVID-19 transmission using model [41] as our framework by taking vaccination to susceptible human class and Holling type II with optimal control into consideration. In order to achieve this, we present SVIDR and this proposed model differs from others in that those susceptible people who received vaccinations joined the class of vaccinated humans, whereas other susceptible people who did not receive vaccinations might contract an infection from an infectious agent with the term βSI/σ+ηI and advance into the infectious class. The rest of this study is organized as follows. We formulated the transmission dynamics of the COVID-19 model in Section 2. The formulated model is analyzed in Section 3. The optimal control problem with some basic properties are discussed in Section 4. We illustrated numerical simulations of the model in Section 5. Finally, conclusion and conclusion remarks along with study suggestions are provided in Section 6.


## 2. Mathematical Model Formulation

In this section, we divided the entire population into five compartments based on the disease status of each individual: susceptible individuals (S), vaccinated individuals (V), individuals infected with COVID-19 (I), quarantined individuals (D), and recovered individuals (R). In the formulation of our model, the following assumptions are taken into account:
(i)
When susceptible individuals receive the COVID-19 vaccination, they will join a class of humans who have received the vaccine at a rate ofϕ(ii)
The class of susceptible people only contract COVID-19 infection when they come into contact with an infected person at a rate ofβ(iii)
The class of vaccine recipients can contract COVID-19 infection only when they come into contact with those who have the disease at a rate ofψ(iv)
The class of COVID-19-infected individuals can join either the quarantine class with the rate ofα or recovered class with the rate of γ(v)
In addition to vaccination, susceptible people take other precautions against COVID-19(vi)
All parameters involved in the model are non-negativeIndividuals were recruited into the susceptible class with the recruitment rateΛ. We assumed that all compartments would experience the same amount of natural death μ. Susceptible individuals decrease due to COVID-19 infection following effective contact with infected individuals at the rate β. The infectious and quarantined human classes decrease due to COVID-19-induced death with the rate of ϵ and π, respectively. With this in mind, the total number of human population at a given time t is given by Nt=St+Vt+It+Dt+Rt. With regard to the above considerations, the compartmental developed COVID-19 model schematic flowchart is shown in Figure 1 and model parameters are described in Table 1.Figure 1
The schematic flowchart of COVID-19 dynamics.Table 1
Description of model parameters.
NotationBiological meaningΛThe recruitment rate of susceptible individualsϕThe progression rate of susceptibles into vaccinated individual classμThe natural death rate of all individualsβThe contact rate of susceptibles with infectious individualsψThe contact rate of vaccinated individuals with infected individualsωRecovery rate of quarantined individuals into recovered individualsϵDisease-induced death rate of infectious individualsπDisease-induced death rate of quarantined individualsγRecovery rate of infected individuals into recovered individuals’ classσHalf saturation constantηPositive constantThe above model assumptions and Figure1 lead to the following autonomous system of nonlinear differential equations:
(1)Ṡt=Λ−βSIσ+ηI−ϕS−μSV̇t=ϕS−μ+ψIVİt=βSIσ+ηI+ψVI−μ+ϵ+γI−αIḊt=αI−μ+π+ωDṘt=γI+ωD−μR,with the initial conditions:
(2)S0=S0>0V0=V0≥0I0=I0≥0D0=D0≥0R0=R0≥0.


## 3. Mathematical Analysis of the Model

In this section, the fundamental properties of the model system (1) are investigated, including the existence and uniqueness of solutions, nonnegativity of solutions, invariant region, basic reproduction number, equilibria, and its stability analysis.Lemma 1 (existence and uniqueness of model solutions).
The model system (1) with initial condition (2) admits the unique solution in Cℝ+,ℝ+5,∀t≥0.Proof.
The model system (1) can be written in the form of ẋ=fx, where
(3)ẋ=ṠtV̇tİtḊtṘt,fx=Λ−βSIσ+ηI−ϕS−μSϕS−μ+ψIVβSIσ+ηI+ψVI−μ+ϵ+γI−αIαI−μ+π+ωDγI+ωD−μR.
All the right-hand side components of the functionfx in the system (3) are continuously differentiable almost everywhere in Cℝ+,ℝ+5, which implies that fx is a family of C1 on ℝ+5 for all t≥0. Hence, by the Picard–Lindelöf theorem [42], the model system (1) admits a unique solution locally in ℝ+5 for all time t≥0. This completes the proof of the theorem.Theorem 2 (nonnegativity of model solutions).
The solution setSt,Vt,It,Dt,Rt of the model system (1) with initial condition (2) in the closed region Ω⊂ℝ+5 is nonnegative for all future time t≥0.Proof.
From the model system (1), we have that
(4)ṠS=0=Λ>0,V̇V=0=ϕS≥0,∀S≥0,İI=0=ϕS≥0,∀S≥0,ḊD=0=αI≥0,∀I≥0,ṘR=0=γI+ωD≥0,∀I,D≥0.
Hence, all the trajectories of model system (1) remain nonnegative for all nonnegative initial conditions.Theorem 3 (invariant region).
The closed regionΩ=S,I,V,D,R∈ℝ+5:0≤Nt≤Λ/μ is positively invariant and absorbing with respect to the set of nonlinear differential system of equations (1).Proof.
DifferentiatingNt with respect to time t and adding all the right-hand sides of model (1) together gives Ṅt=Λ−μN−ϵI−πD. Based on the non-negativity of model (1), we have that Nt=Λ−μN−ϵI−πD≤Λ−μN. Using Gronwall’s inequality, the solution become: Nt≤Λ/μ+N0−Λ/μe−ut. It follows that Nt≤Λ/μ as t⟶∞. Hence, Ω=S,I,V,D,R∈ℝ+5:0≤Nt≤Λ/μ. Thus, the model is biologically meaningful and mathematically well-posed in the region Ω.


### 3.1. Model Equilibria and Basic Reproduction Number

#### 3.1.1. Disease-Free Equilibrium Point of the Model

To obtain the disease-free equilibrium (DFE), we equate the right-hand sides of model system (1) equal to zero at I=D=0. Therefore, the disease-free equilibrium point of model system (1) is denoted by Ef and given by Ef=S∗,V∗,I∗,D∗,R∗=Λ/ϕ+μ,Λϕ/μμ+ϕ,0,0,0.


#### 3.1.2. Basic Reproduction Number of the Model

In this subsection, we obtained the threshold parameter that governs the spread of a COVID-19 disease which is called the basic reproduction number,R0. To obtain the basic reproduction number, we used the next-generation matrix method [43, 44] so that it is the spectral radius of the next-generation matrix. From model system (1), the rate that may cause new infections f and the rate of transfer v are
(5)f=βSIσ+ηI+ψVI0,v=μ+ϵ+γI+αI−αI+μ+π+ωD,respectively. The linearized matrix of f and v at the disease-free equilibrium point Ef is given by F and V as
(6)F=βΛσϕ+μ+ϕΛψμϕ+μ000,V=μ+ϵ+γ+α0−αω+μ+π,respectively. Then, the inverse of V is obtained and given by
(7)V−1=1ϵ+γ+μ+α0απ+ω+μα+γ+μ+ϵ1π+ω+μ.Finally,(8)FV−1=Λβμ+ϕσψμϕ+μα+μ+γ+ϵ000.The two eigenvalues ofFV−1 are
(9)λ1=0,λ2=Λβμ+ϕσψσμϕ+μα+μ+γ+ϵ.Indeed, the spectral radius ofFV−1 is
(10)ρFV−1=R0=maxλ1,λ2=Λβμ+ϕσψσμϕ+μα+μ+γ+ϵ.From (10), 1/ϕ+μ and 1/α+μ+γ+ϵ represent the time spent by susceptible and infectious individuals in states S and I, respectively.


### 3.2. Endemic Equilibrium Point of the Model

Besides the disease-free equilibrium point, the endemic equilibrium point is a positive steady-state solution where the disease persists in the population and it is calculated by setting all differential equations in model (1) equal to zero. Solving the first, second, fourth, and fifth equations in the system (1) gives as
(11)S∗∗=Λσ+ηI∗βI∗+k1σ+ηI∗,V∗∗=ϕΛμ+ψI∗βI∗+k1σ+ηI∗,I∗∗=αI∗k3,R∗∗=γk3+αωI∗μk3,Then, substitutingS∗∗ and V∗∗ into İ, we obtained the quadratic equation and I∗∗ is the positive solution calculated from
(12)AI∗∗2+BI∗∗+C=0,where
(13)k1=ϕ+μ,k2=α+γ+ϵ+μ,k3=π+ω+μ,A=k2βψ+k1ψη,B=k2βμ+μk1η+k1ψσ−βΛψ+ϕψ2Λ,C=−μk1k2σR0−1.It follows thatC<0 whenever R0>1 and the endemic equilibrium exist. Hence, the following result is established.Proposition 4.
The positive endemic equilibrium pointEp=S∗∗,V∗∗,I∗∗,D∗∗,R∗∗ of model system (1) exists when R0>1.


### 3.3. Stability Analysis of Model Equilibria

The subsequent theorems discuss the local and global stability analysis of disease-freeEf and disease-persistent equilibrium Ep points of model system (1).Theorem 5.
The disease-free equilibrium pointEf=S∗,V∗,I∗,D∗,R∗ of model system (1) is locally asymptotically stable whenever R0<1 and unstable whenever R0>1.Proof.
The linearized matrix of the system (3) at the disease-free equilibrium Ef is given by
(14)LEf=−k10−βΛσk100ϕ−μ−ψϕΛμk10000Λβμ+ϕσψσμk1−k20000α−k3000γω−μ,where k1=ϕ+μ,k2=μ+ϵ+α+γ, and k3=μ+π+ω and the characteristic polynomial of LEf is expressed by
(15)Pλ=λ+k1λ+μ2λ+k3λ−Λβμ+ϕσψσμk1−k2=0.
The five eigenvalues of (15) are
(16)λ1=−k1=−ϕ+μ<0,λ2=λ3=−μ<0,λ4=−k3=−α+γ+μ+ϵ<0,λ5=Λβμ+ϕσψσμk1−k2=Λβμ+ϕσψσμμ+ϕ−μ+γ+ϵ+α=μ+γ+ϵ+αΛβμ+ϕσψσμϕ+μμ+γ+ϵ+α−1=μ+γ+ϵ+αR0−1.
It is obvious thatλ1<0,λ2<0,λ3<0, and λ4<0, but λ5=μ+γ+ϵ+αR0−1<0 provided that R0−1<0, which leads to R0<1. Thus, the disease-free equilibrium point Ef of the system (1) is locally asymptotically stable for R0<1. In epidemilogical point of view, this result implies that COVID-19 cannot invade the community for R0<1. Whenever R0>1, the value of λ5>0, which implies that there is one eigenvalue with positive real root and the disease-free equilibrium point Ef is unstable, and this indicates that COVID-19 invasion will occur within the society.Theorem 6.
The disease-free equilibrium pointEf=S∗,V∗,I∗,D∗,R∗ of model system (1) is globally asymptotically stable whenever R0<1, otherwise unstable.Proof.
To show the above result, we follow Castillo-Chavez et al.’s theorem [45]. The model system (1) can be written as
(17)Ẋ=HX,ZŻ=MX,Z,MX,0=0,where X=S,V,R∈R3 denotes the non-infectious human populations, while Z=I,D∈R2 denotes the infectious human populations. According to Castillo-Chavez et al.’s theorem [45], to be the disease-free equilibrium point Ef of the model system (1) globally asymptotically stable, the following two necessary conditions H1 and H2 must be met:
H1 For Ẋ=HX∗,0, X∗ is globally asymptotically stable. Thus, in the system (1) at Ef, we have that limt⟶∞St,Vt,Rt=Λ/μ+ϕ,ϕΛ/μϕ+μ,0. Hence, X∗ is globally asymptotically stable.
H2MX,Z=AZ−M^X,Z, where M^X,Z≥0 in the closed region Ω.
Now, considerMX,Z=AZ−M^X,Z⇒M^X,Z=AZ−MX,Z.
(18)M^X,Z=∂IEf∂I∂IEf∂D∂DEf∂I∂DEf∂DID−dIdtdDdt,M^=βΛμ+ψσϕΛσμϕ+μ−μ+ϵ+α+γ0α−μ+π+ωID−βSIσ+ηI+ψVI−μ+ϵ+γ+αIαI−μ+π+ωD.
FromM^, A is a Metzler matrix (off diagonal elements of A are nonnegative in the closed region Ω in which the model gives biological sense).
Further simplification onM^X,Z, we get that
(19)M^=βμ+ϕσψΛIσμϕ+μ−βSσ+ηI+ψVI0≥0,since at disease-free equilibrium point Ef, the solutions of the first and the second equation in system (1) become
St=Λ/ϕ+μ+S0−Λ/ϕ+μe−ϕ+μt≤Λ/ϕ+μ and Vt=Λϕ/μϕ+μ+1/μS0−Λ/ϕ+μe−ϕ+μt≤Λϕ/μϕ+μ as t⟶∞. This implies that condition H2 is fulfilled.
As a result, the disease-free equilibrium pointEf is globally asymptotically stable.Theorem 7.
The positive endemic equilibrium pointEp=S∗∗,V∗∗,I∗∗,D∗∗,R∗∗ of model (1) is globally asymptotically stable wheneverR0>1.Proof.
To prove this theorem, we employ the Dulac criterion [46] by taking LI,D=1/ID as the candidate of Dual’s function for model system (1) and F=S,V,I,D,R.
SinceLI,D>0,∀I,D>0, we have that
(20)LS,V,I,D,R=dLFdt=dfSdt+dfVdt+dfIdt+dfDdt+dfRdt=∂∂S1IDΛ−βSIσ+ηI−k1S+∂∂V1IDϕS−μ+ψIV+∂∂I1IDβSIσ+ηI+ψVI−k2I+∂∂D1IDαI−k3D+∂∂R1IDγI+ωD−μR=−βDσ+ηI−k1ID−μ+ψIID−αD2−μID−βηSσ+ηI2<0.According to Dual’s criterion [46], the model system (1) does not exhibit any periodic solution in the region Ω and, hence, the disease-persistent equilibrium point Ep is globally asysmptotically stable whenever R0>1. This finding from an epidemiological perspective suggests that the number of infectious agents may fluctuate, making it challenging to allocate resources for disease control.


### 3.4. Sensitivity Analysis of Model Parameters

In this subsection, we investigate sensitivity analysis used to discover the most dominant parameters for the spreading out as well as control of infection in the society. To go through sensitivity analysis, we follow the classical approach used by Chitnis et al. [47].Definition 8.
The normalized forward sensitivity index of variableR0 that depends differentiability on an index parameter v is defined as SvR0=v/R0.dR0/dv [47].Here, for example,SvR0≅+0.1 means increasing (or decreasing) the value of v by 10%, which always results to increase (or decrease) the value of R0 by 10%, whereas SvR0=0.1 means increasing (or decreasing) the value of v by 10%, which always results to decrease (or increase) R0 by 10%.We assessed the sensitivity indices of the system (1) in a similar manner, and the results are as follows:
(21)SβR0=1,SβR0=ββμ+ϕσψ,SϕR0=ϕμδψ−βϕ+μσψϕ+μβ,SαR0=−αα+μ+γ+ε,SϕR0=ϕμδψ−βϕ+μσψϕ+μβ,SδR0=−βμσβμ+ϕψσ,SψR0=ϕμψβμ+ϕψσ.Numerical results of sensitivity analysis from baseline model parameters are arranged in a descending order from Table2. However, in the study of sensitivity analysis, it is not biologically reasonable to increase or decrease human mortality rates to control disease outbreak, and hence we do not consider.Table 2
Baseline values for model parameters of the model system.
ParameterNominal valueReferenceElasticity indexΛ4721.5[31]+1β0.0108[31]+0.847ψ0.01Assumed+0.153ϕ0.0043Assumed-0.0132γ0.1[3]-0.219α0.22Assumed-0.481σ5[48]-1ϵ0.0271[31]Noπ0.4[49]Noη1[48]Noμ0.11[30]Noω0.2AssumedNo


### 3.5. Interpretation of Sensitivity Indices

It is noted from Table2 that the value of R0 increases when the parameter values Λ, β, and ψ increase while other parameter values are kept fixed. This means the spread of COVID-19 increases due to the reason that the basic reproduction number increases as their values increase; it means that the average number of secondary cases of infection increases in the community while the rest of parameters kept fixed. In contrast, parameters σ,α,ϕ, and γ have a positive impact in the reduction of the burden of the disease in the community as their values increase while the others are kept constant. And also, as their values increase, the basic reproduction number decreases. This means the spread of the virus decreases since raising the values of these parameters will consequently decrease the number of total infected population.


## 3.1. Model Equilibria and Basic Reproduction Number

### 3.1.1. Disease-Free Equilibrium Point of the Model

To obtain the disease-free equilibrium (DFE), we equate the right-hand sides of model system (1) equal to zero at I=D=0. Therefore, the disease-free equilibrium point of model system (1) is denoted by Ef and given by Ef=S∗,V∗,I∗,D∗,R∗=Λ/ϕ+μ,Λϕ/μμ+ϕ,0,0,0.


### 3.1.2. Basic Reproduction Number of the Model

In this subsection, we obtained the threshold parameter that governs the spread of a COVID-19 disease which is called the basic reproduction number,R0. To obtain the basic reproduction number, we used the next-generation matrix method [43, 44] so that it is the spectral radius of the next-generation matrix. From model system (1), the rate that may cause new infections f and the rate of transfer v are
(5)f=βSIσ+ηI+ψVI0,v=μ+ϵ+γI+αI−αI+μ+π+ωD,respectively. The linearized matrix of f and v at the disease-free equilibrium point Ef is given by F and V as
(6)F=βΛσϕ+μ+ϕΛψμϕ+μ000,V=μ+ϵ+γ+α0−αω+μ+π,respectively. Then, the inverse of V is obtained and given by
(7)V−1=1ϵ+γ+μ+α0απ+ω+μα+γ+μ+ϵ1π+ω+μ.Finally,(8)FV−1=Λβμ+ϕσψμϕ+μα+μ+γ+ϵ000.The two eigenvalues ofFV−1 are
(9)λ1=0,λ2=Λβμ+ϕσψσμϕ+μα+μ+γ+ϵ.Indeed, the spectral radius ofFV−1 is
(10)ρFV−1=R0=maxλ1,λ2=Λβμ+ϕσψσμϕ+μα+μ+γ+ϵ.From (10), 1/ϕ+μ and 1/α+μ+γ+ϵ represent the time spent by susceptible and infectious individuals in states S and I, respectively.


## 3.1.1. Disease-Free Equilibrium Point of the Model

To obtain the disease-free equilibrium (DFE), we equate the right-hand sides of model system (1) equal to zero at I=D=0. Therefore, the disease-free equilibrium point of model system (1) is denoted by Ef and given by Ef=S∗,V∗,I∗,D∗,R∗=Λ/ϕ+μ,Λϕ/μμ+ϕ,0,0,0.


## 3.1.2. Basic Reproduction Number of the Model

In this subsection, we obtained the threshold parameter that governs the spread of a COVID-19 disease which is called the basic reproduction number,R0. To obtain the basic reproduction number, we used the next-generation matrix method [43, 44] so that it is the spectral radius of the next-generation matrix. From model system (1), the rate that may cause new infections f and the rate of transfer v are
(5)f=βSIσ+ηI+ψVI0,v=μ+ϵ+γI+αI−αI+μ+π+ωD,respectively. The linearized matrix of f and v at the disease-free equilibrium point Ef is given by F and V as
(6)F=βΛσϕ+μ+ϕΛψμϕ+μ000,V=μ+ϵ+γ+α0−αω+μ+π,respectively. Then, the inverse of V is obtained and given by
(7)V−1=1ϵ+γ+μ+α0απ+ω+μα+γ+μ+ϵ1π+ω+μ.Finally,(8)FV−1=Λβμ+ϕσψμϕ+μα+μ+γ+ϵ000.The two eigenvalues ofFV−1 are
(9)λ1=0,λ2=Λβμ+ϕσψσμϕ+μα+μ+γ+ϵ.Indeed, the spectral radius ofFV−1 is
(10)ρFV−1=R0=maxλ1,λ2=Λβμ+ϕσψσμϕ+μα+μ+γ+ϵ.From (10), 1/ϕ+μ and 1/α+μ+γ+ϵ represent the time spent by susceptible and infectious individuals in states S and I, respectively.


## 3.2. Endemic Equilibrium Point of the Model

Besides the disease-free equilibrium point, the endemic equilibrium point is a positive steady-state solution where the disease persists in the population and it is calculated by setting all differential equations in model (1) equal to zero. Solving the first, second, fourth, and fifth equations in the system (1) gives as
(11)S∗∗=Λσ+ηI∗βI∗+k1σ+ηI∗,V∗∗=ϕΛμ+ψI∗βI∗+k1σ+ηI∗,I∗∗=αI∗k3,R∗∗=γk3+αωI∗μk3,Then, substitutingS∗∗ and V∗∗ into İ, we obtained the quadratic equation and I∗∗ is the positive solution calculated from
(12)AI∗∗2+BI∗∗+C=0,where
(13)k1=ϕ+μ,k2=α+γ+ϵ+μ,k3=π+ω+μ,A=k2βψ+k1ψη,B=k2βμ+μk1η+k1ψσ−βΛψ+ϕψ2Λ,C=−μk1k2σR0−1.It follows thatC<0 whenever R0>1 and the endemic equilibrium exist. Hence, the following result is established.Proposition 4.
The positive endemic equilibrium pointEp=S∗∗,V∗∗,I∗∗,D∗∗,R∗∗ of model system (1) exists when R0>1.


## 3.3. Stability Analysis of Model Equilibria

The subsequent theorems discuss the local and global stability analysis of disease-freeEf and disease-persistent equilibrium Ep points of model system (1).Theorem 5.
The disease-free equilibrium pointEf=S∗,V∗,I∗,D∗,R∗ of model system (1) is locally asymptotically stable whenever R0<1 and unstable whenever R0>1.Proof.
The linearized matrix of the system (3) at the disease-free equilibrium Ef is given by
(14)LEf=−k10−βΛσk100ϕ−μ−ψϕΛμk10000Λβμ+ϕσψσμk1−k20000α−k3000γω−μ,where k1=ϕ+μ,k2=μ+ϵ+α+γ, and k3=μ+π+ω and the characteristic polynomial of LEf is expressed by
(15)Pλ=λ+k1λ+μ2λ+k3λ−Λβμ+ϕσψσμk1−k2=0.
The five eigenvalues of (15) are
(16)λ1=−k1=−ϕ+μ<0,λ2=λ3=−μ<0,λ4=−k3=−α+γ+μ+ϵ<0,λ5=Λβμ+ϕσψσμk1−k2=Λβμ+ϕσψσμμ+ϕ−μ+γ+ϵ+α=μ+γ+ϵ+αΛβμ+ϕσψσμϕ+μμ+γ+ϵ+α−1=μ+γ+ϵ+αR0−1.
It is obvious thatλ1<0,λ2<0,λ3<0, and λ4<0, but λ5=μ+γ+ϵ+αR0−1<0 provided that R0−1<0, which leads to R0<1. Thus, the disease-free equilibrium point Ef of the system (1) is locally asymptotically stable for R0<1. In epidemilogical point of view, this result implies that COVID-19 cannot invade the community for R0<1. Whenever R0>1, the value of λ5>0, which implies that there is one eigenvalue with positive real root and the disease-free equilibrium point Ef is unstable, and this indicates that COVID-19 invasion will occur within the society.Theorem 6.
The disease-free equilibrium pointEf=S∗,V∗,I∗,D∗,R∗ of model system (1) is globally asymptotically stable whenever R0<1, otherwise unstable.Proof.
To show the above result, we follow Castillo-Chavez et al.’s theorem [45]. The model system (1) can be written as
(17)Ẋ=HX,ZŻ=MX,Z,MX,0=0,where X=S,V,R∈R3 denotes the non-infectious human populations, while Z=I,D∈R2 denotes the infectious human populations. According to Castillo-Chavez et al.’s theorem [45], to be the disease-free equilibrium point Ef of the model system (1) globally asymptotically stable, the following two necessary conditions H1 and H2 must be met:
H1 For Ẋ=HX∗,0, X∗ is globally asymptotically stable. Thus, in the system (1) at Ef, we have that limt⟶∞St,Vt,Rt=Λ/μ+ϕ,ϕΛ/μϕ+μ,0. Hence, X∗ is globally asymptotically stable.
H2MX,Z=AZ−M^X,Z, where M^X,Z≥0 in the closed region Ω.
Now, considerMX,Z=AZ−M^X,Z⇒M^X,Z=AZ−MX,Z.
(18)M^X,Z=∂IEf∂I∂IEf∂D∂DEf∂I∂DEf∂DID−dIdtdDdt,M^=βΛμ+ψσϕΛσμϕ+μ−μ+ϵ+α+γ0α−μ+π+ωID−βSIσ+ηI+ψVI−μ+ϵ+γ+αIαI−μ+π+ωD.
FromM^, A is a Metzler matrix (off diagonal elements of A are nonnegative in the closed region Ω in which the model gives biological sense).
Further simplification onM^X,Z, we get that
(19)M^=βμ+ϕσψΛIσμϕ+μ−βSσ+ηI+ψVI0≥0,since at disease-free equilibrium point Ef, the solutions of the first and the second equation in system (1) become
St=Λ/ϕ+μ+S0−Λ/ϕ+μe−ϕ+μt≤Λ/ϕ+μ and Vt=Λϕ/μϕ+μ+1/μS0−Λ/ϕ+μe−ϕ+μt≤Λϕ/μϕ+μ as t⟶∞. This implies that condition H2 is fulfilled.
As a result, the disease-free equilibrium pointEf is globally asymptotically stable.Theorem 7.
The positive endemic equilibrium pointEp=S∗∗,V∗∗,I∗∗,D∗∗,R∗∗ of model (1) is globally asymptotically stable wheneverR0>1.Proof.
To prove this theorem, we employ the Dulac criterion [46] by taking LI,D=1/ID as the candidate of Dual’s function for model system (1) and F=S,V,I,D,R.
SinceLI,D>0,∀I,D>0, we have that
(20)LS,V,I,D,R=dLFdt=dfSdt+dfVdt+dfIdt+dfDdt+dfRdt=∂∂S1IDΛ−βSIσ+ηI−k1S+∂∂V1IDϕS−μ+ψIV+∂∂I1IDβSIσ+ηI+ψVI−k2I+∂∂D1IDαI−k3D+∂∂R1IDγI+ωD−μR=−βDσ+ηI−k1ID−μ+ψIID−αD2−μID−βηSσ+ηI2<0.According to Dual’s criterion [46], the model system (1) does not exhibit any periodic solution in the region Ω and, hence, the disease-persistent equilibrium point Ep is globally asysmptotically stable whenever R0>1. This finding from an epidemiological perspective suggests that the number of infectious agents may fluctuate, making it challenging to allocate resources for disease control.


## 3.4. Sensitivity Analysis of Model Parameters

In this subsection, we investigate sensitivity analysis used to discover the most dominant parameters for the spreading out as well as control of infection in the society. To go through sensitivity analysis, we follow the classical approach used by Chitnis et al. [47].Definition 8.
The normalized forward sensitivity index of variableR0 that depends differentiability on an index parameter v is defined as SvR0=v/R0.dR0/dv [47].Here, for example,SvR0≅+0.1 means increasing (or decreasing) the value of v by 10%, which always results to increase (or decrease) the value of R0 by 10%, whereas SvR0=0.1 means increasing (or decreasing) the value of v by 10%, which always results to decrease (or increase) R0 by 10%.We assessed the sensitivity indices of the system (1) in a similar manner, and the results are as follows:
(21)SβR0=1,SβR0=ββμ+ϕσψ,SϕR0=ϕμδψ−βϕ+μσψϕ+μβ,SαR0=−αα+μ+γ+ε,SϕR0=ϕμδψ−βϕ+μσψϕ+μβ,SδR0=−βμσβμ+ϕψσ,SψR0=ϕμψβμ+ϕψσ.Numerical results of sensitivity analysis from baseline model parameters are arranged in a descending order from Table2. However, in the study of sensitivity analysis, it is not biologically reasonable to increase or decrease human mortality rates to control disease outbreak, and hence we do not consider.Table 2
Baseline values for model parameters of the model system.
ParameterNominal valueReferenceElasticity indexΛ4721.5[31]+1β0.0108[31]+0.847ψ0.01Assumed+0.153ϕ0.0043Assumed-0.0132γ0.1[3]-0.219α0.22Assumed-0.481σ5[48]-1ϵ0.0271[31]Noπ0.4[49]Noη1[48]Noμ0.11[30]Noω0.2AssumedNo


## 3.5. Interpretation of Sensitivity Indices

It is noted from Table2 that the value of R0 increases when the parameter values Λ, β, and ψ increase while other parameter values are kept fixed. This means the spread of COVID-19 increases due to the reason that the basic reproduction number increases as their values increase; it means that the average number of secondary cases of infection increases in the community while the rest of parameters kept fixed. In contrast, parameters σ,α,ϕ, and γ have a positive impact in the reduction of the burden of the disease in the community as their values increase while the others are kept constant. And also, as their values increase, the basic reproduction number decreases. This means the spread of the virus decreases since raising the values of these parameters will consequently decrease the number of total infected population.


## 4. Extension of the Model into Optimal Control

In this section, we consider the optimal control of model (1), which describes the interaction of educational campaign with susceptibles and pharmaceutical interventions of the infected individuals. The ultimate goal of optimal control is to minimize the total number of infected individuals with minimum cost. We perform optimal control to the system (1) to investigate the impact of implementing continuous educational awareness creation and pharmaceutical intervention to mitigate the COVID-19 outbreak. To proceed with this, we introduced a set of time-dependent control variables u1t and u2t where
(i)
u1t denotes the efforts intended for susceptible individuals to create awareness through all media outlets(ii)
u2t denotes the efforts intended to encourage continuous treatment usually for quarantined individualsFor simplicity of notation, we representu1t=u1 and u2t=u2. After incorporating the above intervention strategies in model (1), we get an optimal control model described by
(22)Ṡt=Λ−1−u1βSIσ+ηI−ϕS−μS,V̇t=ϕS−μ+ψIV,İt=1−u1βSIσ+ηI+ψVI−μ+ϵ+γI−αI,Ḋt=αI−μ+π+ω+u2D,Ṙt=γI+ω+u2D−μR.The optimal control variablesu1 and u2 minimize the objective function subject to system (23). The objective function is defined as
(23)Ju1,u2=minu1,u2∫0tfA1Iτ+A2Dτ+B1u12τ+B2u22τdτ,where tf is fixed final time and the positive coefficients A1,A2,B1, and B2 are constants that balance the cost and the number of infected individuals at time t. In this paper, since the cost of the controls may not be measured linearly, we use the quadratic form as frequently used in [50, 51]. The ultimate goal is to find the optimal control u1∗,u2∗ such that
(24)Ju1∗,u2∗=minJu1,u2:u1,u2∈U.The controlsu1 and u2 are assumed to be at least Lebesgue measurable on 0,tf [52].


### 4.1. Existence of Optimal Control Problem

Theorem 9.
There exists an optimal controlu1∗,u2∗ that minimizes the objective function Ju1,u2, subject to the control system (22).Proof.
To show the existence of optimal control, we denote the right-hand sides of system (22) by zt,x⟶,u⟶. Then, to prove the existence of optimal control, we follow the same procedure as [53, 54]. To achieve this result, the following conditions must be met:
(1)
z is of class C1, and there exists a constant c such that(25)zt,0,0≤c,zx⟶t,x⟶,u⟶≤c1+u⟶,zu⟶t,x⟶,u⟶≤c(2)
The admissible set of all solutions to system (22) with corresponding control in Uad is nonempty(3)
There exist functionsd1 and d2 such that zt,x⟶,u⟶=d1t,x⟶+d2t,x⟶u⟶(4)
The control setU=0,1×0,1 is closed, convex, and compact(5)
The integrand of the objective function is convex inU
To verify condition (1), rewrite
(26)zt,x⟶,u⟶=Λ−1−u1βSIσ+ηI−μ+ϕSϕS−μ+ψIV1−u1βSIσ+ηI+ψVI−ψ+ϵ+γ+αIαI−u2+π+μ+ωDγI+ω+u2D−μR.
Then, one can see thatzt,x⟶,u⟶ is of class C1 and zt,0,0=Λ.
Moreover, consider(27)zx⟶t,x⟶,u⟶=−1−u1βIσ+ηI−μ+ϕ0−1−u1βσSσ+ηI200ϕ−μ+ψI−ϕV001−u1βSσ+ηI1−u1βσIσ+ηI2+ψV−k100000α−u2+k2000γω+u2−μ,where k1=α+γ+ψ+ϵ,k2=π+μ+ω, and
(28)zu⟶t,x⟶,u⟶=βSIσ+ηI000−βSIσ+ηI00−D0D.
SinceS,V,I,D, and R are bounded, hence, there exists a constant c such that
(29)zt,0,0≤c,zx⟶t,x⟶,u⟶≤c1+u⟶,zu⟶t,x⟶,u⟶≤c.
This shows that condition (1) is satisfied.
According to condition (1), a unique solution for the constant controls exist, which will improve that condition (2), is met.
Beside this,(30)zt,x⟶,u⟶=Λ−βSIσ+ηI−μ+ϕSϕS−μ+ψIVβSIσ+ηI+ψVI−ψ+ϵ+γ+αIαI−π+μ+ωDγI+ωD−μR+βSIσ+ηI0−βSIσ+ηI00−D0Du1u2.
This means condition (3) holds. The control set U is closed and bounded, hence compact. This verifies condition (10). Finally, we try to verify condition (15), (i.e., the convexity of the integrand of the objective function). We must ensure that for any two values u⟶ and k⟶ of the control vector and a constant m∈0,1, the following inequality holds:
(31)1−mft,x⟶,u⟶+mft,x⟶,k⟶≥ft,x⟶,1−mu⟶+mk⟶,where
(32)ft,x⟶,u⟶=A1I+A2D+B1u12+B2u22.
Furthermore,(33)1−mft,x⟶,u⟶+mft,x⟶,k⟶=A1I+A2D+1−mB1u12+B2u22+mB1k12+B2k22,and also,
(34)ft,x⟶,1−mu⟶+mk⟶=A1I+A2D+B11−mu1+mk12+B21−mu2+mk22.
Then,(35)1−mft,x⟶,u⟶+mft,x⟶,k⟶−ft,x⟶,1−mu⟶+mk⟶=1−mB1u21+B2u22+mB1k21+B2k22−B11−mu1+mk12+B21−mu2+mk22=B11−mu21+mk21−1−mu1+mk12+B21−mu22+mk22−1−mu2+mk22=B1m1−mu1−k12+B2m1−mu2−k22≥0.
This shows that condition (15) is fulfilled.
Hence, the proof of the theorem is completed.


### 4.2. Characterization of the Optimal Control

Now, to solve the optimal control problem (23), we have constructed the Lagrangian equation as follows:
(36)Lt,x,u=dJdt=A1It+A2Dt+B1u12t+B2u22t.The necessary conditions that an optimal control should satisfy are derived from Pontryagin’s Maximum Principle [33, 34].To minimize the Lagrangian, we have constructed the Hamiltonian equation by(37)Ht,x,u,ξ=dJdt+∑i=15ξifi=Lt,x,u+∑i=15ξifi,where Lt,x,u=A1I+A2D+B1u12t+B2u22t and fi,i=1,2,3,4,5 are the right-hand side components of (22) and ξ=ξ1,ξ2,ξ3,ξ4,ξ5∈ℝ5 are costate variables associated to the state variables S,L,A,I, and R which can be determined from second partial derivatives of H with respect to state variables as follows:
(38)dξ1dt=−∂H∂S,dξ2dt=−∂H∂V,dξ3dt=−∂H∂I,dξ4dt=−∂H∂D,dξ5dt=−∂H∂R.It follows that(39)ξ1̇=1−u1βIσ+ηI+α1ξ1−ϕξ2−1−u1βIξ3σ+ηI,ξ2̇=μ+ψIξ2−ψIξ3,ξ3̇=−A1+1−u1βδξ1δ+ηI2+ψVξ2+α2−ψV−1−u1βσSσ+ηI2ξ3−αξ4−γξ5,ξ4̇=A2+α3+u2ξ4−ω+u2ξ5,ξ5̇=μξ5,where a1=ϕ+μ,a2=α+μ+γ+ϵ, and a3=μ+π+ω, with transversality conditions:
(40)ξ1tf=ξ2tf=ξ3tf=ξ4tf=ξ5tf=0.Theorem 10.
The optimal control set {u1∗t, u2∗t} that minimizes Ju1,u2 over the region U is obtained and presented in compact form as
(41)u1t=min1,max0,u∆,u2t=min1,max0,u∇,where
(42)u∆=βSIξ3−ξ12B1σ+ηI,u∇=Dξ4−ξ52B2.Proof.
The proof is similar to [55].


## 4.1. Existence of Optimal Control Problem

Theorem 9.
There exists an optimal controlu1∗,u2∗ that minimizes the objective function Ju1,u2, subject to the control system (22).Proof.
To show the existence of optimal control, we denote the right-hand sides of system (22) by zt,x⟶,u⟶. Then, to prove the existence of optimal control, we follow the same procedure as [53, 54]. To achieve this result, the following conditions must be met:
(1)
z is of class C1, and there exists a constant c such that(25)zt,0,0≤c,zx⟶t,x⟶,u⟶≤c1+u⟶,zu⟶t,x⟶,u⟶≤c(2)
The admissible set of all solutions to system (22) with corresponding control in Uad is nonempty(3)
There exist functionsd1 and d2 such that zt,x⟶,u⟶=d1t,x⟶+d2t,x⟶u⟶(4)
The control setU=0,1×0,1 is closed, convex, and compact(5)
The integrand of the objective function is convex inU
To verify condition (1), rewrite
(26)zt,x⟶,u⟶=Λ−1−u1βSIσ+ηI−μ+ϕSϕS−μ+ψIV1−u1βSIσ+ηI+ψVI−ψ+ϵ+γ+αIαI−u2+π+μ+ωDγI+ω+u2D−μR.
Then, one can see thatzt,x⟶,u⟶ is of class C1 and zt,0,0=Λ.
Moreover, consider(27)zx⟶t,x⟶,u⟶=−1−u1βIσ+ηI−μ+ϕ0−1−u1βσSσ+ηI200ϕ−μ+ψI−ϕV001−u1βSσ+ηI1−u1βσIσ+ηI2+ψV−k100000α−u2+k2000γω+u2−μ,where k1=α+γ+ψ+ϵ,k2=π+μ+ω, and
(28)zu⟶t,x⟶,u⟶=βSIσ+ηI000−βSIσ+ηI00−D0D.
SinceS,V,I,D, and R are bounded, hence, there exists a constant c such that
(29)zt,0,0≤c,zx⟶t,x⟶,u⟶≤c1+u⟶,zu⟶t,x⟶,u⟶≤c.
This shows that condition (1) is satisfied.
According to condition (1), a unique solution for the constant controls exist, which will improve that condition (2), is met.
Beside this,(30)zt,x⟶,u⟶=Λ−βSIσ+ηI−μ+ϕSϕS−μ+ψIVβSIσ+ηI+ψVI−ψ+ϵ+γ+αIαI−π+μ+ωDγI+ωD−μR+βSIσ+ηI0−βSIσ+ηI00−D0Du1u2.
This means condition (3) holds. The control set U is closed and bounded, hence compact. This verifies condition (10). Finally, we try to verify condition (15), (i.e., the convexity of the integrand of the objective function). We must ensure that for any two values u⟶ and k⟶ of the control vector and a constant m∈0,1, the following inequality holds:
(31)1−mft,x⟶,u⟶+mft,x⟶,k⟶≥ft,x⟶,1−mu⟶+mk⟶,where
(32)ft,x⟶,u⟶=A1I+A2D+B1u12+B2u22.
Furthermore,(33)1−mft,x⟶,u⟶+mft,x⟶,k⟶=A1I+A2D+1−mB1u12+B2u22+mB1k12+B2k22,and also,
(34)ft,x⟶,1−mu⟶+mk⟶=A1I+A2D+B11−mu1+mk12+B21−mu2+mk22.
Then,(35)1−mft,x⟶,u⟶+mft,x⟶,k⟶−ft,x⟶,1−mu⟶+mk⟶=1−mB1u21+B2u22+mB1k21+B2k22−B11−mu1+mk12+B21−mu2+mk22=B11−mu21+mk21−1−mu1+mk12+B21−mu22+mk22−1−mu2+mk22=B1m1−mu1−k12+B2m1−mu2−k22≥0.
This shows that condition (15) is fulfilled.
Hence, the proof of the theorem is completed.


## 4.2. Characterization of the Optimal Control

Now, to solve the optimal control problem (23), we have constructed the Lagrangian equation as follows:
(36)Lt,x,u=dJdt=A1It+A2Dt+B1u12t+B2u22t.The necessary conditions that an optimal control should satisfy are derived from Pontryagin’s Maximum Principle [33, 34].To minimize the Lagrangian, we have constructed the Hamiltonian equation by(37)Ht,x,u,ξ=dJdt+∑i=15ξifi=Lt,x,u+∑i=15ξifi,where Lt,x,u=A1I+A2D+B1u12t+B2u22t and fi,i=1,2,3,4,5 are the right-hand side components of (22) and ξ=ξ1,ξ2,ξ3,ξ4,ξ5∈ℝ5 are costate variables associated to the state variables S,L,A,I, and R which can be determined from second partial derivatives of H with respect to state variables as follows:
(38)dξ1dt=−∂H∂S,dξ2dt=−∂H∂V,dξ3dt=−∂H∂I,dξ4dt=−∂H∂D,dξ5dt=−∂H∂R.It follows that(39)ξ1̇=1−u1βIσ+ηI+α1ξ1−ϕξ2−1−u1βIξ3σ+ηI,ξ2̇=μ+ψIξ2−ψIξ3,ξ3̇=−A1+1−u1βδξ1δ+ηI2+ψVξ2+α2−ψV−1−u1βσSσ+ηI2ξ3−αξ4−γξ5,ξ4̇=A2+α3+u2ξ4−ω+u2ξ5,ξ5̇=μξ5,where a1=ϕ+μ,a2=α+μ+γ+ϵ, and a3=μ+π+ω, with transversality conditions:
(40)ξ1tf=ξ2tf=ξ3tf=ξ4tf=ξ5tf=0.Theorem 10.
The optimal control set {u1∗t, u2∗t} that minimizes Ju1,u2 over the region U is obtained and presented in compact form as
(41)u1t=min1,max0,u∆,u2t=min1,max0,u∇,where
(42)u∆=βSIξ3−ξ12B1σ+ηI,u∇=Dξ4−ξ52B2.Proof.
The proof is similar to [55].


## 5. Numerical Results and Discussion

In this section, we carried out some numerical results of the model to verify the theoretical results by using MATLAB ode45 solver with the following initial conditions:S0=120,000, V0=40,000, I0=35,000, D0=5000, and R0=40,000. Here, it is important to note that the initial conditions and parameter values used for simulation purpose are taken from the previous published articles and some others are assumed.Figure2 shows that all trajectories converge to their disease-free equilibrium point component when ψ=0.0000007,ϕ=0.043,ε=0.271, and other parameters values are from Table 1. This result supported the analytical conclusion that the disease-free equilibrium point is stable for R0<1 stated from Theorem 5 and the outbreak of COVID-19 may reduce within the society. Figure 3 illustrates that all solution curves approach to their endemic equilibrium point component when parameter values are β=0.408,α=0.32,π=0.34,ϕ=0.063,σ=2, andε=0.0271 and other parameter values are taken from Table 2. This result confirmed with the theoretical result that the COVID-19 persistent equilibrium point is globally asymptotically stable for R0=1.1232>1 stated in Theorem 6 and reinvasion of the COVID-19 outbreak occurs within the population.Figure 2
R0=0.0132<1.Figure 3
R0=1.1232>1.Figures4–7 show a time-serious plot of the total number of susceptible, immunized, infectious, and recovered human populations, with varying initial values, at the disease-free equilibrium point, where the basic reproduction number is less than unity. From Figure 8, it can be seen that large value of σ or ϕ and small value of ψ in the presence of β can lead to a small value of R0. That is to say, if we increase the progression rate of susceptibles into vaccinated and strength of other effective COVID-19 prevention measures, then the transmission of COVID-19 outbreak will end. Figure 9 depicts the physical outlook of the sensitivity analysis when the following values are used: ϕ=0.063,α=0.82,γ=0.51, and ψ=0.999 while other parameter values are from Table 2 of system of equations (1), which means that parameters with positive indices can cause an outbreak of the COVID-19 pandemic, while parameters with negative indices can slow or stop its spread.Figure 4
The time series plot of the susceptible population.Figure 5
The time series plot of the vaccinated population.Figure 6
The time series plot of the infectious population.Figure 7
The time series plot of the recovered population.Figure 8
The effect ofσ,ϕ, and ψ on R0 in the presence of β.Figure 9
The numerical result of sensitivity analysis.From Figure10, it is clear that R0 rises in proportion to the rise in contact rate β and after a certain value of β, R0 becomes larger than 1. This implies that up to a certain value of R0, the disease-free equilibrium point is stable and for R0 greater than 1, the disease-free equilibrium becomes unstable. Figure 11 reveals the variation of R0 with respect to mitigation measures of the COVID-19 outbreak σ. It can easily be observed that R0 decreases with an increase in COVID-19 mitigation rate σ and after a certain value of σ, R0 becomes less than 1. It implies after a certain value of σ, the disease-free equilibrium becomes stable.Figure 10
The effect ofβ on R0.Figure 11
The effect ofσ on R0.Figures12 and 13 present that the progression rate of the number of susceptible individuals into infectious populations due to COVID-19 infection decreases while the number of the infectious individual class decreases when the values of σ and η increase.Figure 12
The effect ofσ and η on the number of susceptible population.Figure 13
The effect ofσ and η on the number of infectious population.The authors find some significant results which are shown in Figures14–17 that the number of susceptible and recovered individuals increases with the optimal control and decreases without optimal control, whereas the number of infectious and quarantined individuals decreases with optimal control and increases without optimal control.Figure 14
Susceptible individuals with and without prevention strategyu1.Figure 15
Infectious individuals with and without prevention strategyu1 and intervention strategy u2.Figure 16
Quarantined individuals with and without prevention strategyu1 and intervention strategy u2.Figure 17
Recovered individuals with and without intervention strategyu2.The profile of optimal controlsu1 and u2 is portrayed in Figure 18. From Figure 18, we observe that to control the spread of COVID-19, the use of prevention,u1, and intervention u2 strategies must be held effectively.Figure 18
The profile of optimal controlsu1 and u2.


## 6. Conclusion and Future Directions

In this study, the authors formulated and analyzed a nonlinear deterministic mathematical model of the COVID-19 outbreak to investigate the impact of vaccination and treatment with optimal control. We have shown the existence and uniqueness, nonnegativity, and invariant region of model solutions under a certain meaningful set. Basic reproduction numberR0 is derived via the help of next-generation approach method, and stability analysis for disease-free equilibrium and endemic equilibrium points is rigorously studied based on this R0.It is observed that the disease-free equilibrium is locally as well as globally asymptotically stable whenever the basic reproduction number is less than unity (i.e., COVID-19 dies out in the population). However, the endemic equilibrium point exists and globally asymptotically stable whenever the basic reproduction number is greater than unity (i.e., COVID-19 persists within the society). Sensitivity analysis of the basic reproduction number with respect to model parameters was demonstrated, and it can be inferred that to minimize the COVID-19 outbreak, parameters with positive indices should be decreased while parameters with negative indices should be increased.The theoretical and numerical studies of the model without control reveal that as the effective mitigations of COVID-19 rates (σ and η) increase, the infected number of individuals decreases in addition to vaccination and treatment. Further, as the contact rate of susceptible and vaccinated individuals with infected individuals decreases, then the number of infected individuals decreases. Then, the model is extended by introducing two time-dependent control variables, educational campaign u1 and continuous treatment u2. The qualitative study proved the existence and characterization of an optimal control solution.The numerical results of the optimal control are illustrated and revealed that the educational campaignu1 and treatment u2 control strategies are the most effective strategies in the control of the COVID-19 pandemic. Therefore, authorities and health practitioners have to give great emphasis for the educational campaign to susceptibles and continuous treatment for quarantined individuals in order to eliminate the COVID-19 pandemic. The impact of asymptomatic of COVID-19 human immigrants and cost-effectiveness analysis of the optimal control strategies would be considered in our future study.



---

*Source: 8570311-2023-03-03.xml*