# Optimal Control and Bifurcation Analysis of HIV Model

Authors: Kumama Regassa Cheneke
Journal: Computational and Mathematical Methods in Medicine
 (2023)

---

## Abstract

In this study, a very crucial stage of HIV extinction and invisibility stages are considered and a modified mathematical model is developed to describe the dynamics of infection. Moreover, the basic reproduction numberR0 is computed using the next-generation matrix method whereas the stability of disease-free equilibrium is investigated using the eigenvalue matrix stability theory. Furthermore, if R0≤1, the disease-free equilibrium is stable both locally and globally whereas if R0>1, based on the forward bifurcation behavior, the endemic equilibrium is locally and globally asymptotically stable. Particularly, at the critical point R0=1, the model exhibits forward bifurcation behavior. On the other hand, the optimal control problem is constructed and Pontryagin’s maximum principle is applied to form an optimality system. Further, forward fourth-order Runge–Kutta’s method is applied to obtain the solution of state variables whereas Runge–Kutta’s fourth-order backward sweep method is applied to obtain solution of adjoint variables. Finally, three control strategies are considered and a cost-effective analysis is performed to identify the better strategies for HIV transmission and progression. In advance, prevention control measure is identified to be the better strategy over treatment control if applied earlier and effectively. Additionally, MATLAB simulations were performed to describe the population’s dynamic behavior.


---

## Body

## 1. Introduction

Human immunodeficiency virus (HIV) is a virus causing HIV infection and results in economic and life devastating crisis if immediate action is not taken to halt further prevalence of the infection [1–11]. Moreover, this infection has no curing medication, but antiretroviral therapy (ART) or its combination is used for halting further progression of infection by inhibiting the virus in human blood, but if left untreated leads to a sever stage called acquired immunodeficiency syndrome (AIDS) [12–23]. HIV infection progress through stages as follows: (i) primary stage (asymptomatic stage): this stage faces human individual where the virus is in the blood cannot be diagnosed with medical instruments; (ii) asymptomatic stage: this stage is symptomless stage of HIV infection but diagnosable with medical test; (iii) symptomatic stage: in this stage, the symptom of HIV infection like tiredness, loss of weight, and extreme loss of water starts to manifest in the life of HIV infected individuals; and (iv) AIDS stage: this is advanced stage of HIV infection where it is difficult for treatment and leads to death soon if special care is not taken [24–32]. Modes of HIV transmission are through unsafe sexual practices with HIV-infected person, through contacts of normal blood with HIV infected blood, mother to child through breast or birth time, and any contacts of HIV contains fluids of human’s with HIV-negative human fluids [33–35]. However, in safe practice, the risk of HIV transmission can be reduced by using principle of abstinence–be faithful–use condom (ABC) [36–38]. Optimal control intervention through public health education, using of condom and treatment benefits a lot for continuing human life in safe [1–3, 7, 11, 15, 19, 25–27, 33, 39–44].A mathematical model is a crucial scientific representation of both biological and physical problems in the form of mathematical equations [45–52]. Moreover, different mathematical models have been developed to describe the dynamics of HIV with optimal controls. Particularly, a SICA mathematical model of HIV with ART effect on HIV patients is described in [9]. But, a mathematical model of HIV infection transmission with undetectable behavior of individuals has not been studied. Hence, in this study, we are motivated to highlight the significance of adherence of ART with optimal control intervention and cost-effectiveness analysis in a modified HIV model. Further, the effect of early starting of ART that leads to an undetectable level of viral load in the blood of humans is considered and the SWIUA mathematical model of HIV with optimal control is developed. Organization of the paper is follows: In Section 2, we formulate the SWIUA model of HIV that describes the stated assumptions, variables, and parameters and in Section 3, a mathematical analysis of the model without control is carried out. Moreover, in Section 4, optimal control problem is discussed. Finally, in Section 5, numerical simulations, results and discussion, and the conclusion are presented.


## 2. Model Formulation

In this study, a deterministic mathematical model is formulated to analyze the transmission dynamics of the human immunodeficiency virus through applying prevention and treatment control strategies. The current model is a modification of the SIA model of HIV discussed in [13]. A base model has three compartments of human population and did not explicitly show the effectiveness of ART that reduces the viral load in the blood to the undetectable level for ART users. The present work is done to fill the gap observed in the previous study. The compartments of the current model are described as follows:
(i)
Susceptible compartment: this compartment is denoted byS and embraces all humans who are free of HIV but have a chance of being infected in the infective environment. All humans in this compartment transfer to the HIV compartment at the transmission rate β, provided that effective contact is done with humans in HIV compartment(ii)
HIV untested compartment: this compartment consists of all individuals who are at pre-AIDS stage and not test for HIV. They transmits virus to others at transmission rateβ1(iii)
HIV tested compartment: it is denoted byI. This compartment includes all pre-AIDS humans who are infected with HIV and tests positive at healthy center and they transmit virus to susceptible individuals at transmission rate β2(iv)
Undetectable compartment: it is denoted byU. This compartment includes all humans who get infected with human immunodeficiency virus, but with undetectable viral load, as a result of effective usage of antiretroviral treatment(v)
AIDS compartment: it is denoted byA. This compartment includes all humans who are at advanced stage of HIV and face loss of life at disease induced death rate δIn the formulated mathematical model of HIV, the following assumptions are stated:(i)
The total size of population is assumed to be nonconstant(ii)
It is assumed that HIV untested individuals get tested for HIV infection before they develop AIDS disease(iii)
The total population size at timet is denoted by Nt is given by(1)Nt=St+Wt+It+Ut+At(iv)
Susceptible humans are recruited to the compartmentSt at some constant rate τ(v)
Susceptible humans get HIV and join the HIV compartment at a constant rateβ(vi)
Individuals transfer from HIV compartment to undetectable compartment at the constant rateθ and transfer to AIDS compartment at a constant rate α(vii)
Individuals in undetectable compartment transfer to the HIV compartment at a constant rate∅(viii)
All categories of human’s compartments face the same natural mortality with a rateμ(ix)
All humans in AIDS compartment suffer disease-induced death at a constant rateδ(x)
u1 is prevention control effort(xi)
u2 is treatment control effort(xii)
All parameters used in the dynamical system are positiveMoreover, the notations and description of model variables are given in Table1, whereas model parameters notations and descriptions are given in Table 2.Table 1
Notations and description of model variables.
VariableDescriptionStSize of susceptible population at timetWtSize of HIV untested population at timetItSize of HIV tested pre-AIDS population with transmittable virus at timet.UtSize of pre-AIDS population with untransmittable virus at timetAtSize of AIDS population at timetTable 2
Model parameter notations and description.
ParameterDescriptionβTransmission rate of HIV to susceptible populationλRecruitment rate of individuals to susceptible classθTransferring rate of population from compartmentI to U as a result of using ART properlyϕProgression rate of population from compartmentU to I as a result of faulty using ART.ξProgression rate of population from compartmentW to IαProgression rate ofI to AκTransferring rate ofI to A due to treatmentμNatural death rate of all populationδDisease induced death rate of AIDS populationBased on flow diagram given in Figure1 and assumptions considered, without control measures, the dynamics of the populations are represented by the subsequent dynamical system:
(2)dSdt=λ−Sβ1W+β2IN−μS,dWdt=Sβ1W+β2IN−ξ+μW,dIdt=ξW+ϕU−η+μI,dUdt=−∅+μU,dAdt=ηI−δ+μA,Figure 1
Flow diagram of model in the study.with nonnegative initial conditions areS0>0,W0≥0,I0≥0,U0≥0,A0≥0.


## 3. Mathematical Analysis of the Model

### 3.1. Invariant Region

Theorem 1.
Ω⊂R+5 is a region in which all solutions of the model (2) are bounded provided that initial conditions are bounded. That is,
(3)Ω=St,Wt,It,Ut,At∈R+5:Nt≤λμ.Proof.
Consider a total population size,Nt, at time t given by
(4)Nt=St+Wt+It+Ut+At.
Now, differentiating both sides of equation (4) with respect to time t, we have
(5)dNdt=dSdt+dWdt+dIdt+dUdt+dAdt
Reduced to,(6)dNdt=λ−μN−δA≤λ−μN⇒dNdt≤λ−μ⇒dNτ−μN≤dt.
Integrating both sides of (6) and using comparison theorem [53], we have,
(7)Nt≤λμ−λμ−N0e−μt.
Here, from (7) it follows that as t⟶∞,Nt⟶λ/μ. That is, for all possibility, the expressions on the right hand side of the inequality (7) either increase to least upper bound or decrease to greatest lower bound as time increase. Hence, Nt≤λ/μ.Therefore, the feasible solution set of model (2) is the invariant region Ω, defined as
(8)Ω=St,Wt,It,Ut,At∈R+5:Nt≤λμ.


### 3.2. Positivity of Model Solutions

Theorem 2.
Solutions of the model (2)are always nonnegative for all t and will remain in ℝ+5.Proof.
The proof follows by showing that each solution variable is nonnegative. Considering the first equation of model (2) we have
(9)dSdt=λ−Sβ1W+β2IN−μS.
Ignoring the positive termλ from the preceding equation we get the following inequality.
(10)dSS≥−β1W+β2INdt.
Solving the forgoing inequality over time interval0,t, we get the following inequality.
(11)St≥S0e−μt−∫0tβ1Wξ+β2Iξ/Nξdξ.In the preceding inequality, we know thate−μt−∫0tβ1Wξ+β2Iξ/Nξdξ is nonnegative exponential expression, and S0 is positive from initial condition. Therefore, St is nonnegative for all time t. Similarly, the remaining solution variables Wt,It,Ut and At are nonnegative. Hence, the solutions of model (2) are nonnegative.Therefore, from the discussions in subsection it can be concluded that the formulated model is mathematically well-posed and solutions of the model are biologically meaningful.


### 3.3. Existence and Uniqueness of Solutions

Theorem 3.
The solutions of model (2) exists and unique.Proof.
Applying Theorem1, all expressions on the right hand of equality of model (2) are continuously differentiable and bounded. Therefore, by Cauchy-Lipchitz theorem we can conclude that the formulated model has unique solution for all positive time t.


### 3.4. Equilibria of Model

#### 3.4.1. Disease-Free Equilibrium (DFE)

The disease-free equilibrium point of model (2) is a steady state point where there is no HIV infection in the population. Hence, to obtain DFE the state variables I,A, and U are set to zero W=I=A=U=0. The computed disease–free equilibrium point E0 of the model (2) is given by
(12)E0=λμ,0,0,0,0.


#### 3.4.2. Endemic Equilibrium (EE)

An endemic equilibrium point is steady state point where disease persists in the population. The endemic equilibrium is obtained by setting the right-hand of model equation to zero and evaluating forS,W,I,U, and A in terms of model parameters. Moreover, adding equations (2) and (4), we obtain I in terms of S. Finally, using equation (4), setting N=S+W+I+U+A, and expressing all variables in terms of S, we obtain the desired solution. Thus, the computed endemic equilibrium E1=S∗,W∗,I∗,U∗,A∗ of model (2) is given by
(13)E1=λ1+d+bdβ1+β2d+μd+bμd−ξ,W∗,dW∗,0,bI∗,where
(14)W∗=λR0−1β1+β2d+μd+bμd−ξ,b=ηδ+μ,d=ξη+μ.


### 3.5. Basic Reproduction Number

The basic reproduction number is the average number of cases generated by one infected person in wholly susceptible population [45–50]. The basic reproductive number R0 can be determined using the next-generation matrix. In this method, R0 is defined as the largest eigenvalue of the next-generation matrix [51]. The formulation of this matrix involves classification of all compartments of the model in to two classes: infected and noninfected. Let f be a matrix of newly infected cases and v be a matrix of transition cases in model (2). Consider model (2).Nowf and v are given, respectively, as,
(15)f=Sβ1W+β2IN000,v=ξ+μW−ξW−ϕU+η+μI∅+μU−ηI+δ+μA,The Jacobian off and v evaluated at disease-free equilibrium point E0 is given by F and V, respectively, as follows,
(16)F=β1β200000000000000,V=ξ+μ000−ξη+μ−ϕ000∅+μ00−η0δ+μ.The next-generation matrix,FV−1, is computed as
(17)FV−1=β1β200000000000000ξ+μ000−ξη+μ−ϕ000∅+μ00−η0δ+μ−1=β1ξ+μ+β2ξξ+μη+μβ2η+μβ2ϕϕ+μη+μ0000000000000.The eigenvalues of next-generation matrix are computed and given by(18)λ1=β1ξ+μ+β2ξξ+μη+μ,λ2=0,λ3=0,λ4=0.Here,λ1 is the largest eigenvalue of next-generation matrix. Therefore, the reproduction number, R0, is given by
(19)R0=β1ξ+μ+β2ξξ+μη+μ.


### 3.6. Bifurcation Analysis

The Hopf bifurcation is studied in [5], and it describes the behavior of the system as system changes its steady state. In this subsection, we verify the possibility of backward and forward bifurcation using the center manifold theory stated in reference [28]. Now, model (2) can be written in the vector form, by renaming the variables as S=x1,W=x2,I=x3,U=x4,A=x5. That is,
(20)dXdt=FX,where, X=x1,x2x3,x4x5,T,FX=f1f2,f3,f4,f5T. Then, model (2) becomes
(21)dx1dt=λ−β1x1x2+β2x1x3N−μx1=f1,dx2dt=β1x1x2+β2x1x3N−ξ+μx2=f2,dx3dt=ξx2+ϕx4−η+μx3=f3,dx4dt=−∅+μx4=f4,dx5dt=ηx3−δ+μx5=f5.Here, from preceding system of nonlinear equation, choosingβ1 as a bifurcation parameter and setting R0=1, we have
(22)β1∗≡β1=ξ+μ−β2ξη+μ.So that the disease-free equilibrium,E0, is locally stable when β1<β1∗ and is unstable when β1>β1∗. The linearized matrix of the system around the disease-free equilibrium E0 and evaluated at β1∗ is given by
(23)JE0,β1∗=−μ−ξ+μ+β2ξη+μ−β2000−β2ξη+μβ2000ξ−η+μϕ0000−∅+μ000η0−δ+μ.The eigenvalues of matrixJE0,β1∗ are λ1=0,λ2=−μ,λ3=−δ−μ. The sign of remaining eigenvalues are determined from characteristic polynomial:
(24)fλ=a0λ2+a1λ+a2,where
(25)a0=1,a1=ϕ2η+ϕ2μ+ϕη2+4ϕημ+3ϕμ2+β2ξϕ+η2μ+3ημ2+2μ3+β2ξμϕη+ϕμ+ημ+μ2,a2=ϕ2η2+β2ξϕ2+2ϕμϕη+η2+2ημ+μ2+β2ξ+ϕ2+η2+2ημ+μ2+β2ξμ2ϕη+ϕμ+ημ+μ2.We observe that coefficients,a0,a1, and a2, are positive, and all eigenvalues of Jacobian matrix are negative while one eigenvalue is zero. Moreover, let w=w1,w2,w3,w4,w5T be the right eigenvector, associated with simple zero eigenvalue and can be obtained by solving equation JE0,β∗w=0. Thus, we have
(26)w1=−ξ+μη+μϕ+μ+θμξμϕ+μw3,w2=η+μξw3,w3=w3,w4=θw3ϕ+μ,w5=ηw3δ+μ.Also, letv=v1,v2,v3,v4,v5 be the left eigenvector corresponding to simple zero eigenvalue, obtained by setting and solving vJE0,β∗=0. This gives, v1=0,v2=η+μ/β2v3,v3=v3,v4=ϕv3/ϕ+μ,v5=0. In order that the required conditions v.w=1, the product gives,
(27)v1w1+v2w2+v3w3+v4w4+v5w5=1.Again substituting the corresponding components in the preceding equation, we have(28)v3w3η+μ2+β2ξβ2ξ=1.Thus, the preceding equality is satisfied if we choose(29)v3=β2ξ,w3=1η+μ2+β2ξ.Next, we compute all second-order partial derivatives of functions on the right hand of (∗) as given in the appendix.Next, we compute bifurcation coefficientsa and b:
(30)a=∑k,i,j=14vkwiwj∂2fk∂xi∂xjE0,β∗=∑i,j=14v2wiwj∂2f2∂xi∂xjE0,β∗=−2β1x1−β1x1−β2x1−β1x1−2β2x1−β2x1<0,b=∑k,i=14vkwi∂2fk∂xi∂β∗E0,β∗=v2w2=η+μ2η+μ2+β2ρ>0.Sincea<0 and b>0, the model exhibits forward bifurcation at R0=1.Next, the following procedures given in (Biswas et al., 2020), we compute the bifurcation coefficientsa and b, to identify the direction of bifurcation at R0=1. Thus, we have
(31)a=∑k,i,j=14vkwiwj∂2fk∂xi∂xjE0,β∗=−2v2μ3w22β∗τ3+2v2w2μβ∗τϕ+μ2,b=∑k,i=14vkwi∂2fk∂xi∂β∗E0,β∗=v2w2.Since all parameters in model (2) are nonnegative and additionally v2 and w2 are positive, we conclude a<0 and b>0. Thus, according to [5], model (2) exhibits a supercritical (forward) bifurcation, when R0 crosses the threshold R0=1. That is, there exist locally asymptotically stable endemic equilibrium point E1=S∗,W∗,I∗,U∗,A∗ for R0>1. Based on the results of the above discussion and [7], the following theorem is stated.Theorem 4.
Thetrans-critical bifurcation of model (2) that occurs at R0=1 is a supercritical (forward) bifurcation. That is, there exists locally asymptotically stable endemic equilibrium point E1=S∗,W∗,I∗,U∗,A∗ for R0>1.Remark 5.
Theorem4 shows that if R0>1, then a few inflow of infectious individuals in fully susceptible population can result in persistence of the HIV in the population.


### 3.7. Stability Analysis of Equilibrium Points

In absence of the infectious disease, model (2) has a unique disease-free steady-state E0, and in the presence of disease, model (2) has unique endemic equilibrium E1.


#### 3.7.1. Local Stability of Disease-Free Equilibrium

Theorem 6.
The DFEE0 of the model (2) is locally asymptotically stable if R0<1 and unstable if R0>1.Proof.
Consider model (2), so that the Jacobian matrix of the system at DFE is given by
(32)J=−μ−β1−β2000β1−ξ+μβ2000ξ−η+μϕ0000−∅+μ000η0−δ+μ.
The eigenvalues of a preceding Jacobian matrixJ is computed and given by
(33)λ1=−μ,λ2=−δ+μ,λ3=−∅+μ,λ4=β1−ξ−η−2μ−β1−ξ2+2ηβ1−2ξη+4β2ξ+η22,λ5=β1−ξ−η−2μ+β1−ξ2+2ηβ1−2ξη+4β2ξ+η22.
Clearly, the first three eigenvalues are negative whereas the fourth and fifth eigenvalues are negative ifR0<1.Further, by a stability analysis of a point using Jacobian matrix, we conclude that the disease-free equilibrium point is locally asymptotically stable ifR0<1 and unstable otherwise.


#### 3.7.2. Local Stability of Endemic Equilibrium Point

Theorem 7.
The EEE1 of the model (2) is locally asymptotically stable if R0>1 and unstable if R0<1. (Proof. Behavior of forward bifurcation.)


#### 3.7.3. Global Stability of Disease–Free Equilibrium Point

To show global stability of disease-free equilibriumE0, we use the technique employed in [4]. Accordingly, let X∈ℝ1 denote individuals in uninfected compartment S and Y∈ℝ4 denotes individuals in infected compartments W,I,U,A. Hence, we write model (2) in the form:
(34)dXdt=HX,Y,dYdt=GX,Y,GX,0=0.Also, the disease-free equilibrium is given by(35)E01=X0,0.Here,X0 is the disease-free equilibrium of the foregoing system.To guarantee global asymptotic stability of the disease-free equilibrium point, the technique we employed must met the following two conditions H1 and H2.H1: fordX/dt=HX,0,X0 is globally asymptotically stableH2:GX,Y=PY−G^X,Y,G^X,Y≥0 for X,Y∈ΩHere,P=DYGX,0 is a Metzler matrix and Ω is a region where solutions of the model are biologically feasible.Theorem 8.
The disease-free equilibrium pointE0 of model (2) is globally asymptotically stable in a region Ω if R0<1 as the fourth and fifth eigenvalues are negative if R0<1 and unstable whenever R0>1 provided that the above two conditions H1 and H2 are satisfied, where Ω is a feasible solution region of model (2) in ℝ+5.Proof.
From model (2), we have
(36)HX,0=λ−μS=HS,0.
PuttingHX,0=0 and solving, we obtain S=λ/μ. Hence, X0=λ/μ,0. Clearly, X0 is globally asymptotically stable equilibrium point of equation:
(37)dXdt=HX,0.From infected compartments of model (2), we have
(38)GX,Y=β1SN−ξ+μβ2SN00ξ−η+μϕ000−∅+μ00η0−δ+μWIUA=β1SNβ2SN00ξ−η+μϕ000−∅+μ00η0−δ+μWIUA−ξ+μ000000000000000WIUA.In the preceding computation, let(39)P=β1SNβ2SN00ξ−η+μϕ000−∅+μ00η0−δ+μ,so that at disease-free equilibrium it reduces to(40)P=β1−ξ+μβ200ξ−η+μϕ000−∅+μ00η0−δ+μ.Again,PY gives
(41)PY=β1−ξ+μβ200ξ−η+μϕ000−∅+μ00η0−δ+μWIUA.Let(42)G^X,Y=PY−G=β1−ξ+μβ200ξ−η+μϕ000−∅+μ00η0−δ+μWIUA−β1SN−ξ+μβ2SN00ξ−η+μϕ000−∅+μ00η0−δ+μWIUA.Therefore,(43)G^X,Y=PY−G=β1W+β2I−Sβ1W+β2IN000.This implies that(44)G^X,Y=β1W+β2I1−SS+W+I+U+A000.SinceS≤S+W+I+U+A, we have β1W+β2I1−S/S+W+I+U+A≥0; the preceding matrix is a nonnegative matrix. Therefore, it can be written as
(45)G^X,Y=β1W+β2I1−SS+W+I+U+A000≥0.Clearly, the second condition is satisfied asGX,Y=PY−G^X,Y, where, G^X,Y≥0,∀X,Y in the invariant region. Comparing the computations, we have
(46)P=β1−ξ+μβ200ξ−η+μϕ000−∅+μ00η0−δ+μ=DYGX,0,which is a Metzler matrix as all off diagonal elements of a matrix are nonnegative.Therefore, the disease-free equilibriumE0 is globally asymptotically stable if R0<1 as the fourth and fifth eigenvalues are negative if R0<1.


### 3.8. Sensitivity Analysis

Sensitivity analysis is used to determine the sensitivity of the variable with respect to the parameters involved in it [8, 13, 28]. The normalized forward sensitivity index of a particular variable, R, with respect to a parameter, p, is defined as
(47)YpR=∂R∂p×pR.It is already shown that the explicit expression ofR0 is given by
(48)R0=β1ξ+μ+β2ξξ+μη+μ.The normalized forward sensitivity indices ofR0 with respect to parameters in it are given by
(49)Yβ1R0=∂R0∂β1β1R0=1ξ+μβ1R0,Yβ2R0=∂R0∂β2β2R0=ξξ+μη+μβ2R0,YξR0=∂R0∂ξ×ξR0=−β1ξ+μ2+β2μξ+μ2η+μξR0,YμR0=∂R0∂μ×μR0=−β1ξ+μ2+−β2ξη+ξ+2μξ+μ2η+μ2μR0,YηR0=∂R0∂η×ηR0=−β2ξξ+μη+μ2ηR0,where, based on parametric values given in Table 3, we obtain
(50)R0=8.23.Table 3
Sensitivity index value ofR0 with respect parameter.
ParameterValueSensitivity index valueDescriptionβ10.9815+0.145β1~R0β20.866+0.854β2~R0ξ0.8-0.121ξ~1R0μ0.02-0.004μ~1R0η0.1-0.701η~1R0From Table3, it can be observed that parameters β1 and β2 have positive sensitivity indices and the values of the remaining two parameters ξ, η, and μ get negative sensitivity indices. Further, these parameters affect the value of a reproduction number that helps in the analysis of virus extinction or persistence of the disease in the population. On the other hand, an increase of positive parameter value will increase the value of R0; this implies that disease transmission increase with significant amount. Also, an increase in the magnitude of negative parameter value will cause a value of reproduction number to decrease in some amount, which means the disease transmission significantly decreases in some amount.


## 3.1. Invariant Region

Theorem 1.
Ω⊂R+5 is a region in which all solutions of the model (2) are bounded provided that initial conditions are bounded. That is,
(3)Ω=St,Wt,It,Ut,At∈R+5:Nt≤λμ.Proof.
Consider a total population size,Nt, at time t given by
(4)Nt=St+Wt+It+Ut+At.
Now, differentiating both sides of equation (4) with respect to time t, we have
(5)dNdt=dSdt+dWdt+dIdt+dUdt+dAdt
Reduced to,(6)dNdt=λ−μN−δA≤λ−μN⇒dNdt≤λ−μ⇒dNτ−μN≤dt.
Integrating both sides of (6) and using comparison theorem [53], we have,
(7)Nt≤λμ−λμ−N0e−μt.
Here, from (7) it follows that as t⟶∞,Nt⟶λ/μ. That is, for all possibility, the expressions on the right hand side of the inequality (7) either increase to least upper bound or decrease to greatest lower bound as time increase. Hence, Nt≤λ/μ.Therefore, the feasible solution set of model (2) is the invariant region Ω, defined as
(8)Ω=St,Wt,It,Ut,At∈R+5:Nt≤λμ.


## 3.2. Positivity of Model Solutions

Theorem 2.
Solutions of the model (2)are always nonnegative for all t and will remain in ℝ+5.Proof.
The proof follows by showing that each solution variable is nonnegative. Considering the first equation of model (2) we have
(9)dSdt=λ−Sβ1W+β2IN−μS.
Ignoring the positive termλ from the preceding equation we get the following inequality.
(10)dSS≥−β1W+β2INdt.
Solving the forgoing inequality over time interval0,t, we get the following inequality.
(11)St≥S0e−μt−∫0tβ1Wξ+β2Iξ/Nξdξ.In the preceding inequality, we know thate−μt−∫0tβ1Wξ+β2Iξ/Nξdξ is nonnegative exponential expression, and S0 is positive from initial condition. Therefore, St is nonnegative for all time t. Similarly, the remaining solution variables Wt,It,Ut and At are nonnegative. Hence, the solutions of model (2) are nonnegative.Therefore, from the discussions in subsection it can be concluded that the formulated model is mathematically well-posed and solutions of the model are biologically meaningful.


## 3.3. Existence and Uniqueness of Solutions

Theorem 3.
The solutions of model (2) exists and unique.Proof.
Applying Theorem1, all expressions on the right hand of equality of model (2) are continuously differentiable and bounded. Therefore, by Cauchy-Lipchitz theorem we can conclude that the formulated model has unique solution for all positive time t.


## 3.4. Equilibria of Model

### 3.4.1. Disease-Free Equilibrium (DFE)

The disease-free equilibrium point of model (2) is a steady state point where there is no HIV infection in the population. Hence, to obtain DFE the state variables I,A, and U are set to zero W=I=A=U=0. The computed disease–free equilibrium point E0 of the model (2) is given by
(12)E0=λμ,0,0,0,0.


### 3.4.2. Endemic Equilibrium (EE)

An endemic equilibrium point is steady state point where disease persists in the population. The endemic equilibrium is obtained by setting the right-hand of model equation to zero and evaluating forS,W,I,U, and A in terms of model parameters. Moreover, adding equations (2) and (4), we obtain I in terms of S. Finally, using equation (4), setting N=S+W+I+U+A, and expressing all variables in terms of S, we obtain the desired solution. Thus, the computed endemic equilibrium E1=S∗,W∗,I∗,U∗,A∗ of model (2) is given by
(13)E1=λ1+d+bdβ1+β2d+μd+bμd−ξ,W∗,dW∗,0,bI∗,where
(14)W∗=λR0−1β1+β2d+μd+bμd−ξ,b=ηδ+μ,d=ξη+μ.


## 3.4.1. Disease-Free Equilibrium (DFE)

The disease-free equilibrium point of model (2) is a steady state point where there is no HIV infection in the population. Hence, to obtain DFE the state variables I,A, and U are set to zero W=I=A=U=0. The computed disease–free equilibrium point E0 of the model (2) is given by
(12)E0=λμ,0,0,0,0.


## 3.4.2. Endemic Equilibrium (EE)

An endemic equilibrium point is steady state point where disease persists in the population. The endemic equilibrium is obtained by setting the right-hand of model equation to zero and evaluating forS,W,I,U, and A in terms of model parameters. Moreover, adding equations (2) and (4), we obtain I in terms of S. Finally, using equation (4), setting N=S+W+I+U+A, and expressing all variables in terms of S, we obtain the desired solution. Thus, the computed endemic equilibrium E1=S∗,W∗,I∗,U∗,A∗ of model (2) is given by
(13)E1=λ1+d+bdβ1+β2d+μd+bμd−ξ,W∗,dW∗,0,bI∗,where
(14)W∗=λR0−1β1+β2d+μd+bμd−ξ,b=ηδ+μ,d=ξη+μ.


## 3.5. Basic Reproduction Number

The basic reproduction number is the average number of cases generated by one infected person in wholly susceptible population [45–50]. The basic reproductive number R0 can be determined using the next-generation matrix. In this method, R0 is defined as the largest eigenvalue of the next-generation matrix [51]. The formulation of this matrix involves classification of all compartments of the model in to two classes: infected and noninfected. Let f be a matrix of newly infected cases and v be a matrix of transition cases in model (2). Consider model (2).Nowf and v are given, respectively, as,
(15)f=Sβ1W+β2IN000,v=ξ+μW−ξW−ϕU+η+μI∅+μU−ηI+δ+μA,The Jacobian off and v evaluated at disease-free equilibrium point E0 is given by F and V, respectively, as follows,
(16)F=β1β200000000000000,V=ξ+μ000−ξη+μ−ϕ000∅+μ00−η0δ+μ.The next-generation matrix,FV−1, is computed as
(17)FV−1=β1β200000000000000ξ+μ000−ξη+μ−ϕ000∅+μ00−η0δ+μ−1=β1ξ+μ+β2ξξ+μη+μβ2η+μβ2ϕϕ+μη+μ0000000000000.The eigenvalues of next-generation matrix are computed and given by(18)λ1=β1ξ+μ+β2ξξ+μη+μ,λ2=0,λ3=0,λ4=0.Here,λ1 is the largest eigenvalue of next-generation matrix. Therefore, the reproduction number, R0, is given by
(19)R0=β1ξ+μ+β2ξξ+μη+μ.


## 3.6. Bifurcation Analysis

The Hopf bifurcation is studied in [5], and it describes the behavior of the system as system changes its steady state. In this subsection, we verify the possibility of backward and forward bifurcation using the center manifold theory stated in reference [28]. Now, model (2) can be written in the vector form, by renaming the variables as S=x1,W=x2,I=x3,U=x4,A=x5. That is,
(20)dXdt=FX,where, X=x1,x2x3,x4x5,T,FX=f1f2,f3,f4,f5T. Then, model (2) becomes
(21)dx1dt=λ−β1x1x2+β2x1x3N−μx1=f1,dx2dt=β1x1x2+β2x1x3N−ξ+μx2=f2,dx3dt=ξx2+ϕx4−η+μx3=f3,dx4dt=−∅+μx4=f4,dx5dt=ηx3−δ+μx5=f5.Here, from preceding system of nonlinear equation, choosingβ1 as a bifurcation parameter and setting R0=1, we have
(22)β1∗≡β1=ξ+μ−β2ξη+μ.So that the disease-free equilibrium,E0, is locally stable when β1<β1∗ and is unstable when β1>β1∗. The linearized matrix of the system around the disease-free equilibrium E0 and evaluated at β1∗ is given by
(23)JE0,β1∗=−μ−ξ+μ+β2ξη+μ−β2000−β2ξη+μβ2000ξ−η+μϕ0000−∅+μ000η0−δ+μ.The eigenvalues of matrixJE0,β1∗ are λ1=0,λ2=−μ,λ3=−δ−μ. The sign of remaining eigenvalues are determined from characteristic polynomial:
(24)fλ=a0λ2+a1λ+a2,where
(25)a0=1,a1=ϕ2η+ϕ2μ+ϕη2+4ϕημ+3ϕμ2+β2ξϕ+η2μ+3ημ2+2μ3+β2ξμϕη+ϕμ+ημ+μ2,a2=ϕ2η2+β2ξϕ2+2ϕμϕη+η2+2ημ+μ2+β2ξ+ϕ2+η2+2ημ+μ2+β2ξμ2ϕη+ϕμ+ημ+μ2.We observe that coefficients,a0,a1, and a2, are positive, and all eigenvalues of Jacobian matrix are negative while one eigenvalue is zero. Moreover, let w=w1,w2,w3,w4,w5T be the right eigenvector, associated with simple zero eigenvalue and can be obtained by solving equation JE0,β∗w=0. Thus, we have
(26)w1=−ξ+μη+μϕ+μ+θμξμϕ+μw3,w2=η+μξw3,w3=w3,w4=θw3ϕ+μ,w5=ηw3δ+μ.Also, letv=v1,v2,v3,v4,v5 be the left eigenvector corresponding to simple zero eigenvalue, obtained by setting and solving vJE0,β∗=0. This gives, v1=0,v2=η+μ/β2v3,v3=v3,v4=ϕv3/ϕ+μ,v5=0. In order that the required conditions v.w=1, the product gives,
(27)v1w1+v2w2+v3w3+v4w4+v5w5=1.Again substituting the corresponding components in the preceding equation, we have(28)v3w3η+μ2+β2ξβ2ξ=1.Thus, the preceding equality is satisfied if we choose(29)v3=β2ξ,w3=1η+μ2+β2ξ.Next, we compute all second-order partial derivatives of functions on the right hand of (∗) as given in the appendix.Next, we compute bifurcation coefficientsa and b:
(30)a=∑k,i,j=14vkwiwj∂2fk∂xi∂xjE0,β∗=∑i,j=14v2wiwj∂2f2∂xi∂xjE0,β∗=−2β1x1−β1x1−β2x1−β1x1−2β2x1−β2x1<0,b=∑k,i=14vkwi∂2fk∂xi∂β∗E0,β∗=v2w2=η+μ2η+μ2+β2ρ>0.Sincea<0 and b>0, the model exhibits forward bifurcation at R0=1.Next, the following procedures given in (Biswas et al., 2020), we compute the bifurcation coefficientsa and b, to identify the direction of bifurcation at R0=1. Thus, we have
(31)a=∑k,i,j=14vkwiwj∂2fk∂xi∂xjE0,β∗=−2v2μ3w22β∗τ3+2v2w2μβ∗τϕ+μ2,b=∑k,i=14vkwi∂2fk∂xi∂β∗E0,β∗=v2w2.Since all parameters in model (2) are nonnegative and additionally v2 and w2 are positive, we conclude a<0 and b>0. Thus, according to [5], model (2) exhibits a supercritical (forward) bifurcation, when R0 crosses the threshold R0=1. That is, there exist locally asymptotically stable endemic equilibrium point E1=S∗,W∗,I∗,U∗,A∗ for R0>1. Based on the results of the above discussion and [7], the following theorem is stated.Theorem 4.
Thetrans-critical bifurcation of model (2) that occurs at R0=1 is a supercritical (forward) bifurcation. That is, there exists locally asymptotically stable endemic equilibrium point E1=S∗,W∗,I∗,U∗,A∗ for R0>1.Remark 5.
Theorem4 shows that if R0>1, then a few inflow of infectious individuals in fully susceptible population can result in persistence of the HIV in the population.


## 3.7. Stability Analysis of Equilibrium Points

In absence of the infectious disease, model (2) has a unique disease-free steady-state E0, and in the presence of disease, model (2) has unique endemic equilibrium E1.


### 3.7.1. Local Stability of Disease-Free Equilibrium

Theorem 6.
The DFEE0 of the model (2) is locally asymptotically stable if R0<1 and unstable if R0>1.Proof.
Consider model (2), so that the Jacobian matrix of the system at DFE is given by
(32)J=−μ−β1−β2000β1−ξ+μβ2000ξ−η+μϕ0000−∅+μ000η0−δ+μ.
The eigenvalues of a preceding Jacobian matrixJ is computed and given by
(33)λ1=−μ,λ2=−δ+μ,λ3=−∅+μ,λ4=β1−ξ−η−2μ−β1−ξ2+2ηβ1−2ξη+4β2ξ+η22,λ5=β1−ξ−η−2μ+β1−ξ2+2ηβ1−2ξη+4β2ξ+η22.
Clearly, the first three eigenvalues are negative whereas the fourth and fifth eigenvalues are negative ifR0<1.Further, by a stability analysis of a point using Jacobian matrix, we conclude that the disease-free equilibrium point is locally asymptotically stable ifR0<1 and unstable otherwise.


### 3.7.2. Local Stability of Endemic Equilibrium Point

Theorem 7.
The EEE1 of the model (2) is locally asymptotically stable if R0>1 and unstable if R0<1. (Proof. Behavior of forward bifurcation.)


### 3.7.3. Global Stability of Disease–Free Equilibrium Point

To show global stability of disease-free equilibriumE0, we use the technique employed in [4]. Accordingly, let X∈ℝ1 denote individuals in uninfected compartment S and Y∈ℝ4 denotes individuals in infected compartments W,I,U,A. Hence, we write model (2) in the form:
(34)dXdt=HX,Y,dYdt=GX,Y,GX,0=0.Also, the disease-free equilibrium is given by(35)E01=X0,0.Here,X0 is the disease-free equilibrium of the foregoing system.To guarantee global asymptotic stability of the disease-free equilibrium point, the technique we employed must met the following two conditions H1 and H2.H1: fordX/dt=HX,0,X0 is globally asymptotically stableH2:GX,Y=PY−G^X,Y,G^X,Y≥0 for X,Y∈ΩHere,P=DYGX,0 is a Metzler matrix and Ω is a region where solutions of the model are biologically feasible.Theorem 8.
The disease-free equilibrium pointE0 of model (2) is globally asymptotically stable in a region Ω if R0<1 as the fourth and fifth eigenvalues are negative if R0<1 and unstable whenever R0>1 provided that the above two conditions H1 and H2 are satisfied, where Ω is a feasible solution region of model (2) in ℝ+5.Proof.
From model (2), we have
(36)HX,0=λ−μS=HS,0.
PuttingHX,0=0 and solving, we obtain S=λ/μ. Hence, X0=λ/μ,0. Clearly, X0 is globally asymptotically stable equilibrium point of equation:
(37)dXdt=HX,0.From infected compartments of model (2), we have
(38)GX,Y=β1SN−ξ+μβ2SN00ξ−η+μϕ000−∅+μ00η0−δ+μWIUA=β1SNβ2SN00ξ−η+μϕ000−∅+μ00η0−δ+μWIUA−ξ+μ000000000000000WIUA.In the preceding computation, let(39)P=β1SNβ2SN00ξ−η+μϕ000−∅+μ00η0−δ+μ,so that at disease-free equilibrium it reduces to(40)P=β1−ξ+μβ200ξ−η+μϕ000−∅+μ00η0−δ+μ.Again,PY gives
(41)PY=β1−ξ+μβ200ξ−η+μϕ000−∅+μ00η0−δ+μWIUA.Let(42)G^X,Y=PY−G=β1−ξ+μβ200ξ−η+μϕ000−∅+μ00η0−δ+μWIUA−β1SN−ξ+μβ2SN00ξ−η+μϕ000−∅+μ00η0−δ+μWIUA.Therefore,(43)G^X,Y=PY−G=β1W+β2I−Sβ1W+β2IN000.This implies that(44)G^X,Y=β1W+β2I1−SS+W+I+U+A000.SinceS≤S+W+I+U+A, we have β1W+β2I1−S/S+W+I+U+A≥0; the preceding matrix is a nonnegative matrix. Therefore, it can be written as
(45)G^X,Y=β1W+β2I1−SS+W+I+U+A000≥0.Clearly, the second condition is satisfied asGX,Y=PY−G^X,Y, where, G^X,Y≥0,∀X,Y in the invariant region. Comparing the computations, we have
(46)P=β1−ξ+μβ200ξ−η+μϕ000−∅+μ00η0−δ+μ=DYGX,0,which is a Metzler matrix as all off diagonal elements of a matrix are nonnegative.Therefore, the disease-free equilibriumE0 is globally asymptotically stable if R0<1 as the fourth and fifth eigenvalues are negative if R0<1.


## 3.7.1. Local Stability of Disease-Free Equilibrium

Theorem 6.
The DFEE0 of the model (2) is locally asymptotically stable if R0<1 and unstable if R0>1.Proof.
Consider model (2), so that the Jacobian matrix of the system at DFE is given by
(32)J=−μ−β1−β2000β1−ξ+μβ2000ξ−η+μϕ0000−∅+μ000η0−δ+μ.
The eigenvalues of a preceding Jacobian matrixJ is computed and given by
(33)λ1=−μ,λ2=−δ+μ,λ3=−∅+μ,λ4=β1−ξ−η−2μ−β1−ξ2+2ηβ1−2ξη+4β2ξ+η22,λ5=β1−ξ−η−2μ+β1−ξ2+2ηβ1−2ξη+4β2ξ+η22.
Clearly, the first three eigenvalues are negative whereas the fourth and fifth eigenvalues are negative ifR0<1.Further, by a stability analysis of a point using Jacobian matrix, we conclude that the disease-free equilibrium point is locally asymptotically stable ifR0<1 and unstable otherwise.


## 3.7.2. Local Stability of Endemic Equilibrium Point

Theorem 7.
The EEE1 of the model (2) is locally asymptotically stable if R0>1 and unstable if R0<1. (Proof. Behavior of forward bifurcation.)


## 3.7.3. Global Stability of Disease–Free Equilibrium Point

To show global stability of disease-free equilibriumE0, we use the technique employed in [4]. Accordingly, let X∈ℝ1 denote individuals in uninfected compartment S and Y∈ℝ4 denotes individuals in infected compartments W,I,U,A. Hence, we write model (2) in the form:
(34)dXdt=HX,Y,dYdt=GX,Y,GX,0=0.Also, the disease-free equilibrium is given by(35)E01=X0,0.Here,X0 is the disease-free equilibrium of the foregoing system.To guarantee global asymptotic stability of the disease-free equilibrium point, the technique we employed must met the following two conditions H1 and H2.H1: fordX/dt=HX,0,X0 is globally asymptotically stableH2:GX,Y=PY−G^X,Y,G^X,Y≥0 for X,Y∈ΩHere,P=DYGX,0 is a Metzler matrix and Ω is a region where solutions of the model are biologically feasible.Theorem 8.
The disease-free equilibrium pointE0 of model (2) is globally asymptotically stable in a region Ω if R0<1 as the fourth and fifth eigenvalues are negative if R0<1 and unstable whenever R0>1 provided that the above two conditions H1 and H2 are satisfied, where Ω is a feasible solution region of model (2) in ℝ+5.Proof.
From model (2), we have
(36)HX,0=λ−μS=HS,0.
PuttingHX,0=0 and solving, we obtain S=λ/μ. Hence, X0=λ/μ,0. Clearly, X0 is globally asymptotically stable equilibrium point of equation:
(37)dXdt=HX,0.From infected compartments of model (2), we have
(38)GX,Y=β1SN−ξ+μβ2SN00ξ−η+μϕ000−∅+μ00η0−δ+μWIUA=β1SNβ2SN00ξ−η+μϕ000−∅+μ00η0−δ+μWIUA−ξ+μ000000000000000WIUA.In the preceding computation, let(39)P=β1SNβ2SN00ξ−η+μϕ000−∅+μ00η0−δ+μ,so that at disease-free equilibrium it reduces to(40)P=β1−ξ+μβ200ξ−η+μϕ000−∅+μ00η0−δ+μ.Again,PY gives
(41)PY=β1−ξ+μβ200ξ−η+μϕ000−∅+μ00η0−δ+μWIUA.Let(42)G^X,Y=PY−G=β1−ξ+μβ200ξ−η+μϕ000−∅+μ00η0−δ+μWIUA−β1SN−ξ+μβ2SN00ξ−η+μϕ000−∅+μ00η0−δ+μWIUA.Therefore,(43)G^X,Y=PY−G=β1W+β2I−Sβ1W+β2IN000.This implies that(44)G^X,Y=β1W+β2I1−SS+W+I+U+A000.SinceS≤S+W+I+U+A, we have β1W+β2I1−S/S+W+I+U+A≥0; the preceding matrix is a nonnegative matrix. Therefore, it can be written as
(45)G^X,Y=β1W+β2I1−SS+W+I+U+A000≥0.Clearly, the second condition is satisfied asGX,Y=PY−G^X,Y, where, G^X,Y≥0,∀X,Y in the invariant region. Comparing the computations, we have
(46)P=β1−ξ+μβ200ξ−η+μϕ000−∅+μ00η0−δ+μ=DYGX,0,which is a Metzler matrix as all off diagonal elements of a matrix are nonnegative.Therefore, the disease-free equilibriumE0 is globally asymptotically stable if R0<1 as the fourth and fifth eigenvalues are negative if R0<1.


## 3.8. Sensitivity Analysis

Sensitivity analysis is used to determine the sensitivity of the variable with respect to the parameters involved in it [8, 13, 28]. The normalized forward sensitivity index of a particular variable, R, with respect to a parameter, p, is defined as
(47)YpR=∂R∂p×pR.It is already shown that the explicit expression ofR0 is given by
(48)R0=β1ξ+μ+β2ξξ+μη+μ.The normalized forward sensitivity indices ofR0 with respect to parameters in it are given by
(49)Yβ1R0=∂R0∂β1β1R0=1ξ+μβ1R0,Yβ2R0=∂R0∂β2β2R0=ξξ+μη+μβ2R0,YξR0=∂R0∂ξ×ξR0=−β1ξ+μ2+β2μξ+μ2η+μξR0,YμR0=∂R0∂μ×μR0=−β1ξ+μ2+−β2ξη+ξ+2μξ+μ2η+μ2μR0,YηR0=∂R0∂η×ηR0=−β2ξξ+μη+μ2ηR0,where, based on parametric values given in Table 3, we obtain
(50)R0=8.23.Table 3
Sensitivity index value ofR0 with respect parameter.
ParameterValueSensitivity index valueDescriptionβ10.9815+0.145β1~R0β20.866+0.854β2~R0ξ0.8-0.121ξ~1R0μ0.02-0.004μ~1R0η0.1-0.701η~1R0From Table3, it can be observed that parameters β1 and β2 have positive sensitivity indices and the values of the remaining two parameters ξ, η, and μ get negative sensitivity indices. Further, these parameters affect the value of a reproduction number that helps in the analysis of virus extinction or persistence of the disease in the population. On the other hand, an increase of positive parameter value will increase the value of R0; this implies that disease transmission increase with significant amount. Also, an increase in the magnitude of negative parameter value will cause a value of reproduction number to decrease in some amount, which means the disease transmission significantly decreases in some amount.


## 4. Optimal Control Problem

In this section, we extend SWIUA model by considering control strategies. The optimal control analysis assists to identify the best control strategies to eradicate or control the disease in the community at a specified period of time [3, 15, 19]. In the optimal control problem, the following three control measures are used.
(i)
Preventive control measureu1 that protects susceptible population from getting the disease(ii)
Treatment control measureu2 that is used by patients to reduce viral load in the body and slow/stop the progression of the virusNow, including these control measuresu1and u2 in model (2), we get the following optimal control model:
(51)dSdt=λ−1−u1Sβ1W+β2IN−μS,dWdt=1−u1Sβ1W+β2IN−ξ+μW,dIdt=ξW+u2κA+ϕU−1−u2η+u2θ+μI,dUdt=u2θI−∅+μU,dAdt=1−u2ηI−u2kA+δ+μA,with nonnegative initial conditionsTo study optimal levels of the controls, we define the Lebesgue measurable control setU as follows:
(52)U=u1t,u2t:0≤u1≤1,0≤u2≤1,0≤t≤tf.Our goal is to find the optimal controlsu1∗, and u2∗ optimal solutions S∗,W∗,I∗,U∗, and A∗ by fixing the terminal time tf that minimize the objective functional J given by
(53)Ju1,u2,u3=minu1,u2∫0tfc1W+c2I+c3A+12w1u12+w2u22dt,wherec1,c2,w1, and w2 are constants. The expressions 0.5w1u12,0.5w2u22, and 0.5w3u32 are costs associated with controls. The form of cost is quadratic because we assumed it to be nonlinear in nature [3]. Our goal is to minimize the population size of Wt,It, and At through intervention with control measures u1t, and u2t along with the costs associated with them. For two optimal controls u1∗, and u2∗, we have
(54)Ju1∗,u2∗,u3∗=minJu1,u2:u1,u2∈U,where U=u1,u2:0≤u1≤1,0≤u2≤1, u1,u2, and u3 are measurable controls.


### 4.1. Existence of the Optimal Control

To show the existence of optimal control, we use the approach used in [3]. We have already proved that the HIV model (2) is bounded, so this result can be used to prove the existence of optimal control over finite time interval as applied in [3, 52]. To ensure the existence of optimal control, we need to check if the following conditions are satisfied:
(i)
The set of controls and state variables be nonempty(ii)
The control setU is convex and closed(iii)
The right hand side of the state system is bounded by a linear function in the state and control variables(iv)
The integrand of objective functional is convex onU(v)
The integrand of objective functional is bounded below byk2−k1u12+u22k/2.k1,k2>0 and k>1An existence of the state system with bounded coefficients used to give condition (i). The control set is convex and closed by definition. The right hand side of the state system satisfies (iii). The state solutions are already bounded (iv). The integrand in the objective functionalc1W+c2I+c3A+1/2w1u12+w2u22 is clearly convex on U. (v) Further, from restriction on control measures, we have
(55)12wiui2≤12wi,ui∈0,1.Also, considering the preceding inequality, the integrand can be written as(56)c1W+c2I+c3A+12w1u12+w2u22≥k1u12+u22k/2−k2.wherek1=minw1/2,w2/2,k2=w2/2,k=2.Therefore, there exists optimal control measuresu1 and u2 that minimize the objective functional Ju1,u2,.


### 4.2. The Hamiltonian and Optimality System

The Pontryagin maximum principle stated the necessary conditions which are satisfied by optimal pair. Hence, by this principle, we obtained the Hamiltonian functionH defined as
(57)HS,W,I,U,A,=c1W+c2I+c3A+12w1u12+w2u22+λ1dSdt+λ2dWdt+λ3dIdt+λ4dUdt+λ5dAdt,where, λi,i=1,2,3,4,5 are the adjoint variables corresponding to state variables S,W,I,U, and A, respectively, and to be determined using Pontryagin’s maximal principle for the existence of optimal pairs.Theorem 9.
LetS,W,I,U, and A are optimal state solutions with associated optimal control variables u1,u2, and u3 for the optimal control model, there exist co-state variables λ1,λ2,λ3,λ4, and λ5 that satisfy
(58)dλ1dt=−∂H∂S,dλ2dt=−∂H∂W,dλ3dt=−∂H∂I,dλ4dt=−∂H∂U,dλ5dt=−∂H∂A.
With transversality or final time conditions,λ1tf=λ2tf=λ3tf=λ4tf=λ5tf=0 and where H is Hamiltonian function given in (∗). Furthermore, the optimal controls u1∗,u2∗, and u3∗ are
(59)u1∗=min1,maxβSIλ2−λ1w1N,0,u2∗=min1,maxηIλ5−λ3+kAλ5−λ3+θIλ3−λ4w2,0.Over the constraints,(60)0≤u1≤1,0≤u2≤1.Proof.
Pontryagin’s maximum principle gives the standard form of adjoint equation with transversality conditions [40]. Now, differentiating the Hamiltonian function with respect to state variables S,W,I,U and A, respectively, the adjoint equations can be written as
(61)dλ1dt=−∂H∂S=1−u1β1W+β2IN−Sβ1W+β2IN2λ1−λ2+μλ1,dλ2dt=−∂H∂W=−c1+1−u1β1SN−Sβ1W+β2IN2λ1−λ2+ξλ2−λ3+μλ2,dλ3dt=−∂H∂I=−c2+1−u1β2SN−Sβ1W+β2IN2λ1−λ2+1−u2ηλ3−λ5+u2θλ3−λ4+μλ3,dλ4dt=−∂H∂U=ϕλ4−λ3+μλ4,dλ5dt=−∂H∂A=−c3+u2κλ5−λ3+δ+μλ5.Further, the characterization of optimal controlsu1∗, and u2∗ shows that
(62)∂H∂u1=∂H∂u2=0.It follows that the optimal solution subject to constraints0≤u1≤1,0≤u2≤1 is
(63)u1∗=u1=Sβ1W+β2Iw1Nλ2−λ1,u2∗=u2=ηIλ5−λ3+kAλ5−λ3+θIλ3−λ4w2.Therefore, considering the bounds of the control, the optimal control variables are given by(64)u1∗=Sβ1W+β2Iw1Nλ2−λ1,if0<Sβ1W+β2Iw1Nλ2−λ1<1,0,ifSβ1W+β2Iw1Nλ2−λ1≤0,1,if1≤Sβ1W+β2Iw1Nλ2−λ1,u2∗=ηIλ5−λ3+kAλ5−λ3+θIλ3−λ4w2,if0<ηIλ5−λ3+kAλ5−λ3+θIλ3−λ4w2<1,0,ifηIλ5−λ3+kAλ5−λ3+θIλ3−λ4w2≤0,1,if1≤ηIλ5−λ3+kAλ5−λ3+θIλ3−λ4w2.In compact form, the optimal controls can be written as(65)u1∗=min1,maxSβ1W+β2Iw1Nλ2−λ1,0,u2∗=min1,maxηIλ5−λ3+kAλ5−λ3+θIλ3−λ4w2,0.Next, we write the optimality system using a state variables system of equations with initial conditions, a costate variables system of equations with final time conditions and optimal control solution.(66)dSdt=λ−1−u1Sβ1W+β2IN−μS,dWdt=1−u1Sβ1W+β2IN−ξ+μW,dIdt=ξW+u2κA+ϕU−1−u2η+u2θ+μI,dUdt=u2θI−∅+μU,dAdt=1−u2ηI−u2kA+δ+μA,dλ1dt=−∂H∂S=1−u1β1W+β2IN−Sβ1W+β2IN2λ1−λ2+μλ1,dλ2dt=−∂H∂W=−c1+1−u1β1SN−Sβ1W+β2IN2λ1−λ2+ξλ2−λ3+μλ2,dλ3dt=−∂H∂I=−c2+1−u1β2SN−Sβ1W+β2IN2λ1−λ2+1−u2ηλ3−λ5+θλ3−λ4+μλ3,dλ4dt=−∂H∂U=ϕλ4−λ3+μλ4,dλ5dt=−∂H∂A=−c3+u3κλ5−λ3+δ+μλ5,with conditionsλ1tf=λ2tf=λ3tf=λ4tf=λ5tf=0,S0=S0,W0=W0,I0=I0,U0=U0,andA0=A0.


## 4.1. Existence of the Optimal Control

To show the existence of optimal control, we use the approach used in [3]. We have already proved that the HIV model (2) is bounded, so this result can be used to prove the existence of optimal control over finite time interval as applied in [3, 52]. To ensure the existence of optimal control, we need to check if the following conditions are satisfied:
(i)
The set of controls and state variables be nonempty(ii)
The control setU is convex and closed(iii)
The right hand side of the state system is bounded by a linear function in the state and control variables(iv)
The integrand of objective functional is convex onU(v)
The integrand of objective functional is bounded below byk2−k1u12+u22k/2.k1,k2>0 and k>1An existence of the state system with bounded coefficients used to give condition (i). The control set is convex and closed by definition. The right hand side of the state system satisfies (iii). The state solutions are already bounded (iv). The integrand in the objective functionalc1W+c2I+c3A+1/2w1u12+w2u22 is clearly convex on U. (v) Further, from restriction on control measures, we have
(55)12wiui2≤12wi,ui∈0,1.Also, considering the preceding inequality, the integrand can be written as(56)c1W+c2I+c3A+12w1u12+w2u22≥k1u12+u22k/2−k2.wherek1=minw1/2,w2/2,k2=w2/2,k=2.Therefore, there exists optimal control measuresu1 and u2 that minimize the objective functional Ju1,u2,.


## 4.2. The Hamiltonian and Optimality System

The Pontryagin maximum principle stated the necessary conditions which are satisfied by optimal pair. Hence, by this principle, we obtained the Hamiltonian functionH defined as
(57)HS,W,I,U,A,=c1W+c2I+c3A+12w1u12+w2u22+λ1dSdt+λ2dWdt+λ3dIdt+λ4dUdt+λ5dAdt,where, λi,i=1,2,3,4,5 are the adjoint variables corresponding to state variables S,W,I,U, and A, respectively, and to be determined using Pontryagin’s maximal principle for the existence of optimal pairs.Theorem 9.
LetS,W,I,U, and A are optimal state solutions with associated optimal control variables u1,u2, and u3 for the optimal control model, there exist co-state variables λ1,λ2,λ3,λ4, and λ5 that satisfy
(58)dλ1dt=−∂H∂S,dλ2dt=−∂H∂W,dλ3dt=−∂H∂I,dλ4dt=−∂H∂U,dλ5dt=−∂H∂A.
With transversality or final time conditions,λ1tf=λ2tf=λ3tf=λ4tf=λ5tf=0 and where H is Hamiltonian function given in (∗). Furthermore, the optimal controls u1∗,u2∗, and u3∗ are
(59)u1∗=min1,maxβSIλ2−λ1w1N,0,u2∗=min1,maxηIλ5−λ3+kAλ5−λ3+θIλ3−λ4w2,0.Over the constraints,(60)0≤u1≤1,0≤u2≤1.Proof.
Pontryagin’s maximum principle gives the standard form of adjoint equation with transversality conditions [40]. Now, differentiating the Hamiltonian function with respect to state variables S,W,I,U and A, respectively, the adjoint equations can be written as
(61)dλ1dt=−∂H∂S=1−u1β1W+β2IN−Sβ1W+β2IN2λ1−λ2+μλ1,dλ2dt=−∂H∂W=−c1+1−u1β1SN−Sβ1W+β2IN2λ1−λ2+ξλ2−λ3+μλ2,dλ3dt=−∂H∂I=−c2+1−u1β2SN−Sβ1W+β2IN2λ1−λ2+1−u2ηλ3−λ5+u2θλ3−λ4+μλ3,dλ4dt=−∂H∂U=ϕλ4−λ3+μλ4,dλ5dt=−∂H∂A=−c3+u2κλ5−λ3+δ+μλ5.Further, the characterization of optimal controlsu1∗, and u2∗ shows that
(62)∂H∂u1=∂H∂u2=0.It follows that the optimal solution subject to constraints0≤u1≤1,0≤u2≤1 is
(63)u1∗=u1=Sβ1W+β2Iw1Nλ2−λ1,u2∗=u2=ηIλ5−λ3+kAλ5−λ3+θIλ3−λ4w2.Therefore, considering the bounds of the control, the optimal control variables are given by(64)u1∗=Sβ1W+β2Iw1Nλ2−λ1,if0<Sβ1W+β2Iw1Nλ2−λ1<1,0,ifSβ1W+β2Iw1Nλ2−λ1≤0,1,if1≤Sβ1W+β2Iw1Nλ2−λ1,u2∗=ηIλ5−λ3+kAλ5−λ3+θIλ3−λ4w2,if0<ηIλ5−λ3+kAλ5−λ3+θIλ3−λ4w2<1,0,ifηIλ5−λ3+kAλ5−λ3+θIλ3−λ4w2≤0,1,if1≤ηIλ5−λ3+kAλ5−λ3+θIλ3−λ4w2.In compact form, the optimal controls can be written as(65)u1∗=min1,maxSβ1W+β2Iw1Nλ2−λ1,0,u2∗=min1,maxηIλ5−λ3+kAλ5−λ3+θIλ3−λ4w2,0.Next, we write the optimality system using a state variables system of equations with initial conditions, a costate variables system of equations with final time conditions and optimal control solution.(66)dSdt=λ−1−u1Sβ1W+β2IN−μS,dWdt=1−u1Sβ1W+β2IN−ξ+μW,dIdt=ξW+u2κA+ϕU−1−u2η+u2θ+μI,dUdt=u2θI−∅+μU,dAdt=1−u2ηI−u2kA+δ+μA,dλ1dt=−∂H∂S=1−u1β1W+β2IN−Sβ1W+β2IN2λ1−λ2+μλ1,dλ2dt=−∂H∂W=−c1+1−u1β1SN−Sβ1W+β2IN2λ1−λ2+ξλ2−λ3+μλ2,dλ3dt=−∂H∂I=−c2+1−u1β2SN−Sβ1W+β2IN2λ1−λ2+1−u2ηλ3−λ5+θλ3−λ4+μλ3,dλ4dt=−∂H∂U=ϕλ4−λ3+μλ4,dλ5dt=−∂H∂A=−c3+u3κλ5−λ3+δ+μλ5,with conditionsλ1tf=λ2tf=λ3tf=λ4tf=λ5tf=0,S0=S0,W0=W0,I0=I0,U0=U0,andA0=A0.


## 5. Numerical Simulations

Next, we investigate qualitatively the effect of optimal control strategies on the spread of HIV in a population. Hence, we categorized strategies as follows:(i)
Strategy 1: only prevention control effortu1(ii)
Strategy 2: only treatment control effortu2(iii)
Strategy 3: both treatment and prevention control effortsu1andu2Solving the optimality system yields the better control. For solving the optimality system, an iterative scheme is used. We begin by using the fourth order Runge–Kutta scheme to solve the state equations with a guess for the controls over simulated time. The adjoint equations are solved using the current iterations solutions of the state equation using the backward fourth-order Runge–Kutta scheme because of the transversality conditions. The controls are then updated by combining the previous controls with the value from the characterizations. In this section, numerical simulations are done to illustrate the analytical results obtained in the above analysis. The initial values for variables of model (2) are S0=6000,W0=200,I0=1000,U0=300,A0=100,c1=20,c2=15,c3=25,w1=1,andw2=1. Moreover, the value and the source of parameters used in the simulations are given in the Table 4.Table 4
Value and source of parameters used in the simulation.
ParameterValueSourceλ200/yearAssumedβ10.9815/yearAssumedβ20.866/year[9]μ0.02/year[36]ξ0.8/yearAssumedη0.1 /year[7]ϕ0.1/year[7]θ0.5/yearAssumedκ0.1/year[30]δ1/year[30]In Figure2, the endemic equilibrium changes its stability from unstable to stable and disease-free equilibrium changes its stability from stable to unstable at the bifurcation point R0=1. Moreover, the point R0=1 is the critical point where forward bifurcation behavior of model is exhibited. Epidemiologically, interpretation gives the disease persists in the population if R0>1 and extinct if R0≤1Figure 2
Forward bifurcation diagram of model (2).In Figure3, the susceptible population dynamics is simulated considering the prevention and treatment control interventions. The susceptible population size increases due to intervention with prevention control. However, without control, more people get infected and the size of the susceptible population decreases.Figure 3
Comparison of susceptible population without and with control.In Figure4, the asymptomatic population dynamics is simulated considering the prevention and treatment control interventions. The asymptomatic population size decrease due to intervention with prevention control. However, without control, more people get infected and the size of the asymptomatic population increases with time.Figure 4
Comparison of asymptomatic population without and with control.In Figure5, the symptomatic HIV population dynamics is simulated considering the prevention and treatment control interventions. The symptomatic HIV population size decreases due to intervention with control measures. However, without control, more people get infected and the size of the HIV population increases with time.Figure 5
Comparison of symptomatic HIV population without and with control.In Figure6, the undetectable population dynamics is simulated considering the prevention and treatment control interventions. The undetectable population size decreases without intervention with controls due to replication of virus without control. However, with effective adherence to controls, the size of the undetectable population increases with time.Figure 6
Comparison of susceptible population without and with control.In Figure7, the AIDS population dynamics is simulated considering the prevention and treatment control interventions. The AIDS population size decreases due to intervention with controls. However, without controls, more people severely attacked with HIV and progress to advanced stage so that the size of the AIDS population increases.Figure 7
Comparison of AIDS population without and with control.In Figure8, the susceptible population dynamics is simulated considering only the prevention control. The susceptible population size increase due to intervention with prevention control. However, without prevention control, people get infected and the size of the susceptible population decreases with time.Figure 8
Comparison of susceptible population without treatment.In Figure9, the asymptomatic population dynamics is simulated considering the prevention control intervention. The asymptomatic population size decreases due to intervention with prevention control as less number of individuals get infected in the presence of prevention control. However, without control, more people get infected and the size of the asymptomatic population increases with time.Figure 9
Comparison of asymptomatic population without treatment.In Figure10, the symptomatic HIV population dynamics is simulated considering the prevention control intervention. The symptomatic HIV population size increases without intervention with prevention control. However, the prevention control reduces the number of individuals infected whose impact reduces the number of symptomatic HIV population.Figure 10
Comparison of HIV population without treatment.In Figure11, the undetectable population dynamics is simulated considering the only prevention control intervention. The undetectable population size remains the same without or with presence of prevention control due to absence of treatment control intervention.Figure 11
Comparison of undetectable population without treatment.In Figure12, the AIDS population dynamics is simulated considering the prevention control intervention. The AIDS population size increases without treatment control intervention control. However, with control, less people progress to the advanced stage of HIV and the size of the AIDS population decreases with time.Figure 12
Comparison of AIDS population without treatment.In Figure13, the susceptible population dynamics is simulated considering only treatment control intervention. The HIV population size increases with treatment control intervention control. However, without treatment control, more people get infected and the size of susceptible population decreases.Figure 13
Comparison of susceptible population without prevention control.Figure14 shows that the absence of the prevention control measure increases the number of HIV-infected individuals in the asymptomatic stage whereas the intervention with only treatment reduces the number HIV-infected individuals in the asymptomatic stage.Figure 14
Comparison of asymptomatic population without prevention control.In Figure15, the symptomatic HIV-infected individuals increase the HIV population without control measures, but the early intervention with treatment reduces the number of individuals that get serious with virus.Figure 15
Comparison of HIV population without prevention control.In Figure16, the dynamics of undetectable stage HIV-infected individuals are described. The simulation shows that the number of undetectable individuals increases continuously with continuous using of treatment, but in the absence or stopping the medicine causes the number of individuals at the undetectable individuals to become detectable.Figure 16
Comparison of undetectable population without prevention control.Figure17 shows that the number of individuals at AIDS stage increases with time in the absence of the medicine that inhibits the replication of the viruses in the body and prevention control measures. However, early intervention with treatment reduces the number of individuals at the AIDS stage effectively.Figure 17
Comparison of AIDS population without prevention control.In Figure18, the prevention and the treatment control profiles are simulated. The simulation results indicate that early and continuous application of the prevention and treatment control measures is effective in the controlling the HIV infection.Figure 18
Simulation of prevention and treatment control profiles.In Figure19, the simulation of adjoint variables are performed to show the conditions required in the analysis of optimal control problems. Moreover, at the final time, the adjoint variables yield zero values.Figure 19
Simulation of adjoint variables for the given time interval.In Figure20, the simulation of St versus Wt is performed; the size of St decreases to 2000 as the size of Wt increases to 1600. However, for size of St≤2000 and Wt≤1600, both population sizes increase (decrease) together.Figure 20
Simulation ofSt and Wt dynamics.In Figure21, the simulation of St versus It is performed. Initially, the size of St decrease approximately to 300 as the size of Wt increase to 4700; however, for size of St≤300 and Wt≤4700, both population sizes increase (decrease) together.Figure 21
Simulation ofSt and It dynamics.In Figure22, the simulation of St versus At is performed. Initially, the size of St decrease approximately to 300 as the size of At increase to 430; however, for size of St≤300 and At≤430 both population sizes increase (decrease) together.Figure 22
Simulation ofSt and It dynamics.In Figure23, the phase portrait of St,Wt, and It is performed. The simulation shows the initial occurrence of HIV prevalence increases. However, as time increases, less number of individuals becomes susceptible with the presence of HIV-infected individuals.Figure 23
Simulation ofSt and It dynamics.


### 5.1. Cost-Effectiveness Analysis

Incremental cost-effectiveness ratio (ICER) used to compare the differences between the costs and health outcomes of two alternative intervention strategies that compete for the same resources and is generally described as the additional cost per additional health outcome [54, 55].
(67)ICER=DifferencesbetweenthecostsDifferencebetweenhealthoutcomes.In ICER, when comparing two competing intervention strategies incrementally, one intervention should be compared with the next less effective alternative.The ICER is computed as follows:(68)ICER2=32187115950=0.2776,ICER1=6437.4−32187274020−115950=−0.1629.Comparing the computations, theICER1<ICER2. This implies to discard strategy 2 and construct a table to compare strategies 1 and 3.The ICER is computed as follows:(69)ICER1=6437.4274020=0.0235,ICER3=38624−6437.4281350−274020=4.3911.Comparing the computations, theICER1<ICER3. This implies that strategy (2) is more effective and economical over strategy (6). Therefore, discard strategy 3 and select strategy 1 as a better intervention for optimal control of HIV transmission dynamics. Therefore, based on the results obtained from Tables 5–7, the prevention control strategy is more cost-effective strategy than other control strategies to overcome the burden of HIV on human life.Table 5
Rank of strategies based on the number of the infections averted.
StrategiesTotal infections avertedTotal costsICERNo strategy00—Strategy 1 (prevention)274020321870Strategy 2 (treatment)1159506437.40Strategy 3 (prevention and treatment)2813502574.90Table 6
Comparison between strategies 2 and 1.
StrategiesTotal avertedTotal costsICERNo strategy00—Strategy 2 (treatment)115950321870.2776Strategy 1 (prevention)2740206437.4−0.1629Table 7
Comparison between strategies 1 and 3.
StrategiesTotal avertedTotal costsICERNo strategy00—Strategy 1 (prevention)2740206437.40.0235Strategy 3 (prevention and treatment)281350386244.3911


### 5.2. Result and Discussion

In this section, we present the result and discuss the results obtained from analytical analysis and numerical simulations of model (2). Particularly, we discuss the impact of ART on HIV transmission as a result of varying the value of transfer rate parameter. In Figure 1, the dynamics of population is described by schematic diagram and the model is developed using it. In Figure 2, we observe that the endemic equilibrium exhibits unstable state behavior if R0<1. But, it exhibits global asymptotic stable state behavior if R0>1. On the other hand, DFE is stable for R0<1 and unstable if R0>1. Moreover, it can be observed that the system shifts from the stable DFE to unstable at the bifurcation point R0=1. That is, there is forward bifurcation at the critical point. Further, in Figures 3–17, the impact of optimal controls is described clearly. Moreover, the presence of optimal control measures increases the number of susceptible individuals. However, the presence of optimal control measures decreases the number of individuals in the W,I, and A classes. In Figure 18, the prevention and the treatment control profiles are simulated. The simulations results indicate early and continuous application of the prevention and treatment control measures are effective in the controlling the HIV infection, whereas in Figure 19 the simulation of adjoint variables are performed to show the conditions required in the analysis of optimal control problems. Moreover, at the final time, the adjoint variables yield zero values.Further, in Figure20, the simulation of St versus Wt is performed; the size of St decreases to 2000 as the size of Wt increases to 1600. However, for size of St≤2000 and Wt≤1600, both population sizes increase (decrease) together, and in Figure 21, the simulation of St versus It is performed. Initially, the size of St decreases approximately to 300 as the size of Wt increases to 4700; however, for size of St≤300 and Wt≤4700, both population sizes increase (decrease) together; in Figure 22, the simulation of St versus At is performed. Initially, the size of St decreases approximately to 300 as the size of At increases to 430; however, for size of St≤300 and At≤430, both population sizes increase (decrease) together. In Figure 23, the phase portrait of St,Wt, and It is performed. The simulation shows the initial occurrence of HIV prevalence increases. However, as time increases, less number of individuals becomes susceptible with presence of HIV-infected individuals.Also, the investigation of cost-effectiveness analysis with possible combination with prevention and treatment control measures for HIV infection shows that early applying of prevention control measures is a better strategy than applying only ART or combined strategy with prevention and treatment control measures.


### 5.3. Conclusion

In this study, we have modified the HIV model by including an undetectable compartment that emphasizes the impact of using ART properly and continuously in the entire period of HIV infection. Also, the study identified the impact of faulty using ART by varying the parameter values that control right and faulty of using ART. From sensitivity analysis, it is concluded that increasing proper and continuous use of ART help significantly to reduce the number of reproduction number. Moreover, the analytical and numerical simulation results show that the disease-free equilibrium is both locally and globally stable forR0≤1. But, the endemic equilibrium is both locally and globally stable for R0>1. Furthermore, at the threshold point R0=1, the system experiences a forward branch. Analysis results show that HIV infection can be controlled with proper administration and use of ART. Particularly, with prevention and treatment as the optimal control measures, it reduces the number of individuals living with HIV and the progression of the HIV patient to the AIDS stage. Also, examining cost-benefit analyses of possible combinations of preventive and therapeutic management measures for HIV infection shows that preventive control measure is a better strategy if applied earlier effectively than using ART alone or combined.


## 5.1. Cost-Effectiveness Analysis

Incremental cost-effectiveness ratio (ICER) used to compare the differences between the costs and health outcomes of two alternative intervention strategies that compete for the same resources and is generally described as the additional cost per additional health outcome [54, 55].
(67)ICER=DifferencesbetweenthecostsDifferencebetweenhealthoutcomes.In ICER, when comparing two competing intervention strategies incrementally, one intervention should be compared with the next less effective alternative.The ICER is computed as follows:(68)ICER2=32187115950=0.2776,ICER1=6437.4−32187274020−115950=−0.1629.Comparing the computations, theICER1<ICER2. This implies to discard strategy 2 and construct a table to compare strategies 1 and 3.The ICER is computed as follows:(69)ICER1=6437.4274020=0.0235,ICER3=38624−6437.4281350−274020=4.3911.Comparing the computations, theICER1<ICER3. This implies that strategy (2) is more effective and economical over strategy (6). Therefore, discard strategy 3 and select strategy 1 as a better intervention for optimal control of HIV transmission dynamics. Therefore, based on the results obtained from Tables 5–7, the prevention control strategy is more cost-effective strategy than other control strategies to overcome the burden of HIV on human life.Table 5
Rank of strategies based on the number of the infections averted.
StrategiesTotal infections avertedTotal costsICERNo strategy00—Strategy 1 (prevention)274020321870Strategy 2 (treatment)1159506437.40Strategy 3 (prevention and treatment)2813502574.90Table 6
Comparison between strategies 2 and 1.
StrategiesTotal avertedTotal costsICERNo strategy00—Strategy 2 (treatment)115950321870.2776Strategy 1 (prevention)2740206437.4−0.1629Table 7
Comparison between strategies 1 and 3.
StrategiesTotal avertedTotal costsICERNo strategy00—Strategy 1 (prevention)2740206437.40.0235Strategy 3 (prevention and treatment)281350386244.3911


## 5.2. Result and Discussion

In this section, we present the result and discuss the results obtained from analytical analysis and numerical simulations of model (2). Particularly, we discuss the impact of ART on HIV transmission as a result of varying the value of transfer rate parameter. In Figure 1, the dynamics of population is described by schematic diagram and the model is developed using it. In Figure 2, we observe that the endemic equilibrium exhibits unstable state behavior if R0<1. But, it exhibits global asymptotic stable state behavior if R0>1. On the other hand, DFE is stable for R0<1 and unstable if R0>1. Moreover, it can be observed that the system shifts from the stable DFE to unstable at the bifurcation point R0=1. That is, there is forward bifurcation at the critical point. Further, in Figures 3–17, the impact of optimal controls is described clearly. Moreover, the presence of optimal control measures increases the number of susceptible individuals. However, the presence of optimal control measures decreases the number of individuals in the W,I, and A classes. In Figure 18, the prevention and the treatment control profiles are simulated. The simulations results indicate early and continuous application of the prevention and treatment control measures are effective in the controlling the HIV infection, whereas in Figure 19 the simulation of adjoint variables are performed to show the conditions required in the analysis of optimal control problems. Moreover, at the final time, the adjoint variables yield zero values.Further, in Figure20, the simulation of St versus Wt is performed; the size of St decreases to 2000 as the size of Wt increases to 1600. However, for size of St≤2000 and Wt≤1600, both population sizes increase (decrease) together, and in Figure 21, the simulation of St versus It is performed. Initially, the size of St decreases approximately to 300 as the size of Wt increases to 4700; however, for size of St≤300 and Wt≤4700, both population sizes increase (decrease) together; in Figure 22, the simulation of St versus At is performed. Initially, the size of St decreases approximately to 300 as the size of At increases to 430; however, for size of St≤300 and At≤430, both population sizes increase (decrease) together. In Figure 23, the phase portrait of St,Wt, and It is performed. The simulation shows the initial occurrence of HIV prevalence increases. However, as time increases, less number of individuals becomes susceptible with presence of HIV-infected individuals.Also, the investigation of cost-effectiveness analysis with possible combination with prevention and treatment control measures for HIV infection shows that early applying of prevention control measures is a better strategy than applying only ART or combined strategy with prevention and treatment control measures.


## 5.3. Conclusion

In this study, we have modified the HIV model by including an undetectable compartment that emphasizes the impact of using ART properly and continuously in the entire period of HIV infection. Also, the study identified the impact of faulty using ART by varying the parameter values that control right and faulty of using ART. From sensitivity analysis, it is concluded that increasing proper and continuous use of ART help significantly to reduce the number of reproduction number. Moreover, the analytical and numerical simulation results show that the disease-free equilibrium is both locally and globally stable forR0≤1. But, the endemic equilibrium is both locally and globally stable for R0>1. Furthermore, at the threshold point R0=1, the system experiences a forward branch. Analysis results show that HIV infection can be controlled with proper administration and use of ART. Particularly, with prevention and treatment as the optimal control measures, it reduces the number of individuals living with HIV and the progression of the HIV patient to the AIDS stage. Also, examining cost-benefit analyses of possible combinations of preventive and therapeutic management measures for HIV infection shows that preventive control measure is a better strategy if applied earlier effectively than using ART alone or combined.



---

*Source: 4754426-2023-02-06.xml*