# Threshold Dynamics and Competitive Exclusion in a Virus Infection Model with General Incidence Function and Density-Dependent Diffusion

Authors: Xiaosong Tang; Zhiwei Wang; Jianping Yang
Journal: Complexity
 (2020)

---

## Abstract

In this paper, we investigate single-strain and multistrain viral infection models with general incidence function and density-dependent diffusion subject to the homogeneous Neumann boundary conditions. For the single-strain viral infection model, by using the linearization method and constructing appropriate Lyapunov functionals, we obtain that the global threshold dynamics of the model is determined by the reproductive numbers for viral infectionℛ0. For the multistrain viral infection model, we have discussed the competitive exclusion problem. If the reproduction number ℛi for strain i is maximal and larger than one, the steady state Ei corresponding to the strain i is globally stable. Thus, competitive exclusion happens and all other strains die out except strain i. Meanwhile, we can prove that the single-strain and multistrain viral infection models are well posed. Furthermore, numerical simulations are also carried out to illustrate the theoretical results, which is seldom seen in the relevant known literatures.


---

## Body

## 1. Introduction

Nowadays, more and more people in the world are dying of various diseases such as AIDS, avian influenza, cholera, Ebola, and Zika virus. To explore the mechanisms of these diseases, scientists have proposed plenty of mathematical models describing the transmission of disease such as the infectious disease model (SI, SIR, SEI, and SIS) and the within-host virus model (HBV, HCV, HIV, Ebola, and Zika). Since Nowak et al. [1, 2] has established the basic within-host virus dynamics model, more and more scholars have devoted to investigate the within-host virus dynamics model and obtained so many significant dynamical results [3–16].Assuming that the motion of the virus belongs to the Fickian diffusion [17], Wang et al. [5] have firstly investigated a delayed diffusive hepatitis B virus model. Here, we have to mention that the equation with Fickian diffusion implies that random walkers in the system have no interactions. In fact, random walkers in the system may have interactions and the density-dependent diffusion happens [17] and we may establish the equation with density-dependent diffusion. So, parabolic or elliptic systems with density-dependent diffusion have been widely studied by many researchers [18–21]. More recently, Wang et al. have investigated a virus infection model with density-dependent diffusion and Holling II type [6] or Beddington–DeAngelis type incidence function [7] subjected to the homogeneous Neumann boundary conditions. And, they have obtained the significant results of threshold dynamics and competitive exclusion.It is well known that there are many factors which affect the dynamics of within-host models. One crucial component is incidence function, which can be classified into many different types, such as bilinear type [2], Holling II type [22], Beddington–DeAngelis type [23], and Crowley–Martin type [24]. Many researchers have concentrated on these within-host models or epidemic models and have obtained some interesting results such as global threshold dynamics [8, 9], delay-induced Hopf bifurcation [10, 11], competitive exclusion [6, 7], traveling wave [12, 13], and pattern formation [14–16]. Recently, Korobeinikov [25] studied an ordinary differential equation model with general incidence function and dose-dependent parasite reproduction and virulence. Huang et al. [26] investigated a viral infection model involving time delay and general incidence function. McCluskey and Yang [27] considered a diffusive viral infection model with general incidence function and time delay. Hattaf and Yousfi [28] studied a diffusive HBV model with two delays and general incidence function. Miao et al. [29] considered the global stability of a diffusive humoral immunity virus infection model with time delay and general incidence function and so on.Based on the aforementioned works, we firstly propose a single-strain viral infection model with density-dependent diffusion and general incidence function as follows:(1)∂ux,t∂t=∇⋅Duu∇u+λ−dux,t−gux,t,vx,t,∂wx,t∂t=∇⋅Dww∇w+gux,t,vx,t−δwx,t,∂vx,t∂t=∇⋅Dvv∇v+pwx,t−cvx,t,where Duu=um1, Dww=wm2, and Dvv=vm3 are the density-dependent diffusion coefficients for mi>0, i=1,2,3; ux,t represents the densities of uninfected cells; wx,t represents the densities of infected cells; vx,t represents the densities of free virus at position x and at time t, respectively; constant λ represents the production rate of the uninfected cells; the death rates of uninfected cells, infected cells, and free virus are represented by the parameters d, δ, and c, respectively; the infected cells produce the free virions at a rate pw; further, we assume that the incidence function in this paper is a general incidence function given by gu,v, which represents the contacts between target cells and viruses; and the general incidence function gu,v satisfies the following conditions [27].H
1
g
:
ℝ
+
2
⟶
ℝ
+ is differentiable; there exists η>0 such that gu,v≤ηuv for all u,v≥0; gu,v>0 for all u,v>0; gu,0=g0,v=0 for all u,v≥0; ∂gu,v/∂u>0 for all u≥0 and v>0; ∂gu,v/∂v≥0 for all u,v≥0; v∂gu,v/∂vgu,v≤0 for all u,v≥0.These assumptions are biologically motivated, see [27] in detail. Moreover, for a,b,c,β>0, we choose the incidence functions such as gu,v=βuv/1+bv [22], gu,v=βuv/1+au+bv [23], and gu,v=βuv/1+au+bv+cuv [24], which satisfy H1 clearly.Model (1) we considered in this paper satisfies the following Neumann boundary conditions:(2)Duu∂u∂ν=0,Dww∂w∂ν=0,Dvv∂v∂ν=0,on∂Ω×0,+∞,and initial conditions(3)ux,0=ϕ1x≥0,wx,0=ϕ2x≥0,vx,0=ϕ3x≥0,x∈Ω¯.Here, we give out a bounded domainΩ in ℝn and its boundary ∂Ω is smooth, and the outward normal derivative on ∂Ω is denoted by ∂/∂ν. The nonnegative functions ϕix,0 i=1,2,3 are Hölder continuous and satisfy ∂ϕi/∂ν=0 i=1,2,3 on ∂Ω.For (1), we can obtain its steady-state system,(4)−∇⋅Duu∇u=λ−dux,t−gux,t,vx,t=f1u,v,w,−∇⋅Dww∇w=gux,t,vx,t−δwx,t=f2u,v,w,−∇⋅Dvv∇v=pwx,t−cvx,t=f3u,v,w,which satisfies boundary conditions (2). We denote(5)fu,v,w=f1u,v,w,f2u,v,w,f3u,v,w,and let fiu,v,wi=1,2,3 satisfy the following conditions.H
2 The following relations(6)γ1um1+∂f1∂u≥0,γ2wm2+∂f2∂w≥0,γ3vm3+∂f3∂v≥0,hold for constants γi≥0 i=1,2,3.Since Gause has discussed the principles that govern species competitive exclusion [30], competitive exclusion has been considered in such settings such as within-host models, ecological models, and epidemiological models. Here, we have to point out that Bremermann and Thieme [31] have given out the first rigorous proof of the competitive exclusion principle in the epidemiological model setting. For the multistrain SIR model, they have proved that the strain with the largest reproduction number persists in the population, while the remaining strains die out. This principle was also used in within-host models and has been studied by many scholars [32–37]. For instance, De Leenheer and Pilyugin [32] proved the global competitive exclusion in an ODE version with a multistrain virus. Feng and Velasco-Hernández [33] discussed competitive exclusion in a vector-host model for the dengue fever. Iggidr et al. [34] analyzed globally a new malaria within-host model through a competitive exclusion principle. Wang et al. investigated the competitive exclusion in the virus infection model with density-dependent diffusion and two types of incidence functions [6, 7]. Duan et al. [36] considered the competitive exclusion in the diffusive virus infection model with age of infection and general incidence function, and so on. So, in this paper, we further investigate the competitive exclusion problem for the multistrain virus infection model. Setting i=1,2,…,n, we can obtain the following multistrain viral infection model with density-dependent diffusion and general incidence function:(7)∂ux,t∂t=∇⋅Duu∇u+λ−dux,t−∑i=1ngiux,t,vix,t,∂wix,t∂t=∇⋅Dwiwi∇wi+giux,t,vix,t−δiwix,t,∂vix,t∂t=∇⋅Dvivi∇vi+piwix,t−civix,t,i=1,2,…,n.Here, we have to stress that the Neumann boundary conditions, initial conditions, and steady-state system of Model (7) are described in Section 3. However, to the best of our knowledge, there are no results on the above models (1) and (7) such as global threshold dynamics and competitive exclusion. So, in this paper, we shall focus on the global threshold dynamics for the single-strain viral infection model (1) and competitive exclusion for the multi-strain viral infection model (7). It is worthwhile to mention that the incidence function discussed in [6, 7] is a special incidence function. But, the opposite case is more difficult and complicated. Through constructing an appropriate Lyapunov functional, Lemma 2, and mathematic analysis skills, we overcome the difficulties successfully and obtain the global threshold dynamics for the single-strain viral infection model (1) and competitive exclusion for the multi-strain viral infection model (7).The rest of this paper is organized as follows: in Section2, for the single-strain viral infection model (1), we shall prove the existence, uniqueness, and boundedness of solutions to problems (1)–(4) and consider the global threshold dynamics of the model through constructing appropriate Lyapunov functionals, which is determined by the reproductive numbers for viral infection ℛ0. In Section 3, for the multistrain viral infection model (7), we also prove the existence, uniqueness, and boundedness of solutions to problems (7) and (42) and discuss the competitive exclusion problem. If the reproduction number ℛi for strain i is maximal and larger than one, the steady state Ei corresponding to the strain i is globally stable. Thus, competitive exclusion happens and all other strains die out except strain i. In Section 4, we illustrate our results with numerical simulations, which support and extend the theoretical results. The paper ends with a conclusion.


## 2. Single-Strain Viral Infection Model

In this section, we shall prove firstly that the single-strain viral infection model (1) is well posed by the method used in [18]. Then, we continue to discuss the existence of equilibria for model (1). Finally, through constructing appropriate Lyapunov functionals, we investigate the global threshold dynamics of model (1), which is determined by the reproductive numbers for viral infection ℛ0.


### 2.1. Well Posedness

In this section, we give out some definitions and results first, which can help us to prove that the solutions to problems (1)–(4) exist and are unique and bounded. Similar to the definition method in [6, 7], we can definite that two pairs of functions U˜=u˜,w˜,v˜, U^=u^,w^,v^ in CΩ¯×0,T∩C1,2Ω×0,T and U˜s=u˜s,w˜s,v˜s, U^s=u^s,w^s,v^s in CΩ¯×0,T∩C1,2Ω×0,T are called coupled upper and lower solutions to system (1)–(3) and to system (4), respectively. Then, we can easily conclude that each pair of coupled upper and lower solutions U˜s, U^s of problems (1)–(3) is ones of (4). In what follows, we shall search a pair of constant solutions of (4), which is called upper and lower solutions of (4).Lemma 1.
For constantl1>0 and L1>0 and any sufficiently large constant Lj>0 and small constant lj>0j=2,3, the coupled upper and lower solutions to (4) is a constant pair L=L1,L2,L3 and l=l1,l2,l3.Proof.
Obviously, the constant pairL=L1,L2,L3 and l=l1,l2,l3 satisfy the homogeneous Neumann boundary conditions (2). Further, by means of condition H1, we can confirm that the differential inequalities in (4) can also be satisfied by the constant pair L and l if(8)λ−dL1−gL1,l3≤λ−dL1≤0,gL1,L3−δL2≤ηL1L3−δL2≤0,pL2−cL3≤0;λ−dl1−gl1,L3≥λ−dl1−ηl1L3=λ−d+ηL3l1≥0,gl1,l3−δl2≥gvl1,l3l3−δl2≥0,pl2−cl3≥0,and then, the above inequalities can be fulfilled for constants l1 and L1, sufficiently large Lj, and small ljj=2,3 if(9)λd≤L1,ληL3dδ≤L2≤cL3p,l1≤λd+ηL3,l2≤gvl1,l3l3δ,l3≤pl2c.
To obtain the upper solution of (4), we only need to take L1=λ/d and sufficiently large L2 and L3 satisfying ληL3/dδ≤L2≤cL3/p. Meanwhile, choosing sufficiently small l3, l2=cl3/p, and l1 satisfies cδ/p≤gvl1,l3 and l1≤λ/d+ηL3, and we can also obtain the lower solution of (4). This means that the coupled upper and lower solutions of (4) are a constant pair L and l. This completes the proof of Lemma 1.
According to the results of Lemma1 and method in [6, 7], we can conclude that problems (1)–(4) is well-posed similarly. That is, the solutions to problems (1)–(4) exist and are unique and bounded.


### 2.2. Equilibria

To describe that a single infected cell generates the average number of newly infected cells when the disease has just entered the body, we define the following basic reproduction number:(10)ℛ0=pcδ∂gλ/d,0∂v,which is independent of spatial variety. Obviously, E0=λ/d,0,0 is always an infection-free equilibrium of model (1). Similar to the analysis and discussion in [27], together with the basic reproduction number ℛ0, it is not hard to see that model (1) has only a positive equilibrium E∗=u∗,w∗,v∗ if ℛ0>1. So, we have the results on the equilibria of model (1) as follows.Theorem 1.
The infection-free equilibriumE0=λ/d,0,0 of model (1) always exists. If ℛ0>1, model (1) has only an infection equilibrium E∗=u∗,w∗,v∗.


### 2.3. Stability Properties of Model (1)

#### 2.3.1. The Local Stability of Model (1)

In this section, we shall investigate the local stability of model (1). Under the homogeneous Neumann boundary conditions, assume that 0=ν0<ν1<⋅⋅⋅<νq<⋅⋅⋅ be the eigenvalues of −Δ on Ω and for q=0,1,2,…, and Eνq be the eigenfunctions space corresponding to νq in C1Ω. Suppose that Eνq has an orthonormal basis φqm:m=1,2,…,dimEνq, X=C1Ω3, and Xqm=cφqm:c∈ℝ3. Thus, we obtain(11)X=⊕q=1∞Xq,Xq=⊕m=1dimEνqXqm.Linearizing system (1) at the equilibrium u˜,w˜,v˜ is as follows:(12)∂ux,t∂t=Duu˜Δu−dux,t−∂gu˜,v˜∂uux,t−∂gu˜,v˜∂vvx,t,∂wx,t∂t=Dww˜Δw+∂gu˜,v˜∂uux,t+∂gu˜,v˜∂vvx,t−δwx,t,∂vx,t∂t=Dvv˜Δv+pwx,t−cvx,t.Note that, under the linearization,Xq is an invariant for each q=0,1,2,…. We let the expression of solution be Qx,t=eμtψx, where ψx∈Xq satisfies −Δψ=νqψ. Thus, the characteristic equation of (12) is(13)μ+Duu˜νq+d+∂gu˜,v˜∂uμ+Dww˜νq+δμ+Dvv˜νq+c−p∂gu˜,v˜∂v+p∂gu˜,v˜∂u∂gu˜,v˜∂v=0.When the equilibrium of model (1) is E0=λ/d,0,0, the characteristic equation (13) is changed into(14)μ+Duλdνq+dμ+δμ+c−p∂gλ/d,0∂v=0.Obviously, equation (14) has a negative real root σ=−Duλ/dνq−d. We now discuss the transcendental equation as follows:(15)μ+δμ+c−p∂gλ/d,0∂v=0.Set(16)Hμ,q=μ+δμ+c−p∂gλ/d,0∂v.Case 1.
ℛ
0
<
1. By direct computation, we have(17)H0,q=cδ−p∂gλ/d,0∂v=cδ1−ℛ0>0,and for μ≥0,(18)∂Hμ,q∂μ=2μ+δ+c>0.
So, there is only a negative real root for equation (14).
Assume that there be a complex rootρ+iϖ with ρ≥0, then we can conclude that(19)μδ+1>1,μc+1>1,pcδ∂gλ/d,0∂v=ℛ0<1,which implies(20)μδ+1μc+1>ℛ0=pcδ∂gλ/d,0∂v.
This contradicts to equation (15). Then, all roots of equation (14) have no positive real parts yet. Thus, the infection-free equilibrium E0 is locally asymptotically stable if ℛ0<1.Case 2.
ℛ
0
>
1. Obviously, we need to only consider the space X0 corresponding to ν0=0 when q=0. By direct computation, we get(21)H0,0=cδ1−ℛ0<0,limμ⟶+∞Hμ,0=+∞.
This implies that the equationHμ,0=0 must hold for some μ0>0. Thus, there exists at least one positive root of equation (14). So, the infection-free equilibrium E0 is not stable for ℛ0>1. Then, we have the following results.Theorem 2.
The infection-free equilibriumE0 of model (1) is locally asymptotically stable if ℛ0<1. However, E0 is unstable if ℛ0>1.
When the equilibrium of model (1) is E∗=u∗,w∗,v∗, the characteristic equation (13) is changed into(22)μ3+Aμ2+Bμ+C=0,where(23)A=Duu∗νq+Dww∗νq+Dvv∗νq+d+δ+c+∂gu∗,v∗∂u,B=Duu∗νq+d+∂gu∗,v∗∂uDww∗νq+Dvv∗νq+δ+c+Dww∗νq+δDvv∗νq+c−p∂gu∗,v∗∂v,C=Duu∗νq+d+∂gu∗,v∗∂uDww∗νq+δDvv∗νq+c−p∂gu∗,v∗∂v+p∂gu∗,v∗∂u∂gu∗,v∗∂v.
It is easy to see thatA>0. From the condition H1, together with gu∗,v∗=δw∗ and pw∗=cv∗, we get(24)p∂gu∗,v∗∂v≤pgu∗,v∗v∗=pδw∗pw∗/c=cδ,which yields B>0 and C>0. By direct computation, we have that(25)AB−C=Duu∗νq+Dww∗νq+Dvv∗νq+d+δ+∂gu∗,v∗∂uDuu∗νq+d+∂gu∗,v∗∂uDww∗νq+Dvv∗νq+δ+c+cDuu∗νq+dc+Dww∗νq+Dvv∗νq+δ+c∂gu∗,v∗∂uDww∗νq+Dvv∗νq+c+∂gu∗,v∗∂ucδ−p∂gu∗,v∗∂v>0.
So, we can see that all eigenvalues of (22) have no positive real parts by the Routh–Hurwitz criterion. Therefore, E∗ is locally asymptotically stable for ℛ0>1. Thus, we obtain the following results.Theorem 3.
The infection equilibriumE∗ of model (1) is locally asymptotically stable if ℛ0>1.


#### 2.3.2. The Global Stability of Model (1)

In this section, by constructing appropriate Lyapunov functionals, we continue to investigate the global stability of model (1). To this end, we need the following lemma.Lemma 2.
Under homogeneous Neumann boundary condition, for a density-dependent diffusion problem, we can obtain the following result:(26)∫Ω1gu∇⋅Du∇udx=∫ΩDug′u∇u2g2udx,where Ω is a bounded domain in Rn and its boundary ∂Ω is smooth.Proof.
According to the homogeneous Neumann boundary condition, we get(27)∫Ω1gu−∇⋅Du∇udx=∫Ω1guD′u∇u2dx+∫Ω1guDuΔudx=∫Ω1guD′u∇u2dx+∫∂Ω1guDu∂u∂νds−∫ΩDu∇1gu⋅∇udx−∫Ω1guD′u∇u2dx=∫ΩDug′u∇u2g2udx.
This completes the proof of Lemma2.Remark 1.
Choosinggu=1, u, and 1/u in Lemma 2, respectively, we can obtain that(i)
∫Ω∇⋅Du∇udx=0;(ii)
∫Ω1/u∇⋅Du∇udx=∫ΩDu∇u2/u2dx;(iii)
∫Ωu∇⋅Du∇udx=−∫ΩDu∇u2dx;which implies that Lemma 2 generalized the results in [6, 7].
Firstly, we discuss the global asymptotic stability of the infection-free equilibriumE0 of model (3).Theorem 4.
Ifℛ0≤1, then the infection-free equilibrium E0 of model (3) is globally asymptotically stable.Proof.
Define the following Lyapunov functional:(28)Vt=∫Ωwx,t+δpvx,tdx.
Through calculating the derivative ofVt about time, we can get(29)dVtdt=∫Ω∂w∂t+δp∂v∂tdx=∫Ω∇⋅Dw∇w+gux,t,vx,t−δwx,t+δp∇⋅Dv∇v+δwx,t−cδpvx,tdx=∫Ωgux,t,vx,t−cδpvx,tdx+∫Ω∇⋅Dw∇wdx+δp∫Ω∇⋅Dv∇vdx.
Note that∂u/∂t≤λ−du, then we have limt⟶∞ux,t≤λ/d, which means that all omega limit points satisfy u≤λ/d. So, we need to only consider solutions for which u≤λ/d. Since ∂g/∂u≥0, we can conclude the inequality gu,v≤gλ/d,v. Under the condition H1, together with the fact that gλ/d,v/v is decreasing in v, we can obtain(30)gux,t,vx,t−cδpvx,t=cδpvpcδgu,vv−1≤cδpvpcδgλ/d,vv−1≤cδpvlimv⟶0+pcδgλ/d,vv−1=cδpvpcδ∂gλ/d,0∂v−1=cδpvℛ0−1.
Meanwhile, from Remark1, we get(31)∫Ω∇⋅Dw∇wdx=0=δp∫Ω∇⋅Dv∇vdx.
Thus, we obtain(32)dVtdt≤cδpℛ0−1∫Ωvdx.
So, we have thatdVt/dt≤0 when ℛ0≤1. Further, it is not hard to believe that dVt/dt≤0 leads to either v=0 or ℛ0=1 and u=λ/d. It is easy to confirm that the singleton E0 is just the largest compact invariant set in u,w,vdVt/dt≤0. When ℛ0≤1, we can know that the infection-free equilibrium E0 of model (3) is globally asymptotically stable by the LaSalle’s invariance principle. The proof of Theorem 4 is complete.
Next, we investigate the global asymptotic stability problem of the infection equilibriumE∗ of model (3). To prove the following theorem, we need the following hypothesis:(33)H3vv∗−gu,vgu,v∗gu,vgu,v∗−1≥0,for allu,v>0.Theorem 5.
SupposeH3 is satisfied. If ℛ0>1, then the infection equilibrium E∗ of model (3) is globally asymptotically stable.Proof.
LetHz=z−1−lnz,z>0, and define(34)Vt=∫Ωu−u∗−∫u∗ux,tgu∗,v∗gs,v∗ds+w∗Hww∗+δpv∗Hvv∗dx,which is called a Lyapunov functional. Obviously, 0 is the minimum value of function Hz at z=1. This leads to Hz≥0. Meanwhile, 0 is also the global minimum of function G:u⟶u−u∗−∫u∗ugu∗,v∗/gs,v∗ds at u=u∗, which implies Gu≥0 for any x>0. Therefore, for all x∈Ω and t≥0, Vt≥0 with equality holding if and only if ux,t/u∗=w/w∗=vx,t/v∗=1.
Through calculating the derivative ofVt about time, we can get(35)dVtdt=∫Ω1−gu∗,v∗gu,v∗∂u∂t+1−w∗w∂w∂t+δp1−v∗v∂v∂tdx=∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇u+λ−dux,t−gux,t,vx,t+1−w∗w∇⋅Dw∇w+gux,t,vx,t−δwx,t+δp1−v∗v∇⋅Dv∇v+pwx,t−cvx,tdx.
Note thatλ=du∗+gu∗,v∗, δw∗=gu∗,v∗, and cv∗=pw∗. So, we have(36)dVtdt=∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇u+du∗−dux,t+gu∗,v∗−gux,t,vx,t+1−w∗w∇⋅Dw∇w+gux,t,vx,t−wx,tw∗gu∗,v∗+δp1−v∗v∇⋅Dv∇v+pwx,t−vx,tv∗w∗dx=∫Ωdu∗1−uu∗1−gu∗,v∗gu,v∗+gu∗,v∗1−gu∗,v∗gu,v∗1−gux,t,vx,tgu∗,v∗+1−w∗wgux,t,vx,tgu∗,v∗−wx,tw∗+1−v∗vwx,tw∗−vx,tv∗dx+∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇udx+∫Ω1−w∗w∇⋅Dw∇wdx+δp∫Ω1−v∗v∇⋅Dv∇vdx=∫Ωdu∗1−uu∗1−gu∗,v∗gu,v∗+gu∗,v∗Hgux,t,vx,tgux,t,v∗−Hgu∗,v∗gux,t,v∗−Hgux,t,vx,tgu∗,v∗+Hgux,t,vx,tgu∗,v∗−Hgux,t,vx,tw∗gu∗,v∗wx,tx,t−Hwx,tw∗+Hwx,tw∗−Hwx,tv∗vx,tw∗−Hvx,tv∗dx+∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇udx+∫Ω1−w∗w∇⋅Dw∇wdx+δp∫Ω1−v∗v∇⋅Dv∇vdx=∫Ωdu∗1−uu∗1−gu∗,v∗gu,v∗+gu∗,v∗Hgux,t,vx,tgux,t,v∗−Hgu∗,v∗gux,t,v∗−Hgux,t,vx,tw∗gu∗,v∗wx,t−Hwx,tv∗vx,tw∗−Hvx,tv∗dx+∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇udx+∫Ω1−w∗w∇⋅Dw∇wdx+δp∫Ω1−v∗v∇⋅Dv∇vdx.
From the hypothesisH1, we can know that gu,v∗ is increasing in u strictly. Hence, we have(37)1−uu∗1−gu∗,v∗gu,v∗≤0,and 1−u/u∗1−gu∗,v∗/gu,v∗ equals 0 only if u=u∗. Based on the assumption H3, we can conclude that gux,t,vx,t/gux,t,v∗ lies between 1 and vx,t/v∗. Thus, by the fact that Hz is a monotone function on each side of 1, we can obtain that Hgux,t,vx,t/gux,t,v∗≤Hvx,t/v∗. Therefore,(38)dVtdt≤−gu∗,v∗∫ΩHgu∗,v∗gu,v∗+Hwx,tv∗vx,tw∗+Hgux,t,vx,tw∗gu∗,v∗wx,tdx+∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇udx+∫Ω1−w∗w∇⋅Dw∇wdx+δp∫Ω1−v∗v∇⋅Dv∇vdx.
From Lemma2 and Remark 1, it is easy to conclude that(39)dVtdt≤0.
Hence,E∗ is stable and dVt/dt=0 iff u=u∗, w=w∗, and v=v∗. It is easy to confirm that the singleton E∗ is just the largest compact invariant set in u,w,vdVt/dt=0. According to the LaSalle’s invariance principle, we can know that the infection equilibrium E∗ of model (3) is globally asymptotically stable when ℛ0>1. This completes the proof of Theorem 5.


## 2.1. Well Posedness

In this section, we give out some definitions and results first, which can help us to prove that the solutions to problems (1)–(4) exist and are unique and bounded. Similar to the definition method in [6, 7], we can definite that two pairs of functions U˜=u˜,w˜,v˜, U^=u^,w^,v^ in CΩ¯×0,T∩C1,2Ω×0,T and U˜s=u˜s,w˜s,v˜s, U^s=u^s,w^s,v^s in CΩ¯×0,T∩C1,2Ω×0,T are called coupled upper and lower solutions to system (1)–(3) and to system (4), respectively. Then, we can easily conclude that each pair of coupled upper and lower solutions U˜s, U^s of problems (1)–(3) is ones of (4). In what follows, we shall search a pair of constant solutions of (4), which is called upper and lower solutions of (4).Lemma 1.
For constantl1>0 and L1>0 and any sufficiently large constant Lj>0 and small constant lj>0j=2,3, the coupled upper and lower solutions to (4) is a constant pair L=L1,L2,L3 and l=l1,l2,l3.Proof.
Obviously, the constant pairL=L1,L2,L3 and l=l1,l2,l3 satisfy the homogeneous Neumann boundary conditions (2). Further, by means of condition H1, we can confirm that the differential inequalities in (4) can also be satisfied by the constant pair L and l if(8)λ−dL1−gL1,l3≤λ−dL1≤0,gL1,L3−δL2≤ηL1L3−δL2≤0,pL2−cL3≤0;λ−dl1−gl1,L3≥λ−dl1−ηl1L3=λ−d+ηL3l1≥0,gl1,l3−δl2≥gvl1,l3l3−δl2≥0,pl2−cl3≥0,and then, the above inequalities can be fulfilled for constants l1 and L1, sufficiently large Lj, and small ljj=2,3 if(9)λd≤L1,ληL3dδ≤L2≤cL3p,l1≤λd+ηL3,l2≤gvl1,l3l3δ,l3≤pl2c.
To obtain the upper solution of (4), we only need to take L1=λ/d and sufficiently large L2 and L3 satisfying ληL3/dδ≤L2≤cL3/p. Meanwhile, choosing sufficiently small l3, l2=cl3/p, and l1 satisfies cδ/p≤gvl1,l3 and l1≤λ/d+ηL3, and we can also obtain the lower solution of (4). This means that the coupled upper and lower solutions of (4) are a constant pair L and l. This completes the proof of Lemma 1.
According to the results of Lemma1 and method in [6, 7], we can conclude that problems (1)–(4) is well-posed similarly. That is, the solutions to problems (1)–(4) exist and are unique and bounded.


## 2.2. Equilibria

To describe that a single infected cell generates the average number of newly infected cells when the disease has just entered the body, we define the following basic reproduction number:(10)ℛ0=pcδ∂gλ/d,0∂v,which is independent of spatial variety. Obviously, E0=λ/d,0,0 is always an infection-free equilibrium of model (1). Similar to the analysis and discussion in [27], together with the basic reproduction number ℛ0, it is not hard to see that model (1) has only a positive equilibrium E∗=u∗,w∗,v∗ if ℛ0>1. So, we have the results on the equilibria of model (1) as follows.Theorem 1.
The infection-free equilibriumE0=λ/d,0,0 of model (1) always exists. If ℛ0>1, model (1) has only an infection equilibrium E∗=u∗,w∗,v∗.


## 2.3. Stability Properties of Model (1)

### 2.3.1. The Local Stability of Model (1)

In this section, we shall investigate the local stability of model (1). Under the homogeneous Neumann boundary conditions, assume that 0=ν0<ν1<⋅⋅⋅<νq<⋅⋅⋅ be the eigenvalues of −Δ on Ω and for q=0,1,2,…, and Eνq be the eigenfunctions space corresponding to νq in C1Ω. Suppose that Eνq has an orthonormal basis φqm:m=1,2,…,dimEνq, X=C1Ω3, and Xqm=cφqm:c∈ℝ3. Thus, we obtain(11)X=⊕q=1∞Xq,Xq=⊕m=1dimEνqXqm.Linearizing system (1) at the equilibrium u˜,w˜,v˜ is as follows:(12)∂ux,t∂t=Duu˜Δu−dux,t−∂gu˜,v˜∂uux,t−∂gu˜,v˜∂vvx,t,∂wx,t∂t=Dww˜Δw+∂gu˜,v˜∂uux,t+∂gu˜,v˜∂vvx,t−δwx,t,∂vx,t∂t=Dvv˜Δv+pwx,t−cvx,t.Note that, under the linearization,Xq is an invariant for each q=0,1,2,…. We let the expression of solution be Qx,t=eμtψx, where ψx∈Xq satisfies −Δψ=νqψ. Thus, the characteristic equation of (12) is(13)μ+Duu˜νq+d+∂gu˜,v˜∂uμ+Dww˜νq+δμ+Dvv˜νq+c−p∂gu˜,v˜∂v+p∂gu˜,v˜∂u∂gu˜,v˜∂v=0.When the equilibrium of model (1) is E0=λ/d,0,0, the characteristic equation (13) is changed into(14)μ+Duλdνq+dμ+δμ+c−p∂gλ/d,0∂v=0.Obviously, equation (14) has a negative real root σ=−Duλ/dνq−d. We now discuss the transcendental equation as follows:(15)μ+δμ+c−p∂gλ/d,0∂v=0.Set(16)Hμ,q=μ+δμ+c−p∂gλ/d,0∂v.Case 1.
ℛ
0
<
1. By direct computation, we have(17)H0,q=cδ−p∂gλ/d,0∂v=cδ1−ℛ0>0,and for μ≥0,(18)∂Hμ,q∂μ=2μ+δ+c>0.
So, there is only a negative real root for equation (14).
Assume that there be a complex rootρ+iϖ with ρ≥0, then we can conclude that(19)μδ+1>1,μc+1>1,pcδ∂gλ/d,0∂v=ℛ0<1,which implies(20)μδ+1μc+1>ℛ0=pcδ∂gλ/d,0∂v.
This contradicts to equation (15). Then, all roots of equation (14) have no positive real parts yet. Thus, the infection-free equilibrium E0 is locally asymptotically stable if ℛ0<1.Case 2.
ℛ
0
>
1. Obviously, we need to only consider the space X0 corresponding to ν0=0 when q=0. By direct computation, we get(21)H0,0=cδ1−ℛ0<0,limμ⟶+∞Hμ,0=+∞.
This implies that the equationHμ,0=0 must hold for some μ0>0. Thus, there exists at least one positive root of equation (14). So, the infection-free equilibrium E0 is not stable for ℛ0>1. Then, we have the following results.Theorem 2.
The infection-free equilibriumE0 of model (1) is locally asymptotically stable if ℛ0<1. However, E0 is unstable if ℛ0>1.
When the equilibrium of model (1) is E∗=u∗,w∗,v∗, the characteristic equation (13) is changed into(22)μ3+Aμ2+Bμ+C=0,where(23)A=Duu∗νq+Dww∗νq+Dvv∗νq+d+δ+c+∂gu∗,v∗∂u,B=Duu∗νq+d+∂gu∗,v∗∂uDww∗νq+Dvv∗νq+δ+c+Dww∗νq+δDvv∗νq+c−p∂gu∗,v∗∂v,C=Duu∗νq+d+∂gu∗,v∗∂uDww∗νq+δDvv∗νq+c−p∂gu∗,v∗∂v+p∂gu∗,v∗∂u∂gu∗,v∗∂v.
It is easy to see thatA>0. From the condition H1, together with gu∗,v∗=δw∗ and pw∗=cv∗, we get(24)p∂gu∗,v∗∂v≤pgu∗,v∗v∗=pδw∗pw∗/c=cδ,which yields B>0 and C>0. By direct computation, we have that(25)AB−C=Duu∗νq+Dww∗νq+Dvv∗νq+d+δ+∂gu∗,v∗∂uDuu∗νq+d+∂gu∗,v∗∂uDww∗νq+Dvv∗νq+δ+c+cDuu∗νq+dc+Dww∗νq+Dvv∗νq+δ+c∂gu∗,v∗∂uDww∗νq+Dvv∗νq+c+∂gu∗,v∗∂ucδ−p∂gu∗,v∗∂v>0.
So, we can see that all eigenvalues of (22) have no positive real parts by the Routh–Hurwitz criterion. Therefore, E∗ is locally asymptotically stable for ℛ0>1. Thus, we obtain the following results.Theorem 3.
The infection equilibriumE∗ of model (1) is locally asymptotically stable if ℛ0>1.


### 2.3.2. The Global Stability of Model (1)

In this section, by constructing appropriate Lyapunov functionals, we continue to investigate the global stability of model (1). To this end, we need the following lemma.Lemma 2.
Under homogeneous Neumann boundary condition, for a density-dependent diffusion problem, we can obtain the following result:(26)∫Ω1gu∇⋅Du∇udx=∫ΩDug′u∇u2g2udx,where Ω is a bounded domain in Rn and its boundary ∂Ω is smooth.Proof.
According to the homogeneous Neumann boundary condition, we get(27)∫Ω1gu−∇⋅Du∇udx=∫Ω1guD′u∇u2dx+∫Ω1guDuΔudx=∫Ω1guD′u∇u2dx+∫∂Ω1guDu∂u∂νds−∫ΩDu∇1gu⋅∇udx−∫Ω1guD′u∇u2dx=∫ΩDug′u∇u2g2udx.
This completes the proof of Lemma2.Remark 1.
Choosinggu=1, u, and 1/u in Lemma 2, respectively, we can obtain that(i)
∫Ω∇⋅Du∇udx=0;(ii)
∫Ω1/u∇⋅Du∇udx=∫ΩDu∇u2/u2dx;(iii)
∫Ωu∇⋅Du∇udx=−∫ΩDu∇u2dx;which implies that Lemma 2 generalized the results in [6, 7].
Firstly, we discuss the global asymptotic stability of the infection-free equilibriumE0 of model (3).Theorem 4.
Ifℛ0≤1, then the infection-free equilibrium E0 of model (3) is globally asymptotically stable.Proof.
Define the following Lyapunov functional:(28)Vt=∫Ωwx,t+δpvx,tdx.
Through calculating the derivative ofVt about time, we can get(29)dVtdt=∫Ω∂w∂t+δp∂v∂tdx=∫Ω∇⋅Dw∇w+gux,t,vx,t−δwx,t+δp∇⋅Dv∇v+δwx,t−cδpvx,tdx=∫Ωgux,t,vx,t−cδpvx,tdx+∫Ω∇⋅Dw∇wdx+δp∫Ω∇⋅Dv∇vdx.
Note that∂u/∂t≤λ−du, then we have limt⟶∞ux,t≤λ/d, which means that all omega limit points satisfy u≤λ/d. So, we need to only consider solutions for which u≤λ/d. Since ∂g/∂u≥0, we can conclude the inequality gu,v≤gλ/d,v. Under the condition H1, together with the fact that gλ/d,v/v is decreasing in v, we can obtain(30)gux,t,vx,t−cδpvx,t=cδpvpcδgu,vv−1≤cδpvpcδgλ/d,vv−1≤cδpvlimv⟶0+pcδgλ/d,vv−1=cδpvpcδ∂gλ/d,0∂v−1=cδpvℛ0−1.
Meanwhile, from Remark1, we get(31)∫Ω∇⋅Dw∇wdx=0=δp∫Ω∇⋅Dv∇vdx.
Thus, we obtain(32)dVtdt≤cδpℛ0−1∫Ωvdx.
So, we have thatdVt/dt≤0 when ℛ0≤1. Further, it is not hard to believe that dVt/dt≤0 leads to either v=0 or ℛ0=1 and u=λ/d. It is easy to confirm that the singleton E0 is just the largest compact invariant set in u,w,vdVt/dt≤0. When ℛ0≤1, we can know that the infection-free equilibrium E0 of model (3) is globally asymptotically stable by the LaSalle’s invariance principle. The proof of Theorem 4 is complete.
Next, we investigate the global asymptotic stability problem of the infection equilibriumE∗ of model (3). To prove the following theorem, we need the following hypothesis:(33)H3vv∗−gu,vgu,v∗gu,vgu,v∗−1≥0,for allu,v>0.Theorem 5.
SupposeH3 is satisfied. If ℛ0>1, then the infection equilibrium E∗ of model (3) is globally asymptotically stable.Proof.
LetHz=z−1−lnz,z>0, and define(34)Vt=∫Ωu−u∗−∫u∗ux,tgu∗,v∗gs,v∗ds+w∗Hww∗+δpv∗Hvv∗dx,which is called a Lyapunov functional. Obviously, 0 is the minimum value of function Hz at z=1. This leads to Hz≥0. Meanwhile, 0 is also the global minimum of function G:u⟶u−u∗−∫u∗ugu∗,v∗/gs,v∗ds at u=u∗, which implies Gu≥0 for any x>0. Therefore, for all x∈Ω and t≥0, Vt≥0 with equality holding if and only if ux,t/u∗=w/w∗=vx,t/v∗=1.
Through calculating the derivative ofVt about time, we can get(35)dVtdt=∫Ω1−gu∗,v∗gu,v∗∂u∂t+1−w∗w∂w∂t+δp1−v∗v∂v∂tdx=∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇u+λ−dux,t−gux,t,vx,t+1−w∗w∇⋅Dw∇w+gux,t,vx,t−δwx,t+δp1−v∗v∇⋅Dv∇v+pwx,t−cvx,tdx.
Note thatλ=du∗+gu∗,v∗, δw∗=gu∗,v∗, and cv∗=pw∗. So, we have(36)dVtdt=∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇u+du∗−dux,t+gu∗,v∗−gux,t,vx,t+1−w∗w∇⋅Dw∇w+gux,t,vx,t−wx,tw∗gu∗,v∗+δp1−v∗v∇⋅Dv∇v+pwx,t−vx,tv∗w∗dx=∫Ωdu∗1−uu∗1−gu∗,v∗gu,v∗+gu∗,v∗1−gu∗,v∗gu,v∗1−gux,t,vx,tgu∗,v∗+1−w∗wgux,t,vx,tgu∗,v∗−wx,tw∗+1−v∗vwx,tw∗−vx,tv∗dx+∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇udx+∫Ω1−w∗w∇⋅Dw∇wdx+δp∫Ω1−v∗v∇⋅Dv∇vdx=∫Ωdu∗1−uu∗1−gu∗,v∗gu,v∗+gu∗,v∗Hgux,t,vx,tgux,t,v∗−Hgu∗,v∗gux,t,v∗−Hgux,t,vx,tgu∗,v∗+Hgux,t,vx,tgu∗,v∗−Hgux,t,vx,tw∗gu∗,v∗wx,tx,t−Hwx,tw∗+Hwx,tw∗−Hwx,tv∗vx,tw∗−Hvx,tv∗dx+∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇udx+∫Ω1−w∗w∇⋅Dw∇wdx+δp∫Ω1−v∗v∇⋅Dv∇vdx=∫Ωdu∗1−uu∗1−gu∗,v∗gu,v∗+gu∗,v∗Hgux,t,vx,tgux,t,v∗−Hgu∗,v∗gux,t,v∗−Hgux,t,vx,tw∗gu∗,v∗wx,t−Hwx,tv∗vx,tw∗−Hvx,tv∗dx+∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇udx+∫Ω1−w∗w∇⋅Dw∇wdx+δp∫Ω1−v∗v∇⋅Dv∇vdx.
From the hypothesisH1, we can know that gu,v∗ is increasing in u strictly. Hence, we have(37)1−uu∗1−gu∗,v∗gu,v∗≤0,and 1−u/u∗1−gu∗,v∗/gu,v∗ equals 0 only if u=u∗. Based on the assumption H3, we can conclude that gux,t,vx,t/gux,t,v∗ lies between 1 and vx,t/v∗. Thus, by the fact that Hz is a monotone function on each side of 1, we can obtain that Hgux,t,vx,t/gux,t,v∗≤Hvx,t/v∗. Therefore,(38)dVtdt≤−gu∗,v∗∫ΩHgu∗,v∗gu,v∗+Hwx,tv∗vx,tw∗+Hgux,t,vx,tw∗gu∗,v∗wx,tdx+∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇udx+∫Ω1−w∗w∇⋅Dw∇wdx+δp∫Ω1−v∗v∇⋅Dv∇vdx.
From Lemma2 and Remark 1, it is easy to conclude that(39)dVtdt≤0.
Hence,E∗ is stable and dVt/dt=0 iff u=u∗, w=w∗, and v=v∗. It is easy to confirm that the singleton E∗ is just the largest compact invariant set in u,w,vdVt/dt=0. According to the LaSalle’s invariance principle, we can know that the infection equilibrium E∗ of model (3) is globally asymptotically stable when ℛ0>1. This completes the proof of Theorem 5.


## 2.3.1. The Local Stability of Model (1)

In this section, we shall investigate the local stability of model (1). Under the homogeneous Neumann boundary conditions, assume that 0=ν0<ν1<⋅⋅⋅<νq<⋅⋅⋅ be the eigenvalues of −Δ on Ω and for q=0,1,2,…, and Eνq be the eigenfunctions space corresponding to νq in C1Ω. Suppose that Eνq has an orthonormal basis φqm:m=1,2,…,dimEνq, X=C1Ω3, and Xqm=cφqm:c∈ℝ3. Thus, we obtain(11)X=⊕q=1∞Xq,Xq=⊕m=1dimEνqXqm.Linearizing system (1) at the equilibrium u˜,w˜,v˜ is as follows:(12)∂ux,t∂t=Duu˜Δu−dux,t−∂gu˜,v˜∂uux,t−∂gu˜,v˜∂vvx,t,∂wx,t∂t=Dww˜Δw+∂gu˜,v˜∂uux,t+∂gu˜,v˜∂vvx,t−δwx,t,∂vx,t∂t=Dvv˜Δv+pwx,t−cvx,t.Note that, under the linearization,Xq is an invariant for each q=0,1,2,…. We let the expression of solution be Qx,t=eμtψx, where ψx∈Xq satisfies −Δψ=νqψ. Thus, the characteristic equation of (12) is(13)μ+Duu˜νq+d+∂gu˜,v˜∂uμ+Dww˜νq+δμ+Dvv˜νq+c−p∂gu˜,v˜∂v+p∂gu˜,v˜∂u∂gu˜,v˜∂v=0.When the equilibrium of model (1) is E0=λ/d,0,0, the characteristic equation (13) is changed into(14)μ+Duλdνq+dμ+δμ+c−p∂gλ/d,0∂v=0.Obviously, equation (14) has a negative real root σ=−Duλ/dνq−d. We now discuss the transcendental equation as follows:(15)μ+δμ+c−p∂gλ/d,0∂v=0.Set(16)Hμ,q=μ+δμ+c−p∂gλ/d,0∂v.Case 1.
ℛ
0
<
1. By direct computation, we have(17)H0,q=cδ−p∂gλ/d,0∂v=cδ1−ℛ0>0,and for μ≥0,(18)∂Hμ,q∂μ=2μ+δ+c>0.
So, there is only a negative real root for equation (14).
Assume that there be a complex rootρ+iϖ with ρ≥0, then we can conclude that(19)μδ+1>1,μc+1>1,pcδ∂gλ/d,0∂v=ℛ0<1,which implies(20)μδ+1μc+1>ℛ0=pcδ∂gλ/d,0∂v.
This contradicts to equation (15). Then, all roots of equation (14) have no positive real parts yet. Thus, the infection-free equilibrium E0 is locally asymptotically stable if ℛ0<1.Case 2.
ℛ
0
>
1. Obviously, we need to only consider the space X0 corresponding to ν0=0 when q=0. By direct computation, we get(21)H0,0=cδ1−ℛ0<0,limμ⟶+∞Hμ,0=+∞.
This implies that the equationHμ,0=0 must hold for some μ0>0. Thus, there exists at least one positive root of equation (14). So, the infection-free equilibrium E0 is not stable for ℛ0>1. Then, we have the following results.Theorem 2.
The infection-free equilibriumE0 of model (1) is locally asymptotically stable if ℛ0<1. However, E0 is unstable if ℛ0>1.
When the equilibrium of model (1) is E∗=u∗,w∗,v∗, the characteristic equation (13) is changed into(22)μ3+Aμ2+Bμ+C=0,where(23)A=Duu∗νq+Dww∗νq+Dvv∗νq+d+δ+c+∂gu∗,v∗∂u,B=Duu∗νq+d+∂gu∗,v∗∂uDww∗νq+Dvv∗νq+δ+c+Dww∗νq+δDvv∗νq+c−p∂gu∗,v∗∂v,C=Duu∗νq+d+∂gu∗,v∗∂uDww∗νq+δDvv∗νq+c−p∂gu∗,v∗∂v+p∂gu∗,v∗∂u∂gu∗,v∗∂v.
It is easy to see thatA>0. From the condition H1, together with gu∗,v∗=δw∗ and pw∗=cv∗, we get(24)p∂gu∗,v∗∂v≤pgu∗,v∗v∗=pδw∗pw∗/c=cδ,which yields B>0 and C>0. By direct computation, we have that(25)AB−C=Duu∗νq+Dww∗νq+Dvv∗νq+d+δ+∂gu∗,v∗∂uDuu∗νq+d+∂gu∗,v∗∂uDww∗νq+Dvv∗νq+δ+c+cDuu∗νq+dc+Dww∗νq+Dvv∗νq+δ+c∂gu∗,v∗∂uDww∗νq+Dvv∗νq+c+∂gu∗,v∗∂ucδ−p∂gu∗,v∗∂v>0.
So, we can see that all eigenvalues of (22) have no positive real parts by the Routh–Hurwitz criterion. Therefore, E∗ is locally asymptotically stable for ℛ0>1. Thus, we obtain the following results.Theorem 3.
The infection equilibriumE∗ of model (1) is locally asymptotically stable if ℛ0>1.


## 2.3.2. The Global Stability of Model (1)

In this section, by constructing appropriate Lyapunov functionals, we continue to investigate the global stability of model (1). To this end, we need the following lemma.Lemma 2.
Under homogeneous Neumann boundary condition, for a density-dependent diffusion problem, we can obtain the following result:(26)∫Ω1gu∇⋅Du∇udx=∫ΩDug′u∇u2g2udx,where Ω is a bounded domain in Rn and its boundary ∂Ω is smooth.Proof.
According to the homogeneous Neumann boundary condition, we get(27)∫Ω1gu−∇⋅Du∇udx=∫Ω1guD′u∇u2dx+∫Ω1guDuΔudx=∫Ω1guD′u∇u2dx+∫∂Ω1guDu∂u∂νds−∫ΩDu∇1gu⋅∇udx−∫Ω1guD′u∇u2dx=∫ΩDug′u∇u2g2udx.
This completes the proof of Lemma2.Remark 1.
Choosinggu=1, u, and 1/u in Lemma 2, respectively, we can obtain that(i)
∫Ω∇⋅Du∇udx=0;(ii)
∫Ω1/u∇⋅Du∇udx=∫ΩDu∇u2/u2dx;(iii)
∫Ωu∇⋅Du∇udx=−∫ΩDu∇u2dx;which implies that Lemma 2 generalized the results in [6, 7].
Firstly, we discuss the global asymptotic stability of the infection-free equilibriumE0 of model (3).Theorem 4.
Ifℛ0≤1, then the infection-free equilibrium E0 of model (3) is globally asymptotically stable.Proof.
Define the following Lyapunov functional:(28)Vt=∫Ωwx,t+δpvx,tdx.
Through calculating the derivative ofVt about time, we can get(29)dVtdt=∫Ω∂w∂t+δp∂v∂tdx=∫Ω∇⋅Dw∇w+gux,t,vx,t−δwx,t+δp∇⋅Dv∇v+δwx,t−cδpvx,tdx=∫Ωgux,t,vx,t−cδpvx,tdx+∫Ω∇⋅Dw∇wdx+δp∫Ω∇⋅Dv∇vdx.
Note that∂u/∂t≤λ−du, then we have limt⟶∞ux,t≤λ/d, which means that all omega limit points satisfy u≤λ/d. So, we need to only consider solutions for which u≤λ/d. Since ∂g/∂u≥0, we can conclude the inequality gu,v≤gλ/d,v. Under the condition H1, together with the fact that gλ/d,v/v is decreasing in v, we can obtain(30)gux,t,vx,t−cδpvx,t=cδpvpcδgu,vv−1≤cδpvpcδgλ/d,vv−1≤cδpvlimv⟶0+pcδgλ/d,vv−1=cδpvpcδ∂gλ/d,0∂v−1=cδpvℛ0−1.
Meanwhile, from Remark1, we get(31)∫Ω∇⋅Dw∇wdx=0=δp∫Ω∇⋅Dv∇vdx.
Thus, we obtain(32)dVtdt≤cδpℛ0−1∫Ωvdx.
So, we have thatdVt/dt≤0 when ℛ0≤1. Further, it is not hard to believe that dVt/dt≤0 leads to either v=0 or ℛ0=1 and u=λ/d. It is easy to confirm that the singleton E0 is just the largest compact invariant set in u,w,vdVt/dt≤0. When ℛ0≤1, we can know that the infection-free equilibrium E0 of model (3) is globally asymptotically stable by the LaSalle’s invariance principle. The proof of Theorem 4 is complete.
Next, we investigate the global asymptotic stability problem of the infection equilibriumE∗ of model (3). To prove the following theorem, we need the following hypothesis:(33)H3vv∗−gu,vgu,v∗gu,vgu,v∗−1≥0,for allu,v>0.Theorem 5.
SupposeH3 is satisfied. If ℛ0>1, then the infection equilibrium E∗ of model (3) is globally asymptotically stable.Proof.
LetHz=z−1−lnz,z>0, and define(34)Vt=∫Ωu−u∗−∫u∗ux,tgu∗,v∗gs,v∗ds+w∗Hww∗+δpv∗Hvv∗dx,which is called a Lyapunov functional. Obviously, 0 is the minimum value of function Hz at z=1. This leads to Hz≥0. Meanwhile, 0 is also the global minimum of function G:u⟶u−u∗−∫u∗ugu∗,v∗/gs,v∗ds at u=u∗, which implies Gu≥0 for any x>0. Therefore, for all x∈Ω and t≥0, Vt≥0 with equality holding if and only if ux,t/u∗=w/w∗=vx,t/v∗=1.
Through calculating the derivative ofVt about time, we can get(35)dVtdt=∫Ω1−gu∗,v∗gu,v∗∂u∂t+1−w∗w∂w∂t+δp1−v∗v∂v∂tdx=∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇u+λ−dux,t−gux,t,vx,t+1−w∗w∇⋅Dw∇w+gux,t,vx,t−δwx,t+δp1−v∗v∇⋅Dv∇v+pwx,t−cvx,tdx.
Note thatλ=du∗+gu∗,v∗, δw∗=gu∗,v∗, and cv∗=pw∗. So, we have(36)dVtdt=∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇u+du∗−dux,t+gu∗,v∗−gux,t,vx,t+1−w∗w∇⋅Dw∇w+gux,t,vx,t−wx,tw∗gu∗,v∗+δp1−v∗v∇⋅Dv∇v+pwx,t−vx,tv∗w∗dx=∫Ωdu∗1−uu∗1−gu∗,v∗gu,v∗+gu∗,v∗1−gu∗,v∗gu,v∗1−gux,t,vx,tgu∗,v∗+1−w∗wgux,t,vx,tgu∗,v∗−wx,tw∗+1−v∗vwx,tw∗−vx,tv∗dx+∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇udx+∫Ω1−w∗w∇⋅Dw∇wdx+δp∫Ω1−v∗v∇⋅Dv∇vdx=∫Ωdu∗1−uu∗1−gu∗,v∗gu,v∗+gu∗,v∗Hgux,t,vx,tgux,t,v∗−Hgu∗,v∗gux,t,v∗−Hgux,t,vx,tgu∗,v∗+Hgux,t,vx,tgu∗,v∗−Hgux,t,vx,tw∗gu∗,v∗wx,tx,t−Hwx,tw∗+Hwx,tw∗−Hwx,tv∗vx,tw∗−Hvx,tv∗dx+∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇udx+∫Ω1−w∗w∇⋅Dw∇wdx+δp∫Ω1−v∗v∇⋅Dv∇vdx=∫Ωdu∗1−uu∗1−gu∗,v∗gu,v∗+gu∗,v∗Hgux,t,vx,tgux,t,v∗−Hgu∗,v∗gux,t,v∗−Hgux,t,vx,tw∗gu∗,v∗wx,t−Hwx,tv∗vx,tw∗−Hvx,tv∗dx+∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇udx+∫Ω1−w∗w∇⋅Dw∇wdx+δp∫Ω1−v∗v∇⋅Dv∇vdx.
From the hypothesisH1, we can know that gu,v∗ is increasing in u strictly. Hence, we have(37)1−uu∗1−gu∗,v∗gu,v∗≤0,and 1−u/u∗1−gu∗,v∗/gu,v∗ equals 0 only if u=u∗. Based on the assumption H3, we can conclude that gux,t,vx,t/gux,t,v∗ lies between 1 and vx,t/v∗. Thus, by the fact that Hz is a monotone function on each side of 1, we can obtain that Hgux,t,vx,t/gux,t,v∗≤Hvx,t/v∗. Therefore,(38)dVtdt≤−gu∗,v∗∫ΩHgu∗,v∗gu,v∗+Hwx,tv∗vx,tw∗+Hgux,t,vx,tw∗gu∗,v∗wx,tdx+∫Ω1−gu∗,v∗gu,v∗∇⋅Du∇udx+∫Ω1−w∗w∇⋅Dw∇wdx+δp∫Ω1−v∗v∇⋅Dv∇vdx.
From Lemma2 and Remark 1, it is easy to conclude that(39)dVtdt≤0.
Hence,E∗ is stable and dVt/dt=0 iff u=u∗, w=w∗, and v=v∗. It is easy to confirm that the singleton E∗ is just the largest compact invariant set in u,w,vdVt/dt=0. According to the LaSalle’s invariance principle, we can know that the infection equilibrium E∗ of model (3) is globally asymptotically stable when ℛ0>1. This completes the proof of Theorem 5.


## 3. Multistrain Viral Infection Model

In this section, we shall begin to study the multi-strain viral infection model (7) with Neumann boundary conditions:(40)Duu∂u∂ν=0,Dwiwi∂wi∂ν=0,Dvivi∂vi∂ν=0,on∂Ω×0,+∞,and initial conditions(41)ux,0=ϕ1x≥0,wix,0=ϕ2,ix≥0,vix,0=ϕ3,ix≥0,x∈Ω¯.SupposeΩ is a bounded domain in ℝn and its boundary ∂Ω is smooth, then the outward normal derivative on ∂Ω is denoted by ∂ν. The nonnegative functions ϕ1x,0 and ϕr,ix,0 r=2,3,i=1,2,…,n are Hölder continuous and satisfy ∂ϕ1/∂ν=0 and ∂ϕr,i/∂ν=0 on ∂Ω, respectively. Moreover, we can get the steady-state system of model (7) as follows:(42)−∇⋅Duu∇u=λ−dux,t−∑i=1ngiux,t,vix,t,−∇⋅Dwiwi∇wi=giux,t,vix,t−δiwix,t,−∇⋅Dvivi∇vi=piwix,t−civix,t,i=1,2,…,n,where Duu=um1, Dwiwi=wim2,i, and Dvivi=vim3,i are the density-dependent diffusion coefficients for m1>0 and mr,i>0, r=2,3,i=1,2,…,n.


### 3.1. Well Posedness

In this section, we give out some definitions and results first, which can help us to prove that the solutions to problems (7) and (42) exist and are unique and bounded. Similar to the definition method in [6, 7], we can definite that two pairs of functions U˜=u˜,w˜i,v˜i, U^=u^,w^i,v^i in CΩ¯×0,T∩C1,2Ω×0,T and U˜s=u˜s,w˜i,s,v˜i,s, U^s=u^s,w^i,s,v^i,s in CΩ¯×0,T∩C1,2Ω×0,T are called coupled upper and lower solutions to system (7) and to system (42), respectively. Then, we can easily conclude that each pair of coupled upper and lower solutions U˜s and U^s of problems (7) is also ones of (42). In what follows, we shall search a pair of constant solutions of (42), which is called upper and lower solutions of (42).Lemma 3.
For constantl1>0 and L1>0 and any sufficiently large constant L2,i>0 and L3,i>0 and small constant l2,i>0 and l3,i>0,i=1,2,…,n, the coupled upper and lower solutions to (42) are a constant pair L=L1,L2,i,L3,i and l=l1,l2,i,l3,i.Proof.
Obviously, the constant pairL=L1,L2,i,L3,i and l=l1,l2,i,l3,i satisfy the homogeneous Neumann boundary conditions (40). Further, by means of the condition H1, we can confirm that the differential inequalities in (42) can also satisfied by a constant pair L and l if(43)λ−dL1−∑i=1ngiL1,l3,i≤λ−dL1≤0,giL1,L3,i−δL2,i≤ηiL1L3,i−δL2,i≤0,piL2,i−ciL3,i≤0;λ−dl1−∑i=1ngil1,L3,i≥λ−dl1−∑i=1nηil1L3,i=λ−d+∑i=1nηiL3,il1≥0,gil1,l3,i−δl2,i≥givl1,l3,il3,i−δl2,i≥0,pil2,i−cil3,i≥0.
Then, the above inequalities can be fulfilled for constantsl1 and L1, sufficiently large L2,i and L3,i, and small l2,i and l3,ii=1,2,…,n if(44)λd≤L1,ληiL3,idδi≤L2,i≤cL3,ip,l1≤λd+∑i=0nηiL3,i,l2,i≤givl1,l3,il3,iδi,l3,i≤pil2,ici.
To obtain the upper solution of (42), we only need to take L1=λ/d and sufficiently large L2 and L3 satisfying ληiL3,i/dδi≤L2,i≤cL3,i/p. Meanwhile, choosing sufficiently small l3,i, l2,i=cil3,i/pi, and l1 satisfies ciδi/pi≤givl1,l3,i and l1≤λ/d+∑i=1nηi, and we can also obtain the lower solution of (42). This means that the coupled upper and lower solutions of (42) are a constant pair L and l. This completes the proof of Lemma 3.
According to the results of Lemma3 and method in [6, 7], we can conclude that problems (7) and (42) is well-posed similarly. That is, the solutions to problems (7) and (42) exist and are unique and bounded.


### 3.2. Competitive Exclusion

To describe a single-straini-infected cell which generates the average number of newly infected cells when the disease has just entered the body, we define(45)ℛi=piciδi∂giλ/d,0∂vi,which is called the basic reproduction number of strain i in model (7). Obviously, E0n=λ/d,0,0 is always an infection-free steady state of system (7), where 0 is an n-dimensional vector of zeros. Similar to the analysis and discussion in [27], together with the basic reproduction number ℛi, it is not hard to see that for each strain i, model (7) has only a corresponding single-strain steady state as follows:(46)Ei=u∗,0,…,0,wi∗,0,…,0,vi∗,0,…,0,iff ℛi>1 and it satisfies(47)λ−dux,t−giux,t,vix,t=0,giux,t,vix,t−δiwix,t=0,piwix,t−civix,t=0..Here, we have to point out that the nonzero componentswi∗ is in position i+1 and vi∗ are in position n+i+1.Now, we begin to investigate the competitive exclusion problem of model (7). From the results in [36], to discuss the competitive exclusion problem, we only need to discuss the globally asymptotical stability of steady state E1 of strain one. Next, we shall verify that the globally asymptotical stability of the infected steady state E1 of strain one by constructing an appropriate Lyapunov functional. Here, we also assume that for i=1,2,…,n, giu,vi satisfies the condition H1 and ηk<η1, k=2,…,n. Moreover, we need the following condition:(48)H4ηkpkckδk≤dλp1c1δ1minu∈0,λ/d∂g1u,v1∗∂v1.Theorem 6.
SupposeH3 and H4 are satisfied. For k=2,…,n, let ℛ1>1 and ℛ1>ℛk hold. Then, the single-strain steady state E1 of model (7) is globally asymptotically stable.Proof.
LetHz=z−1−lnz,z>0 and define(49)Vt=∫Ωu−u1∗−∫u1∗ux,tg1u1∗,v1∗g1s,v1∗ds+w1∗Hw1w1∗+δ1p1v1∗Hv1v1∗+∑k=2nwk+δkpkvkdx,which is call a Lyapunov functional. Obviously, 0 is the minimum value of function Hz at z=1. This leads to Hz≥0. Meanwhile, the 0 is also the global minimum of function G:u⟶u−u1∗−∫u1∗ug1u1∗,v1∗/g1s,v1∗ds at u=u1∗, which implies Gu≥0 for any x>0. Therefore, Vt≥0 with equality holding if and only ux,t/u1∗=w1x,t/w1∗=v1x,t/v1∗=1 and wk=vk=0k=2,…,n for all x∈Ω and t≥0.
Through calculating the derivative ofVt about time, we get(50)dVtdt=∫Ω1−g1u1∗,v1∗g1u,v1∗∂u∂t+1−w1∗w1∂w1∂t+δ1p11−v1∗v1∂v1∂t+∑k=2n∂wk∂t+δkpk∂vk∂tdx=∫Ω1−g1u1∗,v1∗g1u,v1∗∇⋅Du∇u+λ−dux,t−g1ux,t,v1x,t+1−w1∗w1∇⋅Dw1∇w1+g1ux,t,v1x,t−δ1w1x,t+δ1p11−v1∗v1∇⋅Dv1∇v1+p1w1x,t−c1v1x,t−1−g1u1∗,v1∗g1u,v1∗∑k=2ngkux,t,vkx,t+∑k=2n∇⋅Dwk∇wk+gkux,t,vkx,t−δkwkx,t∑k=2nδkpk∇⋅Dvk∇vk+pkwkx,t−ckvkx,tdx.
Note thatλ=du1∗+g1u1∗,v1∗, δ1w1∗=g1u1∗,v1∗, and c1v1∗=p1w1∗. So, we have(51)dVtdt=∫Ω1−g1u1∗,v1∗g1u,v1∗∇⋅Du∇u+du1∗−dux,t+g1u1∗,v1∗−g1ux,t,v1x,t+1−w1∗w1∇⋅Dw1∇w1+g1ux,t,v1x,t−w1x,tw1∗g1u1∗,v1∗+δ1p11−v1∗v1∇⋅Dv1∇v1+p1w1x,t−v1x,tv1∗w1∗+∑k=2n∇⋅Dwk∇wk+g1u1∗,v1∗g1u,v1∗gkux,t,vkx,t−δkwkx,t+∑k=2nδkpk∇⋅Dvk∇vk+pkwkx,t−ckvkx,tdx=∫Ωdu1∗1−uu1∗1−g1u1∗,v1∗g1u,v∗+g1u1∗,v1∗1−g1u1∗,v1∗g1u,v1∗1−g1ux,t,v1x,tg1u1∗,v1∗+1−w1∗w1g1ux,t,v1x,tgu1∗,v1∗−w1x,tw1∗+1−v1∗v1w1x,tw1∗−v1x,tv1∗dx+∑k=2n∫Ωg1u1∗,v1∗g1u,v1∗gkux,t,vkx,t−ckδkpkvkx,tdx+∫Ω1−g1u1∗,v1∗g1u,v1∗∇⋅Du∇udx+∫Ω1−w1∗w1∇⋅Dw1∇w1dx+δ1p1∫Ω1−w1∗w1∇⋅Dv1∇v1dx+∑k=2n∫Ω∇⋅Dwk∇wkdx+∑k=2nδkpk∫Ω∇⋅Dvk∇vkdx=∫Ωdu1∗1−uu1∗1−g1u1∗,v1∗g1u,v1∗+g1u1∗,v1∗Hg1ux,t,v1x,tg1ux,t,v1∗−Hg1u1∗,v1∗g1ux,t,v1∗−Hg1ux,t,v1x,tg1u1∗,v1∗+Hg1ux,t,v1x,tg1u1∗,v1∗−Hg1ux,t,v1x,tw1∗g1u1∗,v1∗wx,t−Hw1x,tw1∗+Hw1x,tw1∗−Hw1x,tv1∗v1x,tw1∗−Hv1x,tv1∗dx+∑k=2n∫Ωg1u1∗,v1∗g1u,v1∗gkux,t,vkx,t−ckδkpkvkx,tdx+∫Ω1−g1u1∗,v1∗g1u,v1∗∇⋅Du∇udx+∫Ω1−w1∗w1∇⋅Dw1∇wdx+δ1p1∫Ω1−v1∗v1∇⋅Dv1∇v1dx+∑k=2n∫Ω∇⋅Dwk∇wkdx+∑k=2nδkpk∫Ω∇⋅Dvk∇vkdx=∫Ωdu1∗1−uu1∗1−g1u1∗,v1∗g1u,v1∗+g1u1∗,v1∗Hg1ux,t,v1x,tg1ux,t,v1∗,
From the hypothesisH1, we can know the function g1u,v1∗ is strictly increasing in u. Hence, we have(52)1−uu1∗1−g1u1∗,v1∗g1u,v1∗≤0,with equality only if u=u1∗. Based on the assumption H3, we can conclude that g1ux,t,v1x,t/g1ux,t,v1∗ lies between 1 and v1x,t/v1∗. Thus, by the fact that Hz is a monotone on each side of 1, we can obtain that Hg1ux,t,v1x,t/g1ux,t,v1∗≤Hv1x,t/v1∗. Meanwhile, for k=2,…,n, we have(53)g1u1∗,v1∗g1u,v1∗gkux,t,vkx,t−ckδkpkvkx,t=ckδkg1u1∗,v1∗pkg1u,v1∗vkpkckδkgku,vkvk−g1u,v1∗g1u1∗,v1∗=ckδkg1u1∗,v1∗pkg1u,v1∗vkgku,vkvk∂gkλ/d,0/∂vkℛk−g1u,v1∗v1∗∂g1λ/d,0/∂v1ℛ1=ckδkg1u1∗,v1∗pkg1u,v1∗vkαkℛk−α1ℛ1,where αk=gku,vk/vk∂gkλ/d,0/∂vk and α1=g1u,v1∗/v1∗∂g1λ/d,0/∂v1.
Under the conditionH1 and H4, we have(54)α1αk=g1u,v1∗/v1∗gku,vk/vk∂gkλ/d,0/∂vk∂g1λ/d,0/∂v1≥∂g1u,v1∗/∂v1gku,vk/vk∂gkλ/d,0/∂vk∂g1λ/d,0/∂v1≥∂g1u,v1∗/∂v1ηkλ/d∂gkλ/d,0/∂vk∂g1λ/d,0/∂v1≥pkckδkc1δ1p1∂gkλ/d,0/∂vk∂g1λ/d,0/∂v1=ℛkℛ1,which yields that(55)∑k=2n∫Ωg1u1∗,v1∗g1u,v1∗gkux,t,vkx,t−ckδkpkvkx,tdx≤0.
Therefore, we obtain(56)dVtdt≤−g1u1∗,v1∗∫ΩHg1u∗,v1∗g1u,v1∗+Hw1x,tv1∗v1x,tw1∗+Hg1ux,t,v1x,tw1∗g1u1∗,v1∗w1x,tdx+∫Ω1−g1u1∗,v1∗g1u,v1∗∇⋅Du∇udx+∫Ω1−w1∗w1∇⋅Dw1∇wdx+δ1p1∫Ω1−v1∗v1∇⋅Dv1∇v1dx+∑k=2n∫Ω∇⋅Dwk∇wkdx+∑k=2nδkpk∫Ω∇⋅Dvk∇vkdx.
It follows from Lemma2 and Remark 1 that(57)dVtdt≤0.
Hence,E1 is stable, and dVt/dt=0 iff u=u1∗, w=w1∗, v=v1∗, and wk=0=vk for k=2,…,n. It is easy to confirm that the singleton E1 is just the largest compact invariant set in u,w,vdVt/dt=0. When ℛ1>1 and ℛ1>ℛk hold for k=2,…,n, we can know that the infection steady state E1 of model (7) is globally asymptotically stable by the LaSalle’s invariance principle. The proof of Theorem 6 is complete.Theorem 7.
Assumeℛi<1 for i=1,…,n. Then, the steady state E0n of model (7) is globally asymptotically stable.Proof.
Define(58)Vt=∫Ω∑i=1nwix,t+δipivix,tdx,which is called the Lyapunov functional. The rest of the proof is similar to Theorem 4. So, we omit it. This completes the proof of Theorem 7.Remark 2.
From Theorem7, when the basic reproductive number of each strain is less than 1, we can see that all viral strains die out.


## 3.1. Well Posedness

In this section, we give out some definitions and results first, which can help us to prove that the solutions to problems (7) and (42) exist and are unique and bounded. Similar to the definition method in [6, 7], we can definite that two pairs of functions U˜=u˜,w˜i,v˜i, U^=u^,w^i,v^i in CΩ¯×0,T∩C1,2Ω×0,T and U˜s=u˜s,w˜i,s,v˜i,s, U^s=u^s,w^i,s,v^i,s in CΩ¯×0,T∩C1,2Ω×0,T are called coupled upper and lower solutions to system (7) and to system (42), respectively. Then, we can easily conclude that each pair of coupled upper and lower solutions U˜s and U^s of problems (7) is also ones of (42). In what follows, we shall search a pair of constant solutions of (42), which is called upper and lower solutions of (42).Lemma 3.
For constantl1>0 and L1>0 and any sufficiently large constant L2,i>0 and L3,i>0 and small constant l2,i>0 and l3,i>0,i=1,2,…,n, the coupled upper and lower solutions to (42) are a constant pair L=L1,L2,i,L3,i and l=l1,l2,i,l3,i.Proof.
Obviously, the constant pairL=L1,L2,i,L3,i and l=l1,l2,i,l3,i satisfy the homogeneous Neumann boundary conditions (40). Further, by means of the condition H1, we can confirm that the differential inequalities in (42) can also satisfied by a constant pair L and l if(43)λ−dL1−∑i=1ngiL1,l3,i≤λ−dL1≤0,giL1,L3,i−δL2,i≤ηiL1L3,i−δL2,i≤0,piL2,i−ciL3,i≤0;λ−dl1−∑i=1ngil1,L3,i≥λ−dl1−∑i=1nηil1L3,i=λ−d+∑i=1nηiL3,il1≥0,gil1,l3,i−δl2,i≥givl1,l3,il3,i−δl2,i≥0,pil2,i−cil3,i≥0.
Then, the above inequalities can be fulfilled for constantsl1 and L1, sufficiently large L2,i and L3,i, and small l2,i and l3,ii=1,2,…,n if(44)λd≤L1,ληiL3,idδi≤L2,i≤cL3,ip,l1≤λd+∑i=0nηiL3,i,l2,i≤givl1,l3,il3,iδi,l3,i≤pil2,ici.
To obtain the upper solution of (42), we only need to take L1=λ/d and sufficiently large L2 and L3 satisfying ληiL3,i/dδi≤L2,i≤cL3,i/p. Meanwhile, choosing sufficiently small l3,i, l2,i=cil3,i/pi, and l1 satisfies ciδi/pi≤givl1,l3,i and l1≤λ/d+∑i=1nηi, and we can also obtain the lower solution of (42). This means that the coupled upper and lower solutions of (42) are a constant pair L and l. This completes the proof of Lemma 3.
According to the results of Lemma3 and method in [6, 7], we can conclude that problems (7) and (42) is well-posed similarly. That is, the solutions to problems (7) and (42) exist and are unique and bounded.


## 3.2. Competitive Exclusion

To describe a single-straini-infected cell which generates the average number of newly infected cells when the disease has just entered the body, we define(45)ℛi=piciδi∂giλ/d,0∂vi,which is called the basic reproduction number of strain i in model (7). Obviously, E0n=λ/d,0,0 is always an infection-free steady state of system (7), where 0 is an n-dimensional vector of zeros. Similar to the analysis and discussion in [27], together with the basic reproduction number ℛi, it is not hard to see that for each strain i, model (7) has only a corresponding single-strain steady state as follows:(46)Ei=u∗,0,…,0,wi∗,0,…,0,vi∗,0,…,0,iff ℛi>1 and it satisfies(47)λ−dux,t−giux,t,vix,t=0,giux,t,vix,t−δiwix,t=0,piwix,t−civix,t=0..Here, we have to point out that the nonzero componentswi∗ is in position i+1 and vi∗ are in position n+i+1.Now, we begin to investigate the competitive exclusion problem of model (7). From the results in [36], to discuss the competitive exclusion problem, we only need to discuss the globally asymptotical stability of steady state E1 of strain one. Next, we shall verify that the globally asymptotical stability of the infected steady state E1 of strain one by constructing an appropriate Lyapunov functional. Here, we also assume that for i=1,2,…,n, giu,vi satisfies the condition H1 and ηk<η1, k=2,…,n. Moreover, we need the following condition:(48)H4ηkpkckδk≤dλp1c1δ1minu∈0,λ/d∂g1u,v1∗∂v1.Theorem 6.
SupposeH3 and H4 are satisfied. For k=2,…,n, let ℛ1>1 and ℛ1>ℛk hold. Then, the single-strain steady state E1 of model (7) is globally asymptotically stable.Proof.
LetHz=z−1−lnz,z>0 and define(49)Vt=∫Ωu−u1∗−∫u1∗ux,tg1u1∗,v1∗g1s,v1∗ds+w1∗Hw1w1∗+δ1p1v1∗Hv1v1∗+∑k=2nwk+δkpkvkdx,which is call a Lyapunov functional. Obviously, 0 is the minimum value of function Hz at z=1. This leads to Hz≥0. Meanwhile, the 0 is also the global minimum of function G:u⟶u−u1∗−∫u1∗ug1u1∗,v1∗/g1s,v1∗ds at u=u1∗, which implies Gu≥0 for any x>0. Therefore, Vt≥0 with equality holding if and only ux,t/u1∗=w1x,t/w1∗=v1x,t/v1∗=1 and wk=vk=0k=2,…,n for all x∈Ω and t≥0.
Through calculating the derivative ofVt about time, we get(50)dVtdt=∫Ω1−g1u1∗,v1∗g1u,v1∗∂u∂t+1−w1∗w1∂w1∂t+δ1p11−v1∗v1∂v1∂t+∑k=2n∂wk∂t+δkpk∂vk∂tdx=∫Ω1−g1u1∗,v1∗g1u,v1∗∇⋅Du∇u+λ−dux,t−g1ux,t,v1x,t+1−w1∗w1∇⋅Dw1∇w1+g1ux,t,v1x,t−δ1w1x,t+δ1p11−v1∗v1∇⋅Dv1∇v1+p1w1x,t−c1v1x,t−1−g1u1∗,v1∗g1u,v1∗∑k=2ngkux,t,vkx,t+∑k=2n∇⋅Dwk∇wk+gkux,t,vkx,t−δkwkx,t∑k=2nδkpk∇⋅Dvk∇vk+pkwkx,t−ckvkx,tdx.
Note thatλ=du1∗+g1u1∗,v1∗, δ1w1∗=g1u1∗,v1∗, and c1v1∗=p1w1∗. So, we have(51)dVtdt=∫Ω1−g1u1∗,v1∗g1u,v1∗∇⋅Du∇u+du1∗−dux,t+g1u1∗,v1∗−g1ux,t,v1x,t+1−w1∗w1∇⋅Dw1∇w1+g1ux,t,v1x,t−w1x,tw1∗g1u1∗,v1∗+δ1p11−v1∗v1∇⋅Dv1∇v1+p1w1x,t−v1x,tv1∗w1∗+∑k=2n∇⋅Dwk∇wk+g1u1∗,v1∗g1u,v1∗gkux,t,vkx,t−δkwkx,t+∑k=2nδkpk∇⋅Dvk∇vk+pkwkx,t−ckvkx,tdx=∫Ωdu1∗1−uu1∗1−g1u1∗,v1∗g1u,v∗+g1u1∗,v1∗1−g1u1∗,v1∗g1u,v1∗1−g1ux,t,v1x,tg1u1∗,v1∗+1−w1∗w1g1ux,t,v1x,tgu1∗,v1∗−w1x,tw1∗+1−v1∗v1w1x,tw1∗−v1x,tv1∗dx+∑k=2n∫Ωg1u1∗,v1∗g1u,v1∗gkux,t,vkx,t−ckδkpkvkx,tdx+∫Ω1−g1u1∗,v1∗g1u,v1∗∇⋅Du∇udx+∫Ω1−w1∗w1∇⋅Dw1∇w1dx+δ1p1∫Ω1−w1∗w1∇⋅Dv1∇v1dx+∑k=2n∫Ω∇⋅Dwk∇wkdx+∑k=2nδkpk∫Ω∇⋅Dvk∇vkdx=∫Ωdu1∗1−uu1∗1−g1u1∗,v1∗g1u,v1∗+g1u1∗,v1∗Hg1ux,t,v1x,tg1ux,t,v1∗−Hg1u1∗,v1∗g1ux,t,v1∗−Hg1ux,t,v1x,tg1u1∗,v1∗+Hg1ux,t,v1x,tg1u1∗,v1∗−Hg1ux,t,v1x,tw1∗g1u1∗,v1∗wx,t−Hw1x,tw1∗+Hw1x,tw1∗−Hw1x,tv1∗v1x,tw1∗−Hv1x,tv1∗dx+∑k=2n∫Ωg1u1∗,v1∗g1u,v1∗gkux,t,vkx,t−ckδkpkvkx,tdx+∫Ω1−g1u1∗,v1∗g1u,v1∗∇⋅Du∇udx+∫Ω1−w1∗w1∇⋅Dw1∇wdx+δ1p1∫Ω1−v1∗v1∇⋅Dv1∇v1dx+∑k=2n∫Ω∇⋅Dwk∇wkdx+∑k=2nδkpk∫Ω∇⋅Dvk∇vkdx=∫Ωdu1∗1−uu1∗1−g1u1∗,v1∗g1u,v1∗+g1u1∗,v1∗Hg1ux,t,v1x,tg1ux,t,v1∗,
From the hypothesisH1, we can know the function g1u,v1∗ is strictly increasing in u. Hence, we have(52)1−uu1∗1−g1u1∗,v1∗g1u,v1∗≤0,with equality only if u=u1∗. Based on the assumption H3, we can conclude that g1ux,t,v1x,t/g1ux,t,v1∗ lies between 1 and v1x,t/v1∗. Thus, by the fact that Hz is a monotone on each side of 1, we can obtain that Hg1ux,t,v1x,t/g1ux,t,v1∗≤Hv1x,t/v1∗. Meanwhile, for k=2,…,n, we have(53)g1u1∗,v1∗g1u,v1∗gkux,t,vkx,t−ckδkpkvkx,t=ckδkg1u1∗,v1∗pkg1u,v1∗vkpkckδkgku,vkvk−g1u,v1∗g1u1∗,v1∗=ckδkg1u1∗,v1∗pkg1u,v1∗vkgku,vkvk∂gkλ/d,0/∂vkℛk−g1u,v1∗v1∗∂g1λ/d,0/∂v1ℛ1=ckδkg1u1∗,v1∗pkg1u,v1∗vkαkℛk−α1ℛ1,where αk=gku,vk/vk∂gkλ/d,0/∂vk and α1=g1u,v1∗/v1∗∂g1λ/d,0/∂v1.
Under the conditionH1 and H4, we have(54)α1αk=g1u,v1∗/v1∗gku,vk/vk∂gkλ/d,0/∂vk∂g1λ/d,0/∂v1≥∂g1u,v1∗/∂v1gku,vk/vk∂gkλ/d,0/∂vk∂g1λ/d,0/∂v1≥∂g1u,v1∗/∂v1ηkλ/d∂gkλ/d,0/∂vk∂g1λ/d,0/∂v1≥pkckδkc1δ1p1∂gkλ/d,0/∂vk∂g1λ/d,0/∂v1=ℛkℛ1,which yields that(55)∑k=2n∫Ωg1u1∗,v1∗g1u,v1∗gkux,t,vkx,t−ckδkpkvkx,tdx≤0.
Therefore, we obtain(56)dVtdt≤−g1u1∗,v1∗∫ΩHg1u∗,v1∗g1u,v1∗+Hw1x,tv1∗v1x,tw1∗+Hg1ux,t,v1x,tw1∗g1u1∗,v1∗w1x,tdx+∫Ω1−g1u1∗,v1∗g1u,v1∗∇⋅Du∇udx+∫Ω1−w1∗w1∇⋅Dw1∇wdx+δ1p1∫Ω1−v1∗v1∇⋅Dv1∇v1dx+∑k=2n∫Ω∇⋅Dwk∇wkdx+∑k=2nδkpk∫Ω∇⋅Dvk∇vkdx.
It follows from Lemma2 and Remark 1 that(57)dVtdt≤0.
Hence,E1 is stable, and dVt/dt=0 iff u=u1∗, w=w1∗, v=v1∗, and wk=0=vk for k=2,…,n. It is easy to confirm that the singleton E1 is just the largest compact invariant set in u,w,vdVt/dt=0. When ℛ1>1 and ℛ1>ℛk hold for k=2,…,n, we can know that the infection steady state E1 of model (7) is globally asymptotically stable by the LaSalle’s invariance principle. The proof of Theorem 6 is complete.Theorem 7.
Assumeℛi<1 for i=1,…,n. Then, the steady state E0n of model (7) is globally asymptotically stable.Proof.
Define(58)Vt=∫Ω∑i=1nwix,t+δipivix,tdx,which is called the Lyapunov functional. The rest of the proof is similar to Theorem 4. So, we omit it. This completes the proof of Theorem 7.Remark 2.
From Theorem7, when the basic reproductive number of each strain is less than 1, we can see that all viral strains die out.


## 4. Numerical Calculations and Simulations

In this section, we perform some numerical simulations that illustrate and supplement the analytic results given in the previous sections. In what follows, we perform the numerical simulations from the two aspects of global threshold dynamics and competitive exclusion, respectively.Firstly, we carry out the numerical simulation for global threshold dynamics. Taking the incidence function asgux,t,vx,t=βux,tvx,t/1+aux,t+bvx,t+abux,tvx,t, we can change model (1) into the following form:(59)∂ux,t∂t=∇⋅Duu∇u+λ−dux,t−βux,tvx,t1+aux,t+bvx,t+abux,tvx,t,∂wx,t∂t=∇⋅Dww∇w+βux,tvx,t1+aux,t+bvx,t+abux,tvx,t−δwx,t,∂vx,t∂t=∇⋅Dvv∇v+pwx,t−cvx,t,with Neumann boundary conditions (2) and initial conditions (3). From (10), together with direct computation, we get the basic reproduction number ℛ0=pλβ/cδd+aλ for system (59). Meanwhile, the equilibria are given out by E0=λ/d,0,0 and E∗=u∗,w∗,v∗, where w∗=λ−du∗/δ, v∗=pλ−du∗/cδ, and(60)u∗=pabλ+acδ−pβ−pbd+pabλ+acδ−pβ−pbd2+4pabdpbλ+cδ2pabd,is a unique positive solution of the following quadratic equation:(61)pabdu2+pβ+pbd−pabλ−acδu−pbλ−cδ=0.It is not hard to verify that ifℛ0>1, then λ−du∗>0. So, E∗=u∗,w∗,v∗ if ℛ0>1 is a unique interior equilibrium of system (59).For simplicity of the numerical illustration, we consider only one spatial dimension with spatial domainx∈Ω=0,6. The numerical simulations are observed for a time duration of 500 or 1000 days. Here, we have to point out that some parameter values are taken from the literature [24]. Then, we set choose the parameter values as m1=m2=m3=2, λ=10, d=0.09, β=0.0005, δ=0.4, c=2.4, p=300, a=0.5, and b=0.4. Through direct computation, we obtain ℛ0=0.3070<1 and E0=111.1111,0,0. From Theorem 4, we can know that infection-free equilibrium E0 of model (59) is globally asymptotically stable, see Figure 1. Further, we choose the parameter values as m1=m2=m3=2, λ=10, d=0.09, β=0.0025, δ=0.4, c=2.4, p=300, a=0.0015, and b=0.4. Through direct computation, we obtain ℛ0=74.4048>1 and E∗=104.9054,1.3963,174.5354. From Theorem 5, we can know that the infection equilibrium E∗ of model (59) is globally asymptotically stable, see Figure 2.Figure 1
Choosingm1=m2=m3=2, λ=10, d=0.09, β=0.0005, δ=0.4, c=2.4, p=300, a=0.5, and b=0.4, we have ℛ0=0.3070<1 and the infection-free equilibrium E0=111.1111,0,0 of model (59) is globally asymptotically stable.Figure 2
Takingm1=m2=m3=2, λ=10, d=0.09, β=0.0025, δ=0.4, c=2.4, p=300, a=0.0015, and b=0.4, we have ℛ0=74.4048>1 and the infection equilibrium E∗=104.9054,1.3963,174.5354 of model (59) is globally asymptotically stable.Secondly, we perform the numerical simulation for competitive exclusion. Taking the incidence function asgiux,t,vx,t=βiux,tvix,t/1+aiux,t+bivix,t+aibiux,tvix,t and n=3, we can change model (7) into the following form:(62)∂ux,t∂t=∇⋅Duu∇u+λ−dux,t−∑i=13βiux,tvix,t1+aiux,t+bivix,t+aibiux,tvix,t,∂wix,t∂t=∇⋅Dwiwi∇wi+βiux,tvix,t1+aiux,t+bivix,t+aibiux,tvix,t−δiwix,t,∂vix,t∂t=∇⋅Dvivi∇vi+piwix,t−civix,t,i=1,2,3,with Neumann boundary conditions (40) and initial conditions (41). By ℛi=pi/ciδi∂giλ/d,0/∂vi, together with direct computation, we can get the basic reproduction number for viral strain i for model (62) as ℛi=piλβi/ciδid+aiλ. Meanwhile, we can obtain the equilibria E03=λ/d,0,0, where 0 is a 3-dimensional vector of zeros, and Ei=u∗,0,wi∗,0,0,vi∗,0 for each strain ii=1,2,3, where wi∗=λ−du∗/δi, vi∗=piλ−du∗/ciδi, and(63)u∗=λpiaibi+aiciδi−piβi−dpibi+λpiaibi+aiciδi−piβi−dpibi2+4dλpiaibipibi+ciδi2dpiaibi.Choosing the parameter values asm1=m2=m3=2, λ=10, and d=0.09 and β1=0.0025, δ1=0.4, c1=2.4, p1=300, a1=0.0015, and b1=0.4 for strain 1; β2=0.0015, δ2=0.6, c2=2.4, p2=250, a2=0.5, and b2=0.2 for strain 2; β3=0.0002, δ3=0.2, c3=3.2, p3=320, a3=0.8, and b3=0.6 for strain 3, we can obtain ℛ1=74.4048>1, ℛ2=0.5147<1, ℛ3=0.1357<1, and E1=104.9054,1.3963,0,0,174.5354,0,0 by direct computation. Moreover, we can get that(64)α1αk=1+aiux,t+bivix,t+aibiux,tvix,t1+a1ux,t+b1v1∗+a1b1ux,tv1∗≥11+a1λ/d+b1v1∗+a1b1λ/dv1∗=0.0121,k=2,3,which yields that α1/α2×ℛ1=0.9003>0.5147=ℛ2 and α1/α3×ℛ1=0.9003>0.1357=ℛ3. From Theorem 6, we can know that the infection equilibrium E1 of model (62) is globally asymptotically stable, which implies that strain one persists, while strains 2 and 3 die out. These facts are numerically confirmed in Figure 3.Figure 3
Takingm1=m2=m3=2, λ=10, and d=0.09 and β1=0.0025, δ1=0.4, c1=2.4, p1=300, a1=0.0015, and b1=0.4 for strain 1; β2=0.0015, δ2=0.6, c2=2.4, p2=250, a2=0.5, and b2=0.2 for strain 2; β3=0.0002, δ3=0.2, c3=3.2, p3=320, a3=0.8, and b3=0.6 for strain 3, we have ℛ0=74.4048>1, ℛ2=0.5147<1, ℛ3=0.1357<1, and the infection equilibrium E1=104.9054,1.3963,0,0,174.5354,0,0 of model (62) is globally asymptotically stable.Further, we take another parameter value asβ1=0.0025, δ1=0.4, c1=2.4, p1=300, a1=0.0015, and b1=0.1 for strain 1; β2=0.0015, δ2=0.6, c2=2.4, p2=250, a2=0.25, and b2=0.2 for strain 2; β3=0.0002, δ3=0.2, c3=3.2, p3=320, a3=0.092, and b3=0.6 for strain 3. By direct computation again, we can also obtain ℛ1=74.4048>1, ℛ2=1.0055>1, ℛ3=1.0204>1, E1=89.5396,4.8536,0,0,606.6994,0,0, and(65)α1αk=1+aiux,t+bivix,t+aibiux,tvix,t1+a1ux,t+b1v1∗+a1b1ux,tv1∗≥11+a1λ/d+b1v1∗+a1b1λ/dv1∗=0.0139,k=2,3,which yields that α1/α2×ℛ1=1.0342>1.0055=ℛ2 and α1/α3×ℛ1=1.0342>1.0204=ℛ3. From Theorem 6, we can know that the infection equilibrium E1 of model (62) is globally asymptotically stable, which implies that strain one persists, while strains 2 and 3 die out. These facts are numerically confirmed in Figure 4.Figure 4
Takingm1=m2=m3=2, λ=10, and d=0.09 and β1=0.0025, δ1=0.4, c1=2.4, p1=300, a1=0.0015, and b1=0.1 for strain 1; β2=0.0015, δ2=0.6, c2=2.4, p2=250, a2=0.25, and b2=0.2 for strain 2; β3=0.0002, δ3=0.2, c3=3.2, p3=320, a3=0.092, and b3=0.6 for strain 3, we have ℛ0=74.4048>1, ℛ2=1.0055>1, ℛ3=1.0204>1, and the infection equilibrium E1=89.5396,4.8536,0,0,606.6994,0,0 of model (62) is globally asymptotically stable.Next, we take the parameter values asm1=m2=m3=2, λ=10, and d=0.09 and β1=0.0005, δ1=0.4, c1=2.4, p1=300, a1=0.5, and b1=0.4 for strain 1; β2=0.0015, δ2=0.6, c2=2.4, p2=250, a2=0.5, and b2=0.2 for strain 2; β3=0.0002, δ3=0.2, c3=3.2, p3=320, a3=0.8, and b3=0.6 for strain 3. By direct computation, we can get ℛ1=0.3070<1, ℛ2=0.5147<1, ℛ3=0.1357<1, and E03=111.1111,0,0,0,0,0,0. From Theorem 7, we can know that infection-free equilibrium E03 of model (62) is globally asymptotically stable, which implies that all viral strains die out. These facts are numerically confirmed in Figure 5.Figure 5
Takingm1=m2=m3=2, λ=10, and d=0.09 and β1=0.0005, δ1=0.4, c1=2.4, p1=300, a1=0.5, and b1=0.4 for strain 1; β2=0.0015, δ2=0.6, c2=2.4, p2=250, a2=0.5, and b2=0.2 for strain 2; β3=0.0002, δ3=0.2, c3=3.2, p3=320, a3=0.8, and b3=0.6 for strain 3, we have ℛ1=0.3070<1, ℛ2=0.5147<1, ℛ3=0.1357<1, and the infection-free equilibrium E03=111.1111,0,0,0,0,0,0 of model (62) is globally asymptotically stable.


## 5. Conclusion

In this paper, we investigate the global threshold dynamics in single-strain viral infection models (1) and competitive exclusion in multistrain viral infection models (7) under the homogeneous Neumann boundary conditions.For model (1), by constructing a Lyapunov functional, we obtain that the infection-free equilibrium E0 of model (1) is globally asymptotically stable if the basic reproduction number ℛ0<1, and the infection equilibrium E∗ of model (1) is globally asymptotically stable if the basic reproduction number ℛ0>1. These facts are numerically confirmed in Figures 1 and 2, respectively. For model (7), by constructing a Lyapunov functional, we also obtain that the steady state E1 of model (7) is globally asymptotically stable if ℛ1>1 and ℛ1>ℛk for k=2,…,n, which implies that strain one persists, while strains k=2,…,n die out. These facts are numerically confirmed in Figures 3 and 4. Meanwhile, the steady state E0n of model (7) is globally asymptotically stable if ℛi<1 for i=1,…,n, which implies that all viral strains die out. These facts are numerically confirmed in Figure 5. Here, we have to point out that the incidence function in the present paper is general and performed the numerical simulation for global threshold dynamics and competitive exclusion. So, our results in the present paper are new and can also be applied to other reaction-diffusion systems. Of course, we hope that our work could be instructive to study the dynamics of the viral infection model.



---

*Source: 4923856-2020-04-09.xml*