




































©2023 Ada Academica https://adac.eeEur. J. Math. Anal. 3 (2023) 1doi: 10.28924/ada/ma.3.1
Global Analysis of a Spatiotemporal Cellular Model for the Transmission of Hepatitis C Virus

With Hattaf-Yousfi Functional Response

Alexis Nangue1,∗, Bruno Nde Tchiffo2

1University of Maroua, Higher Teachers’ Training College, Department of Mathematics, P.O.Box 55
Maroua, Cameroon

alexnanga02@gmail.com
2University of Maroua, Faculty of Science, Department of Mathematics and Computer Science, P.O. Box

814 Maroua, Cameroon
ndebruno2@gmail.com

∗Correspondence: alexnanga02@gmail.com

Abstract. This paper carries out a mathematical analysis of the global dynamics of a partial differen-tial equation viral infection cellular model. We study the dynamics of a hepatitis C virus (HCV) model,under therapy, that considers both absorption phenomenon and diffusion of virions, infected and un-infected hepatocytes in the liver. Firstly, we prove the boundedness of the potential solutions, globalexistence, uniqueness, and positivity of the obtained initial value and boundary problem solution.Then, the dynamical behaviour of the model is entirely determined by a threshold parameter calledthe basic reproduction number denoted R0. We show that the uninfected spatially homogeneousequilibrium of the model is globally asymptotically stable if R0 ≤ 1 by using the direct Lyapunovmethod. The latter means that the HCV infection is cleared, and the disease dies out. Also, the globalasymptotical properties stability of the infected spatially homogeneous equilibrium of the model arestudied via a skilful construction of a suitable Lyapunov functional. It means that the HCV infectionpersists in the host, and the infection becomes chronic. Finally, numerical simulations are performedto support the obtained theoretical results.

1. Introduction
The dynamics of viruses, in particular the dynamics of the hepatitis C virus, remains a veryactive field of research in the world of sciences. Moreover, the 2020 Nobel Prize in Medicinewas awarded to three researchers, namely the British Michael Hougton and the Americans HarveyAlter and Charles Rice. They were awarded this Nobel Prize for their very advanced research workon the hepatitis C virus. According to World Health Organization(WHO) [41], 71 million personswere living with chronic hepatitis C virus (HCV) infection worldwide and 399 000 persons had
Received: 29 Apr 2022.
Key words and phrases. reaction-diffusion model; HCV infection; Hattaf-Yousfi functional response; semigroup; globalstability; variational method; cure rate. 1

https://adac.ee
https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 2died from cirrhosis or hepatocellular carcinoma following a survey done in 2015. Aside from theburden of HCV infection secondary to liver-related sequelae, HCV causes an additional burdenthrough comorbidities among persons with HCV infection, including depression, diabetes mellitusand chronic renal disease. In May 2016, the World Health Assembly endorsed the Global HealthSector Strategy for 2016-2021 on viral hepatitis (HBV and HCV infection), which proposes toeliminate viral hepatitis as a public health threat by 2030. Elimination is defined as a 90%reduction in new chronic infections and a 65% reduction in mortality compared with the 2015baseline. Mathematicians cannot stay aside from this disastrous situation decried by WHO. Inview of the vital importance of the liver and the aforementioned facts, any contribution to a betterunderstanding of HCV infection process and strategy to eradicate this infection is of great interest.Mathematical models have been developed to help understand and control the dynamics of HCVwithin an infected host such as in [6, 7, 14, 35]. The dynamics of viral infections such as theEbola virus disease(EVD), the human immunodeficiency virus (HIV) infection, the hepatitis B virus(HBV) infection, the hepatitis C virus (HCV) infection and, new corona virus infection have beenmodeled mathematically in a host. One of the earliest temporal models was the within-host basicviral infection model proposed in [31] to study HIV infection, and later adopted to HBV [8, 32].Particularly, numerous mathematical models describing the temporal dynamics of HCV have beeninitially proposed by Neumann and al [30] using the classical viral infection cellular model, andlater have been extended in [6, 10, 14, 35]. Motivated by what has been done in [8, 30, 32], Chongand al. [7] formulated the basic HCV temporal intra-host model with therapy as a system of threedifferential equations : 

dH(t)

dt
= λ− dH(t)− (1− η)βH(t)V (t),

dI(t)

dt
= (1− η)βH(t)V (t)− αI(t), (1.1)

dV (t)

dt
= (1− ε)kI(t)− µV (t),

where the equations relate the dynamics relationship between, H as the uninfected target cells(hepatocytes), I as the infected cells and V as the viral load (amount of viruses present in theliver). In the system (1.1) the key assumption is that hepatocytes and viruses are well mixed, andneglects the mobility of hepatocytes C viruses, the infected and uninfected target cells. To studythe influences of spatial structures of virus dynamics, Wang and Wang in [39] assuming that themotion of virus follows Fickian diffusion, that is to say, the population flux of virus is proportionalto the concentration gradient and the proportionality constant is taken to be negative [13]. More-over, in model (1.1), the rate of infection is assumed to be bilinear in the virus V and uninfectedhepatocytes T. It is shown in [29] that this bilinear rate of infection could be unrealistic. However,the actual incidence rate is probably not linear over the entire range of T and V. Thus is reason-able to assume that the infection rate is given by a more general one, known as the Hattaf-Yousfi

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 3

functional response [18] of the form βHV
α0+α1H+α2V +α3HV

where α0 > 0, α1 ≥ 0, α2 ≥ 0, α3 ≥ 0are constants. The function βH
α0+α1H+α2V +α3HV

satisfies the hypotheses (H1), (H2) and (H3) ofgeneral incidence rate presented in [16,19–21]. The Hattaf-Yousfi type of functional response wasintroduced by Hattaf and al. [18]. This functional response generalizes many functional responsesand it was used in [34] to describe the dynamics of labour market. Thus, when α0 = 1, the Hattaf-Yousfi functional response is reduced to the specific functional response used by Hattaf and alin [17]. Furthermore, if α3 = α1α2 and α0 = 1, the Hattaf-Yousfi functional response is reducedto Crowley-Martin functional response [9] and was used in [43]. When α3 = 0 et α0 = 1 theHattaf-Yousfi functional response is simplified to Beddington-DeAngelis functional response [5,11],and was used in [25, 26, 38, 42]. When α1 > 0, α2 = α3 = 0 and α0 = 1, the Hattaf-Yousfifunctional response is reduced to Holling type II functional response [28]. And when α1 = α3 = 0,
α2 > 0 and α0 = 1 it expresses a saturation response [36]. Moreover, when α1 = α2 = α3 = 0,and α0 = 1 the Hattaf-Yousfi functional response is reduced to the mass action principle(or Hollingtype I functional response). Also ordinary differential system (1.1) don’t take into consideration thecure of infected hepatocytes. In this work, motivated by the breaches observed in the analysis andthe formulation of system (1.1), we construct and analyze a partial differential equation (PDE)-cellular model system for HCV infection, which derives from system (1.1) by incorporating the space,Hattaf-Yousfi incidence rate, absorption effect and spontaneous cure. It is worth mentioning thatin [7] the authors used mass-action kinetics for viral infection, neglected the cure rate, ignored theabsorption effect and the diffusion of free virions, susceptible cells and infected cells. Thus theobtained model is an extension of the one in the first part of the work done by Chong et al. [7].The work is organized as follows. In section 2, we model the phenomenon described througha reaction-diffusion equations which leads to a initial value and boundary problem. Section 3is devoted to the study of the existence and uniqueness of the global solution of our initial andboundary value problem, and of the properties of this solution, namely positivity and bounded-ness. Section 4 deals with the stability and the analysis of spatially homogeneous equilibria andnumerical simulations in section 5. We conclude our work and provide a discussion in section 6.

2. Formulation of the PDE-cellular model
Let Ω ⊂ R3 be a bounded connected domain representing the liver. Let t ≥ 0 be a given timeand x = (x1, x2, x3) ∈ Ω. Denote respectively by H(x, t), I(x, t) and V (x, t) the concentrations ofhealthy hepatocytes, HCV infected hepatocytes, and free HCV virions at time t and location x . Thedynamics of HCV infection intra-host is the result of the dynamics of each compartment H, I, andV, and the various interactions between them. We now describe the evolution of each compartment.

2.1. Fluctuation of healthy hepatocytes. Let ν be an elementary volume in Ω. The variation ofthe quantity of healthy hepatocytes in ν is described under the following assumptions. Healthy

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 4hepatocytes are produced at constant rate λ from the bone narrow and die at rate dH. Virions infectthe healthy hepatocytes at the rate βHV
α0+α1H+α2V +α3HV

, where β is the rate of transmission of theinfection and αj , j = 0, 1, 2, 3 are positive constants. This generalized incidence function replacesthe mass-action function which has been shown to cause unrealistic conditions for successful chronicHCV infection. ρI is the cure rate of infected hepatocytes either by noncytolytic mechanism orimmunity or treatment. In addition, the therapeutic effect of treatment in this model involved thereduction of new infections, which is described in a fraction as (1 − η). The spatial motion ofhealthy hepatocytes follows the Fickian diffusion law. Thus, the variation of healthy hepatocytesis expressed by the following equation:
∂H

∂t
= D1∆H(x, t) + λ− dH(x, t)−

(1− η)βH(x, t)V (x, t)

α0 + α1H(x, t) + α2V (x, t) + α3H(x, t)V (x, t)
+ ρI(x, t),

where D1 represents the Healthy hepatocytes diffusion coefficient and
∆ =

∂2

∂x2
1

+
∂2

∂x2
2

+
∂2

∂x2
3is the usual Laplacian operator in three-dimensional space.

2.2. Fluctuation of HCV infected cells. The HCV infected cells die at rate α per day so that 1
αis the life-expectancy of HCV infected hepatocytes. Healthy hepatocytes become infected at therate βHV

α0+α1H+α2V +α3HV
. The spatial motion of HCV infected cells follows the Fickian diffusion law.Thus, the variation of infected hepatocytes is expressed by the following equation

∂I

∂t
= D2∆I(x, t) +

(1− η)βH(x, t)V (x, t)

α0 + α1H(x, t) + α2V (x, t) + α3H(x, t)V (x, t)
− (α+ ρ)I(x, t),

where D2 represents the HCV infected cells diffusion coefficient.
2.3. Fluctuation of free HCV virions. The infected hepatocytes produce virus at rate kI , and virusis cleared at the rate µV . Also, the population of virions decreases due to the infection at the rate

u(1−η)βHV
α0+α1H+α2V +α3HV

due to absorption effect, where u ∈ {0, 1}. The spatial motion of virions followsthe Fickian diffusion law. In addition, the therapeutic effect of treatment in this model involvedblocking virions production (referred to as drug effectiveness) which, is described in fraction (1−ε).Thus, the variation of free virions is expressed by the following equation:
∂V

∂t
= D3∆V (x, t)+(1−ε)kI(x, t)−µV (x, t)−u

(1− η)βH(x, t)V (x, t)

α0 + α1H(x, t) + α2V (x, t) + α3H(x, t)V (x, t)
,

where D3 represents the free HCV virions diffusion coefficient.
2.4. The initial boundary value problem associated to PDE-cellular model. In this section, we usethe previous equations describing variables variations to set up a complete PDE system modellingbiological dynamics for HCV infection. Let T > 0 be a fixed time and define

ΩT = Ω× (0, T ).

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 5Therefore, in ΩT the full system of PDE governing the HCV infection becomes :

∂H

∂t
= D1∆H(x, t) + λ− dH −

(1− η)βHV

α0 + α1H + α2V + α3HV
+ ρI,

∂I

∂t
= D2∆I(x, t) +

(1− η)βHV

α0 + α1H + α2V + α3HV
− (α+ ρ)I, (2.1)

∂V

∂t
= D3∆V + (1− ε)kI − µV −

u(1− η)βHV

α0 + α1H + α2V + α3HV
.

We use the Neumann homogeneous boundary conditions:
∂H

∂η
=
∂I

∂η
=
∂V

∂η
= 0 on ∂Ω× [0, T ], (2.2)

where ∂
∂η denotes the outward normal derivative on ∂Ω. The initial conditions are the following :

H(x, 0) = H0, I(x, 0) = I0, V (x, 0) = V0, x ∈ Ω. (2.3)
The boundary conditions in (2.2) imply that the Healthy hepatocytes, the HCV infected cells and freeHCV virions do not move across the boundary ∂Ω. For an epidemiological significance, we assumethat the initial conditions are positive and Hölder continuous, and satisfy ∂H0

∂η = ∂I0
∂η = ∂V0

∂η = 0 on
∂Ω. We then obtain the following initial boundary value problem, denoted IBVP associated to theprevious PDE-cellular model:

∂H

∂t
= D1∆H + λ− dH −

(1− η)βHV

α0 + α1H + α2V + α3HV
+ ρI in ΩT ,

∂I

∂t
= D2∆I +

(1− η)βHV

α0 + α1H + α2V + α3HV
− (α+ ρ)I in ΩT ,

∂V

∂t
= D3∆V + (1− ε)kI − µV −

u(1− η)βHV

α0 + α1H + α2V + α3HV
in ΩT , (2.4)

∂H

∂η
=
∂I

∂η
=
∂V

∂η
= 0 on ∂Ω× [0, T ],

H(x, 0) = H0, I(x, 0) = I0, V (x, 0) = V0, x ∈ Ω,on which our study will focus on.
3. Qualitative and quantitative analysis and some properties of the solutions for IBVP (2.4)

In this section, we provide a thorough study of the dynamics of IBVP (2.4) which yields variousoutcomes. Precisely, we prove existence, uniqueness, positivity and boundedness of solutions forIBVP (2.4). This is done by combining variational method and semigroups techniques to someuseful functional analysis arguments.
3.1. Local existence and uniqueness of solutions for the IBVP (2.4). Set

F (H, I, V ) = (F1(H, I, V ), F2(H, I, V ), F3(H, I, V ))T (3.1)
where

F1(H, I, V ) = λ− dH −
(1− η)βHV

α0 + α1H + α2V + α3HV
+ ρI,

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 6

F2(H, I, V ) =
(1− η)βHV

α0 + α1H + α2V + α3HV
− (α+ ρ)I,and

F3(H, I, V ) = (1− ε)kI − µV − u
(1− η)βHV

α0 + α1H + α2V + α3HV
.We have the following result which guarantees that the right-hand side, without diffusion, of thePDE-model system (2.4) is Lipschitz.

Proposition 3.1. Let T ∈ R∗+ and (H, I, V ) ∈ (C0
B(Ω × [0, T ]))3, where C0

B(Ω × [0, T ]) is the
space of bounded and continuous functions on Ω × [0, T ]. We suppose that F in (3.1) is defined
on L2(Ω× (0, T )). Then F1, F2 and F3 are uniformly Lipschitz continuous on L2(Ω× (0, T )) with
respect to H, I and V.

Proof. Let T ∈ R∗+ and (H1, I1, V1), (H2, I2, V2) ∈ (C0
B(Ω× [0, T ]))3. First, by direct computation,we have :

‖F1(H1, I1, V1)− F1(H2, I2, V2)‖2 ≤ K1
1‖H1 −H2‖2 +K1

2‖I1 − I2‖2 +K1
3‖V1 − V2‖2, (3.2)

with
K1

1 = d + (1− η)β

(
1

α2
+
Vm
α0

)
, K1

2 = ρ, K1
3 = (1− η)β

(
1

α1
+
Hm
α0

)
, (3.3)with Hm and Vm given below.Then

‖F2(H1, I1, V1)− F2(H2, I2, V2)‖2 ≤ K2
1‖H1 −H2‖2 +K2

2‖I1 − I2‖2 +K2
3‖V1 − V2‖2, (3.4)

with
K2

1 = (1− η)β

(
1

α2
+
Vm
α0

)
, K2

2 = (α+ ρ) and K2
3 = (1− η)β

(
1

α1
+
Hm
α0

)
. (3.5)

Finally
‖F3(H1, I1, V1)− F3(H2, I2, V2)‖2 ≤ K3

1‖H1 −H2‖2 +K3
2‖I1 − I2‖2 +K3

3‖V1 − V2‖2, (3.6)
with
K3

1 = u(1− η)β

(
1

α2
+
Vm
α0

)
, K3

2 = k(1− ε), and K3
3 = µ+ u(1− η)β

(
1

α1
+
Hm
α0

)
. (3.7)

This completes the proof of Proposition 3.1. �

Now, consider the following IBVP

∂tH −D1∆H = f (t, H, I, V ) in Ω× (0, T )

∂t I −D2∆I = g(t, H, I, V ) in Ω× (0, T )

∂tV −D3∆V = h(t, H, I, V ) in Ω× (0, T )

∂H
∂η = 0; ∂I

∂η = 0; ∂V
∂η = 0 on ∂Ω× [0, T ]

H = H0, I = I0, V = V0 on Ω× {t = 0}.

(3.8)

In what follows, we will need the following definition and results.

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 7

Definition 3.1. (Sectorial operator, [22]) Let A be a linear operator in a Banach space X andsuppose A is closed and densely defined. If there exist real numbers a, ω ∈ (0, π), M ≥ 1 suchthat
ρ(A) ⊃ Σ = {λ0 ∈ C : ω ≤ arg(λ0 − a) ≤ π, λ0 6= 0} (3.9)and

‖Rλ0
(A)‖ ≤

M

|λ0 − a|
for all λ0 ∈ Σ, (3.10)then we say that A is sectorial.

Remark 3.1. The Neumann realization of the Laplacian A = −∆, with domain
D(A) =

{
ω ∈ H2(Ω) :

∂ω

∂η
= 0
}

is a sectorial operator in L2(Ω). But since C∞0 (Ω) ⊂ D(A), it is densely defined in L2(Ω). For
β ≥ 0 large enough, we define the fractional powers of the Helmholtz operator, H� = −∆ + βI ,with domain D(H�) equipped with graph norm ‖ . ‖D(H�) = ‖ . ‖2 + ‖H�.‖2.

We have the following general results.
Lemma 3.2. [1] Let 1 ≤ p <∞. Then D(H�) ⊂ C∞0 (Ω) with continuous injection for β > n

2p .

Lemma 3.3. [22] D(Hβ) ⊂ C0
B(Ω) with continuous injection for β > n

4 .

Theorem 3.4. [22] If A is sectorial, then −A is the infinitesimal generator of an analytic semigroup,
G(t).
If Rλ0

> a, a ∈ R whenever λ0 ∈ σ, then for any t > 0,

‖G(t)‖ ≤ Ce−at , ‖AG(t)‖ ≤
C

t
e−at

and
d

dt
G(t) = −AG(t) , t > 0.

Corollary 3.5. Let G be the analytic semigroup generated by −A. The following properties hold
for the semigroup G and the fractional powers of the Helmholtz operator Hβ :1) G(t) : L2(Ω)→ D(H�) for all t > 0,2) ‖G(t)ω‖H� ≤ Cβ,2t−β‖ω‖2 for all t > 0, ω ∈ L2,3) G(t)H�ω = H�G(t)ω for all t > 0, ω ∈ D(H�).

Remark 3.2. The following basic hypotheses are assumed to hold :
(H1): D1 > 0, D2 > 0 and D3 > 0,
(H2): H0 ≥ 0, I0 ≥ 0 and V0 ≥ 0 are continuous on Ω, H0, I0, V0 ∈ C0

B(Ω),
(H3): f , g and h are continuously differentiable functions fromR4

+ intoR with f (t, 0, s, z) ≥ 0,
g(t, r, 0, z) ≥ 0 and h(t, r, s, 0) ≥ 0 for all t , r , s , z ≥ 0.

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 8For x ∈ Ω, t ≥ 0, H, I , V ∈ (C0
B(Ω))3, define F ,G and Q on R+ × (C0

B(Ω))3 by :
[F(t, H, I, V )](x) = f (t, H(x), I(x), V (x)), [G(t, H, I, V )](x) = g(t, H(x), I(x), V (x)),

[Q(t, H, I, V )](x) = h(t, H(x), I(x), V (x)).In addition, we let G1, G2 and G3 be the analytical semigroup generated by A1 = D1 × ∆,
A2 = D2 × ∆ and A3 = D3 × ∆ respectively.

In the sequel, we will need the following results.
Lemma 3.6. [22] If H, I and V are continuous from [0, T ] to L2(Ω), then the integrals :

I1(t) =
∫ t

0 G1(t − τ)F(τ,H(τ), I(τ), V (τ))dτ , I2(t) =
∫ t

0 G2(t − τ)G(τ,H(τ), I(τ), V (τ))dτ ,
I3(t) =

∫ t
0 G3(t − τ)Q(τ,H(τ), I(τ), V (τ))dτ ,

exist and I1(t), I2(t) and I3(t) are continues on [0, T [ with I1(t) ∈ D(A1), I2(t) ∈ D(A2),
I3(t) ∈ D(A3) and I1(t)→ 0+ in L2 as t → 0+, I2(t)→ 0+ in L2 as t → 0+ and I3(t)→ 0+ in
L2 as t → 0+.

Lemma 3.7. If the IBVP (3.8) has a classical solution, then H, I and V satisfy the following equalities
:

H(t) = G1(t)H0 +

∫ t

0

G1(t − τ)F(τ,H(τ), I(τ), V (τ))dτ, (3.11)
I(t) = G2(t)I0 +

∫ t

0

G2(t − τ)G(τ,H(τ), I(τ), V (τ))dτ, (3.12)
V (t) = G3(t)V0 +

∫ t

0

G3(t − τ)Q(τ,H(τ), I(τ), V (τ))dτ. (3.13)
Proof. Consider the L2−valued functions θj(τ) = Gj(t − τ)ωj(τ), j = 1, 2, 3 with ω1 = H, ω2 = Iand ω3 = V . Then θj is differentiable since Gj is analytic and ωj is differentiable. Then byTheorem 3.4, we have
dθ1

dτ
=

d

dτ

[
G1(t − τ)

]
H(τ) + G1(t − τ)H

′
(τ),

= −D1 × ∆
(
G1(t − τ)

)
H(τ) + G1(t − τ)

[
D1 × ∆H(τ) + F(τ,H(τ), I(τ), V (τ))

]
,

= −D1 × ∆
(
G1(t − τ)

)
H(τ) +D1 × G1(t − τ)× ∆H(τ) + G1(t − τ)F(τ,H(τ), I(τ), V (τ)),

According to corollary 3.5 with β = 0, we have
D1 × G1(t − τ)∆H(τ) = D1 × ∆G1(t − τ)H(τ).

Therefore
dθ1

dτ
= −D1 × ∆G1(t − τ)H(τ) +D1 × ∆G1(t − τ)H(τ) + G1(t − τ)F(τ,H(τ), I(τ), V (τ)),

= G1(t − τ)F(τ,H(τ), I(τ), V (τ)). (3.14)

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 9In the similar way, we also have :
dθ2

dτ
=

d

dτ

[
G2(t − τ)

]
I(τ) + G2(t − τ)I

′
(τ),

= −D2 × ∆
(
G2(t − τ)

)
I(τ) +D2 × G2(t − τ)× ∆I(τ) + G2(t − τ)G(τ,H(τ), I(τ), V (τ)),

= G2(t − τ)G(τ,H(τ), I(τ), V (τ)). (3.15)
dθ3

dτ
=

d

dτ

[
G3(t − τ)

]
V (τ) + G3(t − τ)V

′
(τ),

= −D × ∆G3(t − τ)V (τ) +D × G3(t − τ)∆V (τ) + G3(t − τ)Q(τ,H(τ), I(τ), V (τ)),

= G3(t − τ)Q(τ,H(τ), I(τ), V (τ)). (3.16)
Integrating equations (3.14), (3.15) and (3.16) with respect to time, we obtain equations (3.11),(3.12) and (3.13) respectively. �

Remark 3.3. Since in this work, n = 3, we take p = 2 so that β > 3
4 and therefore the domain

D(H�) is continuously embedded in C∞0 (Ω) by Lemma 3.2. Now, let H, I and V be continuousfunctions from [0, T ] to D(H�) ↪→ C0
B(Ω) satisfying (3.11), (3.12) and (3.13) respectively. We canthen claim that H, I and V verify system (3.8). The continuity of H, I and V implies continuity of

t 7→ F(t, H(t), I(t), V (t)), t 7→ G(t, H(t), I(t), V (t)) and t 7→ Q(t, H(t), I(t), V (t)).One can then conclude that, the linear Cauchy problem
∂ty1 −D1∆y1 = F(t, H(t), I(t), V (t)),

∂ty2 −D2∆y2 = G(t, H(t), I(t), V (t)),

∂ty3 −D∆y3 = Q(t, H(t), I(t), V (t)),

y1(0) = H0, y2(0) = I0, y3(0) = V0,

has a unique solution, with y1, y2 and y3 given by (3.11), (3.12) and (3.13) respectively.
Following [2], [4], [3], [22], [24], we have the following main result for the local existence of (2.4),based on L2-theory.

Proposition 3.8. If hypotheses (H1), (H2) and (H3) are satisfied, then the initial value and boundary
problem (3.8) admits a unique solution (H, I, V ) ∈ (C0

B(]0, T ], D(Hβ)))3, with H(0) = H0 ∈ C0
B(Ω),

I(0) = I0 ∈ C0
B(Ω) and V (0) = V0 ∈ C0

B(Ω).

The proof of this proposition is given in "Appendix A ".
3.2. Boundedness of the solutions for IBVP (2.4).

Proposition 3.9. Let (H, I, V ) ∈
(
C0
(

Ω× [0, T )
)
∩ C2,1

B (Ω× [0, T ))
)3

be the solution of (2.4) with bounded initial conditions i.e. 0 < H0(x) < Hm, 0 < I0(x) < Hm,

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 10

0 < V0(x) < Vm for all x ∈ Ω, and satisfying the boundary condition ∂H0
∂η = 0, ∂I0∂η = 0, ∂V0

∂η = 0

on ∂Ω. Then,

∀(x, t) ∈ Ω× [0, T ], H(x, t) ≤ Hm, I(x, t) ≤ Hm and V (x, t) ≤ Vm

with

Hm = max

{
λ

δ2
,max
x∈Ω̄
{H(x, 0) + I(x, 0)}

}
and Vm = max

{
(1− ε)kHm

µ
,max
x∈Ω

V0(x)

}
.

Proof. Consider the function S defined for all (x, t) ∈ Ω× [0, T ] by
S(x, t) = H(x, t) + I(x, t).

Adding the first two equations in (2.4), yields
∂S(x, t)

∂t
−D1∆H(x, t)−D2∆I(x, t) = λ− dH(x, t)− αI(x, t).

It follows that
∂S(x, t)

∂t
−max{D1, D2}∆ (H(x, t) + I(x, t)) ≤ λ−min{d, α} (H(x, t) + I(x, t)) ,

we have 

∂S(x,t)
∂t − δ1∆S(x, t) ≤ λ− δ2S(x, t), x ∈ Ω, t ∈ [0, T ]

∂S(x,t)
∂η = 0, x ∈ ∂Ω, t ∈ [0, T ]

S(x, 0) = max
x∈Ω

S0(x),

(3.17)
where S0(x) = {H(x, 0) + I(x, 0)}, δ1 = max{D1, D2} et δ2 = min{d, α}. By using the standardparabolic comparison of the scalar parabolic equations [33], one has

S(x, t) ≤ S̄(t),

where S̄(t) = λ
δ2

(
1− e−δ2t

)
+ max

x∈Ω̄
S0(x)e−δ2t is the solution of the problem

dS̄(t)
dt = λ− δ2S̄(t),

S̄(0) = max
x∈Ω̄

S0(x),

(3.18)
which dominates system (3.17). The general solution of (3.18) is on the form S̄(t) = k(t)e−δ2t . ByLagrange’s method, we have k(t) = λ

δ2
eδ2t + c , c ∈ R. Hence

S̄(t) =

(
λ

δ2
eδ2t + c

)
e−δ2t .

Initial condition yields c = max
x∈Ω

S0(x)− λ
δ2

.Therefore
S̄(t) =

λ

δ2

(
1− e−δ2t

)
+ max

x∈Ω̄
S0(x)e−δ2t .

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 11Then, it follows that :
S(x, t) ≤ S̄(t)

≤
λ

δ2

(
1− e−δ2t

)
+ max

x∈Ω̄
S0(x)e−δ2t

≤ max

{
λ

δ2
,max
x∈Ω̄

S0(x)

}(
1− e−δ2t

)
+ max

{
λ

δ2
,max
x∈Ω̄

S0(x)

}
e−δ2t

≤ max

{
λ

δ2
,max
x∈Ω̄

S0(x)

}
.

Thus,
S(x, t) ≤ max

{ λ
δ2
,max
x∈Ω̄
{H(x, 0) + I(x, 0)}

}
.

Therefore
S(x, t) ≤ Hm = max

{
λ

δ2
,max
x∈Ω̄
{H(x, 0) + I(x, 0)}

}
,∀(x, t) ∈ Ω× [0, Tmax),

where Tmax is the maximal time of existence of the solution of system (2.4), this implies that S isbounded.Hence H and I are bounded since S is bounded. This prove that H and I are bounded.Now, to show that V is bounded, from the third equation of IBVP (2.4), we have
∂V (x, t)

∂t
−D3∆V (x, t) ≤ (1− ε)kI(x, t)− µV (x, t), x ∈ Ω, t ∈ [0, T ]

∂V (x, t)

∂η
= 0, x ∈ ∂Ω, t ∈ [0, T ]

V (x, 0) = max
x∈Ω

V0(x).

It follows from the previous system, inequality

∂V (x,t)
∂t −D∆V (x, t) ≤ (1− ε)kHm − µV (x, t)

∂V (x,t)
∂η = 0

V (x, 0) = max
x∈Ω̄

V0(x),

(3.19)
By using the standard parabolic comparison of the scalar parabolic equations [33], one has

V (x, t) ≤ V (t),

where V (t) = (1−ε)kHm
µ (1− e−µt) + maxx∈Ω̄ V0(x)e−µt is the solution of the problem

dV (t)
dt = (1− ε)kHm − µV (t),

V (0) = max
x∈Ω

V0(x),

(3.20)

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 12which dominates system (3.19). Indeed, the general solution of (3.20) is on the form V (t) =

c(t)e−µt . By the Lagrange’s method, we have c(t) = (1−ε)kHm
µ (eµt − 1) + c0, c0 ∈ R.thus,

V (t) =
[(1− ε)kHm

µ
(eµt − 1) + c0

]
e−µt .

Initial condition yields max
x∈Ω̄

V0(x) = V (0) = c0. It follows from that
V (t) =

(1− ε)kHm
µ

(1− e−µt) + max
x∈Ω̄

V0(x)e−µt .

Therefore
V (x, t) ≤ V (t)

≤ max
{(1− ε)kHm

µ
,max
x∈Ω̄

V0(x)
}

(1− e−µt) + max
{(1− ε)kHm

µ
,max
x∈Ω̄

V0(x)
}
e−µt

≤ max
{(1− ε)kHm

µ
,max
x∈Ω̄

V0(x)
}
.

Since
V (x, t) ≤ V (t) ≤ max

{(1− ε)kHm
µ

,max
x∈Ω

V0(x)
}
, ∀(x, t) ∈ Ω× [0, Tmax);

where Tmax is the maximal time of existence of the solution of system (2.4), this implies that V isbounded.Thus H(x, t), I(x, t) and V (x, t) are bounded on Ω × [0, Tmax). Therefore, it follows from thestandard theory of semi-linear parabolic system in [23] that Tmax = +∞. This completes the proofof proposition 3.9. �

3.3. Global existence, uniqueness and positivity for the IBVP (2.4). We recast the IBVP (2.4) asfollows:

∂w

∂t
−D∆w + q(w)w = f (w) in Ω× [0, T ),

∂w1

∂η
= 0,

∂w2

∂η
= 0,

∂w3

∂η
= 0 on ∂Ω× [0, T ), (3.21)

w(x, 0) = w0(x) in Ω,

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 13where w = (w1, w2, w3)T = (H, I, V )T , D = diag(D1, D2, D3), q(w) =

diag (q1(w), q2(w), q3(w)), f (w) =
(
f1(w), f2(w), f3(w)

)T , with
q1(w) = d +

(1− η)βw3

α0 + α1w1 + α2w3 + α3w1w3
, q2(w) = (α+ ρ),

q3(w) = µ+
u(1− η)βw1

α0 + α1w1 + α2w3 + α3w1w3
, f1(w) = λ+ ρw2,

f2(w) =
(1− η)βw1w3

α0 + α1w1 + α2w3 + α3w1w3
, f3(w) = (1− ε)kw2.

Note that D1, D2, D3 > 0. Denote H = L2(Ω) and E = H1(Ω) and define as in [12] the Hilbertspace
W (0, T, E, E′) =

{
u ∈ L2 ((0, T ), E) :

∂u

∂t
∈ L2

(
(0, T ), E′

)}
,

endowed with the norm
‖u‖2

W = ‖u‖2
L2((0,T ),E) +

∥∥∥∂u
∂t

∥∥∥2

L2((0,T ),E′)and the following hypothesis for initial conditions:
w01 ∈ L∞(Ω), w02, w03 ∈ H and w0i ≥ 0 for i ∈ {1, 2, 3}. (3.22)

Here, we apply Theorem 2.7 of [12]. So, one approaches the solution by a sequence of solutions oflinear equations. For n = 0, w0 denotes the solution of
∂w0

∂t −D∆w0 = 0 in Ω× (0, T ),

w0(0) = w0 in Ω,
∂w0

i
∂η = 0. on ∂Ω

(3.23)
This equation admits a strong solution and w0 ≥ 0.By induction, wn = ((wn1 , w

n
1 , w

n
1 )) denotes the solution of

∂wn

∂t −D∆wn + q(wn−1)wn = f (wn−1) in Ω× (0, T ),

wn(0) = w0 in Ω,

∂wn

∂η = 0. on ∂Ω.

(3.24)
Since (3.24) is a linear equation, qi(wn−1) and fi(wn−1) can replace a0 and f (t) of Corollary 2.10in [12]. Suppose that there exists a unique nonnegative solution wn−1

i . Assuming by induction that
w ji ≥ 0 for 0 ≤ j ≤ n − 1 and that by Proposition 3.9 w ji is bounded for 0 ≤ j ≤ n − 1, one has

0 ≤
u(1− η)βwn−1

1

α0 + α1w
n−1
1 + α2w

n−1
3 + α3w

n−1
1 wn−1

3

≤ u(1− η)β

(3.25)

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 14which implies that
µ ≤ q3(wn−1) ≤ µ+ u(1− η)β. (3.26)

Since w ji are bounded , we have
d ≤ q1(wn−1) ≤ d + (1− η)β.

In addition, q2 is a constant.It then follows that q1(wn−1), q2(wn−1), q3(wn−1) ∈ L∞(Ω× (0, T )). We also have f (wn−1) ≥ 0and f (wn−1) ∈ L2((0, T ), E′). Then, by Corollary 2.10 of [12], there exists a unique solution
wn ∈ W (0, T, E, E′) with wn ≥ 0. Since f1(w) = λ + ρw2, f2(w) = (1−η)βw1w3

α0+α1w1+α2w3+α3w1w3
≤

(1− η)βw3 and f3(w) = (1− ε)kw2, then f1(wn−1) = λ+ ρwn−1
2 , f2(wn−1) ≤ (1− η)βwn−1

3 and
f3(wn−1) = (1 − ε)kwn−1

2 remain bounded in L2(]0, T [, E). We deduce that wn2 and wn3 remainbounded in C0([0, T ],H) and L2((0, T ), E).Now, we deduce that the sequence (wni )n≥0 (one can extract a subsequence (wmi )m≥0) convergesweakly to wi in L2((0, T ), E) and weakly star in L∞((0, T ),H) to wi . Applying Proposition 2.11in [12], it holds that for all n,
wni (t) = Gi(t)w0i +

∫ t

0

Gi(t − s)gni (s)ds, (3.27)
where Gi(t) is the semigroup generated by the unbounded operator Ai = −DiAH, and

gni (s) = −qi(wn−1(s))wni (s) + fi(w
n−1(s)). (3.28)

Then, gni ∈ L2((0, T ), E). Since the sequence (wni )n≥0 is bounded in C0([0, T ],H), the se-quence (gni )n≥0 is bounded in C0([0, T ],H). Now, consider the operator Gi from C0((0, T ),H)into C0((0, T ),H) defined by
Gi(f ) =

∫ t

0

Gi(t − s)f (s)ds. (3.29)
Let us prove that Gi is a compact operator. Considering the triple (L2(Ω), H1(Ω), a) with

a(w, v) =

3∑
j=1

∫
Ω

∂w

∂xj

∂v

∂xj
dx, (3.30)

where Ω is regular and bounded. As in [12], the unbounded variational operator AH associated to
a is a positive symmetric operator with compact resolvent Rλ(AH). It admits a sequence (λk)k ofpositive eigenvalues with lim

k→+∞
λk = +∞ and a Hilbert basis (ek)k ofH consisting of eigenvectorsof AH. If (G(t))t>0 is the semigroup generated by −AH, then for all w0 ∈ H,
G(t)w0 =

+∞∑
k=0

e−tλk (w0, ek)ek . (3.31)
This proves that the operator is compact for all t > 0 since

lim
k→+∞

e−tλk = 0.

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 15We have the same formula for Gi(t), and it suffices to replace λk by Diλk . Setting
GN(t)w =

N∑
k=0

e−tλk (w, ek)ek , (3.32)
one sees that GN(t) is an operator with finite rank which converges to G(t). The following Theoremis relevant in the sequel.
Theorem 3.10. [12] Let t 7→ G(t) be an application from [0,+∞) into L(H).
One assumes that there exists a sequence of operators (GN(t))N≥0 on H verifying the following
properties:1) for all N and all t > 0, GN(t) has finite rank independent of t ,2) t 7→ GN(t) is continuous from [0,+∞) into L(H) for all N ,3) for N → +∞, GN(t) converges to G(t) in L1(]0, T [,L(H)) for all T > 0.

Then the operator G is compact from C0([0, T ],H) to C0([0, T ],H) for all T > 0.

We are now in the position to prove the global existence, uniqueness and positivity of the solutionto the IBVP (2.4).
Theorem 3.11. If the initial condition satisfies (3.22), then the IBVP (3.21) admits a unique non-
negative solution w ∈ (W (0, T, E, E′))3.

The proof of Theorem 3.11 is contained in "appendix B ".
Remark 3.4. It is worth noting that positivity of the solution may be proved by applying the maximumprinciple. Moreover, from the above results and the boundedness of the solution, one has observedthat the solution of IBVP (2.4) enters the region:

Σ =
{

(H, I, V ) ∈ Ω3 ×R3
+ : 0 < H(x, t) ≤ Hm, 0 < I(x, t) ≤ Hm, 0 < V (x, t) ≤ Vm

}
,

where
Hm = max

{ λ
δ2
,max
x∈Ω
{H(x, 0) + I(x, 0)}

} et Vm = max
{(1− ε)kHm

µ
,max
x∈Ω

V (x, 0)
}
.

Hence the region Σ, of biological interest, is positively-invariant under the flow induced by IBVP(2.4).
4. Stability analysis of the spatially homogeneous equilibria

4.1. HCV-spatial homogeneous uninfected equilibrium E0. The spatial homogeneous uninfectedequilibrium of the PDE-model system (2.4) arises when there is no virus within a host i.e., V=0.Easy calculations shows that the HCV-spatial homogeneous uninfected equilibrium for PDE-modelsystem (2.4) is given by
E0 = (Λ, 0, 0)

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 16where
Λ =

λ

d
.

4.2. Basic reproduction number R0. In order to define the Basic reproduction number R0 forsystem (2.4), we first observe that system (2.4) has a spatially homogeneous uninfected equilibrium
E0. It should be noted that one of the main tools in epidemic models is the basic reproductionnumber R0 which is an important threshold parameter to discuss the dynamic behaviour of theepidemic model. It quantifies the infection risk. It measures the expected average number ofnew infected hepatocytes generated by a single virion in a completely healthy hepatocyte. Itshould be also noted that, while a huge number of works deals with the threshold dynamics forODE-models, very few studies are devoted to PDE-models. This is eventually due to the factthat the concept of basic reproduction number has just recently been extended to PDE-modelssuch as reaction-diffusion and reaction-convection-diffusion epidemic models with mixed boundaryconditions [37,40]. The definition of R0 in this work follows the approach developed in [40].

In order to find the basic reproduction number R0 for the system (2.4), we obtain the followinglinear system at E0 for the infected classes:


∂I

∂t
= D2∆I − (α+ ρ)I +

(1− η)βΛ

α0 + α1Λ
V in ΩT ,

∂V

∂t
= D3∆V + (1− ε)kI − µV −

u(1− η)βΛ

α0 + α1Λ
V in ΩT , (4.1)

∂I

∂η
=
∂V

∂η
= 0 on ∂Ω× [0, T ].

Substituting I(x, t) = eλtψ2(x) and V (x, t) = eλtψ3(x) in (4.1), we obtain the following cooper-ative eigenvalue problem:


λψ2(x) = D2∆ψ2(x)− (α+ ρ)ψ2(x) +
(1− η)βΛ

α0 + α1Λ
ψ3(x) in Ω,

λψ3(x) = D3∆ψ3(x) + (1− ε)kψ2(x)− µψ3(x)−
u(1− η)βΛ

α0 + α1Λ
ψ3(x) in Ω, (4.2)

∂ψ2(x)

∂η
=
∂ψ3(x)

∂η
= 0 on ∂Ω.

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 17As in [40], let T: C(Ω̄,R2)→ C(Ω̄,R2) be the solution semigroup of the following reaction-diffusionsystem: 

∂I

∂t
= D2∆I − (α+ ρ)I in ΩT ,

∂V

∂t
= D3∆V + (1− ε)kI − µV −

u(1− η)βΛ

α0 + α1Λ
V in ΩT , (4.3)

I(x, 0) = ψ2(x), V (x, 0) = ψ3(x), in ΩT

∂I

∂η
=
∂V

∂η
= 0 on ∂Ω.

Thus, with initial infection Ψ(x) = (ψ2, ψ3) the distribution of those infections members becomes
T (t)Ψ(x) as time evolves. Therefore, the distribution of total new infections is∫ ∞

0

F (x)T (t)Ψ(x)dt,

then, we define
L(φ)(x) :=

∫ ∞
0

F (x)T (t)Ψ(x)dt = F (x)

∫ ∞
0

T (t)Ψ(x)dt.

L is a positive and continuous operator which maps the initial infection distribution to the distri-bution of the total infective members produced during the infection period. Applying the idea ofnext generation operators [40], we define the spectral radius of L as the basic reproduction number
R0 := ρ(L).

The matrices F and V defined as
F (x) =


0 (1−η)βΛ

α0+α1Λ

0 0

 , V (x) =

 α+ ρ 0

−(1− ε)k
[
µ+ u

(1−η)βΛ
α0+α1Λ

]
 .

Then
FV −1 =

α0 + α1Λ

(α+ ρ) [µ(α0 + α1Λ) + u(1− η)βΛ]


(1−η)(1−ε)kβΛ

α0+α1Λ
(α+ρ)(1−η)βΛ

α0+α1Λ

0 0

 .
By [40] (theorem 3.4), one has

R0 =
(1− η)(1− ε)kβΛ

(α+ ρ) [µ(α0 + α1Λ) + u(1− η)βΛ]
. (4.4)

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 184.3. Existence and uniqueness of HCV-spatial homogeneous infected equilibrium E∗. In thissection, we address the existence and uniqueness of infected spatial homogeneous equilibrium(2.4). The latter denoted as E∗ = (H∗, I∗, V ∗) with H∗ 6= 0, I∗ 6= 0 et V ∗ 6= 0 satisfying thefollowing algebraic system :
λ− dH∗ − (1− η)L(H∗, I∗, V ∗)V ∗ + ρI∗ = 0,

(1− η)L(H∗, I∗, V ∗)V ∗ − (α+ ρ)I∗ = 0,

(1− ε)kI∗ − µV ∗ − u(1− η)L(H∗, I∗, V ∗)V ∗ = 0,

(4.5)
where

L(H, I, V ) =
βH

α0 + α1H + α2V + α3HV
.Adding the first and second equation of (4.5), we have

λ− dH∗ − αI∗ = 0

which yields
I∗ =

λ− dH∗

α
. (4.6)As far as, using the second and third equation of (4.5), one has

−u(α+ ρ)I∗ + (1− ε)kI∗ − µV ∗ = 0,

i.e.,
V ∗ =

(1− ε)k − u(α+ ρ)

µ
I∗,

hence
V ∗ =

(1− ε)k − u(α+ ρ)

µ

λ− dH∗

α
(4.7)

according to (4.6). The substitution of (4.7) in (4.5) yields :
(1− η)L

(
H∗,

λ− dH∗

α
,

(1− ε)k − u(α+ ρ)

µ

λ− dH∗

α

)
(1− ε)k − u(α+ ρ)

µ
I∗ − (α+ ρ)I∗ = 0.

Thus, we have
(1− η)

(
(1− ε)k − u(α+ ρ)

)
L
(
H∗,

λ− dH∗

α
,

(1− ε)k − u(α+ ρ)

µ

λ− dH∗

α

)
= (α+ ρ)µ

since I∗ 6= 0. Furthermore, I∗ ≥ 0, gives λ−dH∗
α ≥ 0. Thus H∗ ≤ λ

d . Hence there is not a biologicalequilibrium when H∗ > λ
d .Let us consider the function ψ defined on [0, λd ] by :

ψ(x) = (1− η)γL
(
x,
λ− dx
α

,
γ(λ− dx)

µα

)
− (α+ ρ)µ,

where
γ = (1− ε)k − u(α+ ρ).

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 19We have
ψ(0) = −(α+ ρ)µ < 0

and
ψ

(
λ

d

)
= (1− η)γL

(
λ

d
, 0, 0

)
− (α+ ρ)µ,

= (1− η) [(1− ε)k − u(α+ ρ)]
βΛ

α0 + α1Λ
− (α+ ρ)µ,

=
(1− η)(1− ε)kβΛ

α0 + α1Λ
−
u(1− η)(α+ ρ)βΛ

α0 + α1Λ
− (α+ ρ)µ,

=
1

α0 + α1Λ
[(1− η)(1− ε)kβΛ− µ(α+ ρ)(α0 + α1Λ)− u(1− η)(α+ ρ)βΛ] ,

=
1

α0 + α1Λ

[
(1− η)(1− ε)kβΛ− (α+ ρ) [µ(α0 + α1Λ) + u(1− η)βΛ]

]
,

=
(α+ ρ) [µ(α0 + α1Λ) + u(1− η)βΛ]

α0 + α1Λ
(R0 − 1),

It follows that
ψ
(λ
d

)
=

(α+ ρ) [µ(α0 + α1Λ) + u(1− η)βΛ]

α0 + α1Λ
(R0 − 1) > 0 if and only if R0 > 1.

Moreover, letting y =
λ− d.x
α

and z =
γ(λ− d.x)

µα
, we have

ψ
′
(x) = (1− η)γ.

d

dx

[
L

(
x,
λ− dx
α

,
γ(λ− dx)

µα

)
− (α+ ρ)µ

]
,

= (1− η)γ

(
∂L

∂x
−
d

α

∂L

∂y
−
γd

µα

∂L

∂z

)
,

= (1− η)γ

(
∂L

∂x
−
d

α

∂L

∂x

∂x

∂y
−
γd

µα

∂L

∂x

∂x

∂z

)
,

= (1− η)γ

(
∂L

∂x
−
d

α

∂L

∂x

(
−
α

d

)
−
γd

µα

∂L

∂x

(
−
µα

γd

))
,

= 3(1− η)γ
∂L

∂x
,

= 3(1− η)γ
βα0 + βα2V

(α0 + (α1 + α3V )x + α2V )2 > 0 if γ > 0.

Therefore, if R0 > 1 there exists a unique spatially homogeneous infected equilibrium E∗ =

(H∗, I∗, V ∗) with H∗ ∈ (0, λd ), I∗ > 0 and V ∗ > 0.The previous investigations can be summarized in the following theorem :
Theorem 4.1. 1) If R0 ≤ 1, then the PDE-system (2.4) admits a unique spatially homogeneous
uninfected equilibrium E0 =

(
λ
d , 0, 0

)
.

2) If R0 > 1 and γ > 0, then the PDE system (2.4) admits a unique spatially homogeneous infected
equilibrium E∗ = (H∗, I∗, V ∗) with H∗ ∈

(
0, λd

)
, I∗ > 0 and V ∗ > 0.

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 20

Remark 4.1. Due to the spatial dependence of the state variables, spatially-inhomogeneous steadystates can exist.Indeed, any spatially-inhomogeneous equilibrium point E = (H, I, V ) of the model (2.4) subjectto the homogeneous Neumann boundary condition must solve the following system.

D1∆H + λ− dH −
(1− η)βHV

α0 + α1H + α2V + α3HV
+ ρI = 0,

D2∆I +
(1− η)βHV

α0 + α1H + α2V + α3HV
− (α+ ρ)I = 0,

D3∆V + (1− ε)kI − µV −
u(1− η)βHV

α0 + α1H + α2V + α3HV
= 0, (4.8)

∂H

∂η
=
∂I

∂η
=
∂V

∂η
= 0.

Investigation of the local stability of such spatially-inhomogeneous equilibria will be the concernof a forthcoming paper via an in-depth analysis of the above system.
4.4. Local stability of HCV-uninfected equilibrium. The objective of this section is to discuss thelocal stability of the spatially homogeneous uninfected equilibrium for the PDE system (2.4). Weaddress local stability by analysing the characteristic equation.
Theorem 4.2. The spatially homogeneous uninfected equilibrium E0 of PDE-model system (2.4) is
locally asymptotically stable if R0 ≤ 1 and it is unstable if R0 > 1.

Proof. Let {µl , ϕl} be an eigenpair of the Laplace operator −∆ on Ω with the homogeneous Neu-mann boundary condition where 0 = µ1 < µ2 < µ3 < · · ·. Let Eµl be the eigenspace correspondingto µl in C1(Ω) and {ϕl j , j = 1, 2, · · ·, dimEµl} be an orthogonal basis of Eµl . Let X = (C1(Ω))3and Xl j = {ϕl jc, / c ∈ R3}.Consider the following direct sum
X =

∞⊕
l=1

Xl with Xl =

dimEµl⊕
j=1

Xl j ,

where Xl j is the eigenspace corresponding to µl . Linearizing (2.4) at the spatially homogeneousuninfected equilibrium E0 we obtain the following linearized system :

∂w1

∂t
= D1∆w1 − dw1 + ρw2 −

(1− η)βΛ

α0 + α1Λ
w3,

∂w2

∂t
= D2∆w2 − (α+ ρ)w2 +

(1− η)βΛ

α0 + α1Λ
w3, (4.9)

∂w3

∂t
= D3∆w3 + (1− ε)kw2 −

[
µ+ u

(1− η)βΛ

α0 + α1Λ

]
w3,

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 21where W = (w1, w2, w3)T = (H, I, V )T .From the previous system (4.9), we obtain
Wt = LW = D∆W +K(E0)W

where
K(E0)W =

 −dw1 + ρw2 − (1−η)βΛ
α0+α1Λ w3

−(α+ ρ)w2 + (1−η)βΛ
α0+α1Λ w3

(1− ε)kw2 −
(
µ+ u

(1−η)βΛ
α0+α1Λ

)
w3

 . (4.10)
For each l ≥ 1, Xl is invariant under the operator L, and λ̃ is an eigenvalue of L if and only if it isan eigenvalue of the matrix −µlD + K(E0) for some l ≥ 1, in which case, there is an eigenvectorin Xl . So, one has
det

(
−µlD +K(E0)− λ̃Id

)
=

∣∣∣∣∣∣∣∣∣∣∣∣∣

−
(
µlD1 + d + λ̃

)
ρ − (1−η)βΛ

α0+α1Λ

0 −
(
µlD2 + (α+ ρ) + λ̃

)
(1−η)βΛ
α0+α1Λ

0 (1− ε)k −
(
µlD3 + µ+ u (1−η)βΛ

α0+α1Λ

)
− λ̃

∣∣∣∣∣∣∣∣∣∣∣∣∣
The characteristic equation of −µlD +K(E0) is
− (µlD1 + d + λ̃)

[[
µlD2 + (α+ ρ) + λ̃

] [(
µlD3 + µ+ u

(1− η)βΛ

α0 + α1Λ

)
+ λ̃

]
− (1− ε)(1− η)kβΛ

α0 + α1Λ

]
= 0, (4.11)

from (4.11), we get
λ̃0 = −µlD1 − d < 0,

and another characteristic eigenvalues are the roots of the following equation :
λ̃2 + Bλ̃+ (µlD2 + α+ ρ)

(
µlD3 + µ+ u

(1− η)βΛ

α0 + α1Λ

)
−

(1− ε)(1− η)kβΛ

α0 + α1Λ
= 0, (4.12)

where
B =

[(
µlD3 + µ+ u

(1− η)βΛ

α0 + α1Λ

)
+ µlD2 + (α+ ρ)

]
.

Let
C = (µlD2 + α+ ρ)

(
µlD3 + µ+ u

(1− η)βΛ

α0 + α1Λ

)
−

(1− ε)(1− η)kβΛ

α0 + α1Λ
.

One has,
C = µlD2

(
µlD3 + µ+ u

(1− η)βΛ

α0 + α1Λ

)
+ (α+ ρ)µlD3 + (α+ ρ)

(
µ+ u

(1− η)βΛ

α0 + α1Λ

)
−

(1− ε)(1− η)kβΛ

α0 + α1Λ
,

= µlD2

(
µlD3 + µ+ u

(1− η)βΛ

α0 + α1Λ

)
+ (α+ ρ)µlD3

+
1

α0 + α1Λ

[
(α+ ρ) [µ(α0 + α1Λ) + u(1− η)βΛ]− (1− ε)(1− η)kβΛ

]
,

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 22

= µlD2

(
µlD3 + µ+ u

(1− η)βΛ

α0 + α1Λ

)
+ (α+ ρ)µlD3

+
(α+ ρ) [µ(α0 + α1Λ) + u(1− η)βΛ]

α0 + α1Λ
(1−R0).

Since B > 0, if R0 ≤ 1 then C is also positive. Hence by virtue of the Routh-Hurwitz criterion,equation (4.12) does not admit solution with positive real part. Thus none characteristic eigenvaluehave positive real part. Therefore if R0 ≤ 1, the spatially homogeneous uninfected equilibrium
E0 =

(
λ
d , 0, 0

) of (2.4) is locally asymptotically stable.Otherwise if R0 > 1, then for l = 1, (in this case µ1 = 0) one has,
C =

(α+ ρ) [µ(α0 + α1Λ) + u(1− η)βΛ]

α0 + α1Λ
(1−R0) < 0.

Hence there is a complex root of equation (4.12) with positive real part in the spectrum of Kaccording to Routh-Hurwitz criterion. Therefore the uninfected equilibrium E0 =
(
λ
d , 0, 0

) of (2.4)is unstable. This completes the proof of Theorem 4.2. �

4.5. Global stability of HCV-Uninfected equilibrium. The objective of this section is to discuss theglobal stability of the spatially homogeneous uninfected equilibrium for the PDE system (2.4). Weaddress global stability by using the construction of Lyapunov functional method. This Lyapunovfunctional is obtained from those of differential equations by applying the method presented in [15].For this purpose, we start by letting
τ0 =

(1− ε)k(1− η)βΛ

µα0(α+ ρ)
.

Then, it is easy to see that
(1− ε)(1− η)kβΛ

(α+ ρ) [µ(α0 + α1Λ) + u(1− η)βΛ]
≤

(1− ε)k(1− η)βΛ

µα0(α+ ρ)
,

i.e.,
R0 ≤ τ0.We state the following result on global stability at E0 as follows :

Theorem 4.3. The spatially homogeneous uninfected equilibrium E0 of PDE-model system (2.4) is
globally asymptotically stable in the positively-invariant region Σ if τ0 < 1.

Proof. Let us consider the following function
G1(t) =

(1− ε)k

α+ ρ
I(t) + V (t).

Then, the differentiation of G1 with respect to t gives
dG1

dt
=

(
(1− ε)k − u(α+ ρ)

µ(α+ ρ)(α0 + α1H + α2V + α3HV )
(1− η)βH − 1

)
µV.

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 23Since H ≤ λ
d = Λ in the positively-invariant region Σ, one has

dG1

dt
≤

[[
(1− ε)k − u(α+ ρ)

]
(1− η)βΛ

µα0(α+ ρ)
− 1

]
µV,

≤
[(1− ε)(1− η)kβΛ

µα0(α+ ρ)
− 1
]
µV,

≤ (τ0 − 1)µV.

Now, we define the Lyapunov function as follows
L1 =

∫
Ω

G1dx.

The computation of the time derivative of L1 along the positive solutions of the PDE-model system(2.4) yields
dL1

dt
=

d

dt

[ ∫
Ω

G1dx
]
,

=

∫
Ω

dG1

dt
dx,

≤
∫

Ω

[
(τ0 − 1)µV

]
dx.

It is clear that the condition τ0 ≤ 1 gives dL1
dt ≤ 0 for all H, I, V > 0. We note that the solutionsof system (2.4) are limited by Υ, the greatest invariant subset of E={(H, I, V ) ∈ Σ|dL1

dt = 0
}.We realize that dL1

dt = 0 if and only if V = 0 and I = 0. Each element of Υ satisfies V = 0 andconsequently I=0. By Lyapunov-LaSalle invariance principle [27], E0 is globally asymptoticallystable if τ0 < 1. So, we obtain a sufficient condition R0 ≤ τ0 which ensures that the HCV spatiallyhomogeneous equilibrium E0 of PDE-model system (2.4) is globally asymptotically stable if τ0 < 1.This completes the proof of theorem 4.3. �

4.6. Local stability of HCV spatially homogeneous infected equilibrium. Let us study the localstability of the unique infected spatially homogeneous equilibrium E∗ of our PDE-model system.Consider the Laplace operator −∆ and let 0 = µ1 < µ2 < µ3 < · · · be its eigenvalues on Ω withthe homogeneous Neumann boundary condition, and Eµl be the eigenspace corresponding to µlin C1(Ω). Let also X = (C1(Ω))3, {ϕl j , j = 1, 2, · · ·, d imEµl} be an orthogonal basis of Eµl and
Xl j = {ϕl jc / c ∈ R3}.Then,

X =

∞⊕
l=1

Xl with Xl =

dimEµl⊕
j=1

Xl j .

Now, let set w1 = H, w2 = I , w3 = V . Further we use the vector notation W = (w1, w2, w3)T =

(H, I, V )T . Then the linearization of the PDE system at E∗ is of the form
Wt = LW = D∆W +K(E∗)W,

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 24where
K(E∗)W =

 − (d + A)w1 + ρw2 − Bw3

Aw1 − (α+ ρ)w2 + Bw3

−uAw1 + (1− ε)kw2 − (µ+ uB)w3

 , (4.13)
with

A =
(1− η)(α0 + α2V

∗)βV ∗

(α0 + α1H∗ + α2V ∗ + α3H∗V ∗)2and
B =

(1− η)(α0 + α1H
∗)βH∗

(α0 + α1H∗ + α2V ∗ + α3H∗V ∗)2
.

For each l ≥ 1, Xl is invariant under the operator L, and λ̃ is an eigenvalue L if and only if it isan eigenvalue of the matrix −µlD + K(E∗) for some l ≥ 1, in which case, there is an eigenvectorin Xl . Therefore we get:
det

(
− µlD+K(E∗)− λ̃Id

)
=

∣∣∣∣∣∣∣∣
−(µlD1 + d + A)− λ̃ ρ −B

A −(µlD2 + α+ ρ)− λ̃ B

−uA (1− ε)k − (µlD3 + µ+ uB)− λ̃

∣∣∣∣∣∣∣∣ .
The characteristic equation of −µlD +K(E∗) is on the form

λ̃3 + a2λ̃
2 + a1λ̃+ a0 = 0 (4.14)

where
a2 = (µlD1 + d + A+ µlD2 + α+ ρ+ µlD3 + µ+ uB) > 0,

a1 = (µlD1 +d+A)(µlD2 +α+ρ+µlD3 +µ+uB) + (µlD2 +α+ρ)(µlD3 +µ+uB)−(1−ε)kB,

a0 = (µlD1 + d + A)(µlD2 + α+ ρ)(µlD3 + µ+ uB)− (µlD1 + d + A)(1− ε)kB.If a1 > 0 and a1a2 > a0 from the above investigations, it then follows from Routh-Hurwitz criterionthat all roots of (4.14) have negative real parts and therefore we have the following result.
Theorem 4.4. If a1 > 0 and a1a2 > a0, then the spatially homogeneous infected equilibrium
E∗ = (H∗, I∗, V ∗) of the PDE-model system (2.4) is locally asymptotically stable when it exists.

4.7. Global stability of HCV-spatially homogeneous infected equilibrium. The objective of thissection is to discuss the global stability of the spatially homogeneous infected equilibrium E∗ forthe PDE system (2.4). We address global stability by using the method of construction of Lyapunovfunctionals. These Lyapunov functional is obtained from those for differential equations by applyingthe method of Hattaf and Yousfi presented in [15]. We address this study with certain assumptionsnamely : u = 0 (i.e., there is no absorption effect), α0 = 1 et α3 = α1α2. Thus we have thefollowing results.
Theorem 4.5. The spatially homogeneous infected equilibrium E∗ of PDE-model system (2.4) is
globally asymptotically stable when it exists.

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 25

Proof. We first define the function
G2(H, I, V ) = H −H∗ −

∫ H

H∗

(α+ ρ)I∗

(1−η)βτV ∗

(1+α1τ)(1+α2V ∗)

dτ + I − I∗ − I∗ ln(
I

I∗
)

+
α+ ρ

(1− ε)k
(α0 + α2V

∗)

V − V ∗ − ∫ V

V ∗

(α+ ρ)I∗

(1−η)βτH∗

(1+α1H∗)(1+α2τ)

dτ

 .
Then, the computation of the derivative of G2 with respect to t yields :

dG2

dt
=

[
λ− dH − αI −

(α+ ρ)µ

(1− ε)k
V

]
− (α+ ρ)I∗

(1 + α1H)(1 + α2V
∗)

(1− η)βHV ∗

[
λ− dH −

(1− η)βHV

(1 + α1H)(1 + α2V )
+ ρI

]
−
I∗

I

[
(1− η)βHV

(1 + α1H)(1 + α2V )
− (α+ ρ)I

]
−

α+ ρ

(1− ε)k

V ∗

V
[(1− ε)kI − µV ] .

Since
(1− η)βH∗V ∗

(1 + α1H∗)(1 + α2V ∗)
= (α+ ρ)I∗,

λ = dH∗ + αI∗,

and
(α+ ρ)µ

(1− ε)k
=

(α+ ρ)I∗

V ∗
,

we have
dG2

dt
=

[
dH∗ + αI∗ − dH − αI − (α+ ρ)I∗

V

V ∗

]
−(α+ ρ)I∗

(1 + α1H)(1 + α2V
∗)

(1− η)βHV ∗

[
dH∗ + αI∗ − dH −

(1− η)βHV

(1 + α1H)(1 + α2V )
+ ρI

]
−
I∗

I

(1− η) (1+α1H
∗)(1+α2V

∗)(α+ρ)I∗

(1−η)H∗V ∗ HV

(1 + α1H)(1 + α2V )
− (α+ ρ)I

− (α+ ρ)I
V ∗

V
+

(α+ ρ)µ

(1− ε)k
V ∗,

=

[
dH∗ + (α+ ρ)I∗ − ρI∗ − dH − αI − (α+ ρ)I∗

V

V ∗

]
−
[H∗
H

1 + α1H

1 + α1H∗
dH∗ +

H∗

H

1 + α1H

1 + α1H∗
αI∗ −

1 + α1H

1 + α1H∗
dH∗ −

V

V ∗
1 + α2V

∗

1 + α2V
(α+ ρ)I∗

+
H∗

H

1 + α1H

1 + α1H∗
ρI
]

+ (α+ ρ)I∗
[

1−
HI∗V (1 + α1H

∗)(1 + α2V
∗)

H∗IV ∗(1 + α1H)(1 + α2V )

]
+ (α+ ρ)I∗

(
1−

I

I∗
V ∗

V

)
,

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 26

=dH∗
[

1−
H

H∗
−
H∗

H

1 + α1H

1 + α1H∗
+

1 + α1H

1 + α1H∗

]
+ (α+ ρ)I∗

[
1−

HI∗V (1 + α1H
∗)(1 + α2V

∗)

H∗IV ∗(1 + α1H)(1 + α2V )
+
V

V ∗
1 + α2V

∗

1 + α2V

]
+ (α+ ρ)I∗

(
2−

V

V ∗
−
I

I∗
V ∗

V

)
− αI∗

(
ρ

α
+
I

I∗

)
−
H∗

H

1 + α1H

1 + α1H∗
αI∗ −

H∗

H

1 + α1H

1 + α1H∗
ρI,

=−
d(H −H∗)2

H(1 + α1H∗)
+ (α+ ρ)I∗

[
−1−

V

V ∗
+
V

V ∗
1 + α2V

∗

1 + α2V
+

1 + α2V

1 + α2V ∗

]
+ (α+ ρ)I∗

[
4−

H∗

H

1 + α1H

1 + α1H∗
−
HI∗V (1 + α1H

∗)(1 + α2V
∗)

H∗IV ∗(1 + α1H)(1 + α2V )
−
I

I∗
V ∗

V
−

1 + α2V

1 + α2V ∗

]
− αI∗

(
ρ

α
+
I

I∗

)
−
H∗

H

1 + α1H

1 + α1H∗
αI∗ −

H∗

H

1 + α1H

1 + α1H∗
ρI + (α+ ρ)I∗

H∗

H

1 + α1H

1 + α1H∗
.

Therefore
dG2

dt
=−

d(H −H∗)2

H(1 + α1H∗)
−

α2(α+ ρ)I∗(V − V ∗)2

V ∗(1 + α2V ∗)(1 + α2V )
− αI∗

(
ρ

α
+
I

I∗

)
+ (α+ ρ)I∗

[
4−

H∗

H

1 + α1H

1 + α1H∗
−
HI∗V (1 + α1H

∗)(1 + α2V
∗)

H∗IV ∗(1 + α1H)(1 + α2V )
−
I

I∗
V ∗

V
−

1 + α2V

1 + α2V ∗

]
.

We have :
(α+ ρ)I∗

[
4−

H∗

H

1 + α1H

1 + α1H∗
−
HI∗V (1 + α1H

∗)(1 + α2V
∗)

H∗IV ∗(1 + α1H)(1 + α2V )
−
I

I∗
V ∗

V
−

1 + α2V

1 + α2V ∗

]
≤ 0

since the left side of the latter inequality is the difference between the geometric mean and thearithmetic mean. That is dG2
dt ≤ 0. Otherwise dG2

dt = 0 if and only if H = H∗, I = I∗ et V = V ∗.Thus G2 is a Lyapunov functional of the differential equation associated to the PDE-model system(2.4). Therefore using Lyapunov-LaSalle invariance principle [27] combined to the method presentedin [15], the functional defined by
L2(t) =

∫
Ω

G2(t)dxis a Lyapunov functional of the PDE-model system (2.4) at the spatially homogeneous infectedequilibrium E∗. Therefore E∗ is globally asymptotically stable. This completes the proof of Theo-rem 4.5. �

5. Numerical simulations
In this section, we present the numerical simulations to illustrate our theoretical results. Tosimplify, we consider IBVP (2.4) with Ω = (1) under Neumann boundary condition

∂H

∂ν
= 0,

∂I

∂ν
= 0,

∂V

∂ν
= 0 t > 0, x = 1 (5.1)

and, following initial conditions
H(x, 0) = 5, I(x, 0) = 5, V (x, 0) = 5, (5.2)

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 27and
H(x, 0) = 15, I(x, 0) = 5, V (x, 0) = 5. (5.3)Now we choose the numerical values of the parameters for the PDE-cellular model system (2.4)as follows: λ = 50; d = 5; ρ = 0, 01; α = 0, 05; D1 = D2 = D3 = 0, 1; η = 0, 00004; α3 = 0, 03;

ε = 0, 5; α2 = 0, 02; k = 2; α1 = 0, 1; µ = 20; α0 = 1; β = 0, 24 et u = 1. By calculation wehave R0 = 0.943361. In this case, PDE-cellular model system (2.4) has a spatially homogeneousequilibrium E0 = (10, 0, 0). Hence by Theorem 4.3 E0 is globally asymptotically stable. Numericalsimulation illustrates our result (see FIGURE 1). Otherwise we choose the numerical values of

(a) (b) (c)

Figure 1. Simulations of IBVP (2.4) under Neumann boundary conditions (5.1) andinitial condition (5.2)
the parameters for the PDE-cellular model system (2.4) as follows : λ = 50; d = 5; ρ = 0, 01;
α = 0, 05; D1 = D2 = D3 = 0, 1; η = 0, 00004; ε = 0, 5; α0 = 1; α1 = 0, 1; α2 = 0, 02;
α3 = 0, 03; k = 2; µ = 2; β = 0, 24 et u = 1. By calculation we have R0 = 6.25009. In this case,PDE-cellular model system (2.4) has a spatially homogeneous equilibrium E∗ = (5; 500; 235).

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 28Hence by Theorem 4.5 E∗ is globally asymptotically stable. Numerical simulation illustrates ourresult (see FIGURE 2).

(a) (b) (c)

Figure 2. Simulations of IBVP (2.4) under Neumann boundary conditions (5.1) andinitial condition (5.3)
6. Conclusion

In this work, we addressed the dynamics of a reaction diffusion HCV intra-host infection modelwith the Hattaf-Yousfi incidence rate, which is a generalized non-linear incidence rate. The object ofthis work was to make a mathematical analysis of a cellular model of HCV infection which assumesthat virions diffuse into the liver, which uses the Hattaf-Yousfi functional response generalizingmost of the functional responses that exist. Our model also takes into account the absorption effectwhich is much neglected in the literature. We first showed that the initial value and boundaryproblem (2.4) admits a unique global solution in time. And secondly, we have shown that thisunique solution is positive and uniformly bounded. Then, we determined the expression of basic

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 29reproduction number R0 which is the parameter from which we studied the dynamics of our modelat equilibria whose existence and uniqueness of these have been previously proven. More precisely,we have shown that, if R0 < 1, the unique uninfected equilibrium point is locally and globallyasymptotically stable. This means that, under this condition,infection disappears. Otherwise, theuninfected equilibrium point is unstable; and in this case the infection persists in the host. It hasalso been shown under the hypothesis R0 > 1, that the infected equilibrium point is locally andglobally asymptotically stable. This edifying work ended with numerical simulations, carried outon the Mathematica software, which confirmed our theoretical results. In order to get as close aspossible to complex reality of biological phenomena, we envisage in the future, the mathematicalanalysis of models taking into account cell proliferation, the delays and more generalization of theincident rate function.
Appendix A. Proof of Proposition 3.8

The proof is established by using Banach’s Fixed Point Theorem.Choose β such that 3
4 < β < 1, then the injection I : D(Hβ) → C0

B is continuous by Lemma 3.3.For r0 > 0 and T > 0, the following closed ball is considered :
Br0 (H0, I0, V0) =

{
(H, I, V ) ∈ (C0

B((0, T ], D(Hβ)))3 : ‖H −H0‖β, ‖I − I0‖β, ‖V − V0‖β ≤ r0
}
.(A.1)By Proposition 3.1, we have the local Lipschitz properties:

‖F(t, H1, I1, V1)−F(t, H2, I2, V2)‖2 ≤ K1
1‖H1 −H2‖D(Hβ) +K1

2‖I1 − I2‖D(Hβ) +K1
3‖V1 − V2‖D(Hβ),

‖G(t, H1, I1, V1)− G(t, H2, I2, V2)‖2 ≤ K2
1‖H1 −H2‖D(Hβ) +K2

2‖I1 − I2‖D(Hβ) +K2
3‖V1 − V2‖D(Hβ),

‖Q(t, H1, I1, V1)−Q(t, H2, I2, V2)‖2 ≤ K3
1‖H1 −H2‖D(Hβ) +K3

2‖I1 − I2‖D(Hβ) +K3
3‖V1 − V2‖D(Hβ),(A.2)

for t ∈ [0, T ], (H1, I1, V1), (H2, I2, V2) ∈ Br0 (H0, I0, V0) and Lipschitz-constants K ij > 0, i , j =

1, 2, 3. In addition (y1, y2, y3) ∈ Br0 (H0, I0, V0), define P : [0, T ] → L2(Ω), Q : [0, T ] → L2(Ω)and R : [0, T ]→ L2(Ω) as follows
Py1(t) = G1(t)H0 +

∫ t

0

G1(t − τ)F(t, y1(τ), y2(τ), y3(τ))dτ,

Qy2(t) = G2(t)I0 +

∫ t

0

G2(t − τ)G(t, y1(τ), y2(τ), y3(τ))dτ,

Ry3(t) = G3(t)V0 +

∫ t

0

G3(t − τ)Q(t, y1(τ), y2(τ), y3(τ))dτ.

Finally, set M1 = sup
t∈[0,T ]

‖F(t, y1(0), y2(0), y3(0))‖2, M2 = sup
t∈[0,T ]

‖G(t, y1(0), y2(0), y3(0))‖2,
M3 = sup

t∈[0,T ]
‖Q(t, y1(0), y2(0), y3(0))‖2 and choose T so that:
‖G1(h)H0 −H0‖β, ‖G2(h)I0 − I0‖β, ‖G3(h)V0 − V0‖β ≤

3r0
4
, 0 ≤ h ≤ T, (A.3)

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 30

C iβ,2(Mi + r0K
i
1 + r0K

i
2 + r0K

i
3)

∫ T

0

s−βds ≤
r0
4
, i = 1, 2, 3, (A.4)

where C iβ,2 is the constant in property 2) of Corollary 3.5 for operator A (A is zero for i = 1, 2).Note that such T exists since G1(h), G2(h) and G3(h) converge to Id as h tends to 0+ by definitionof an analytic semigroup and∫ T

0

s−βds =
1

1− βh
1−β → 0, h → 0+ f or β < 1.

Then the proof will continue according to the following two points :
(a): it is shown that (P,Q,R) maps Br0 (H0, I0, V0) into itself,
(b): it is shown that (P,Q,R) is a strict contraction on Br0 (H0, I0, V0), allowing the use ofBanach’s Fixed Point Theorem to get the existence of a unique fixed point in Br0 (H0, I0, V0).Let (H, I, V ) ∈ Br0 (H0, I0, V0). Then, using (A.3) and property 2) of Corollary 3.5, we have :

‖PH(t)−H0‖D(Hβ) =
∥∥∥G1(t)H0 −H0 +

∫ t

0

G1(t − τ)F(τ,H(τ), I(τ), V (τ))dτ
∥∥∥
D(Hβ)

,

≤ ‖G1(t)H0 −H0‖D(Hβ) +

∫ t

0

‖G1(t − τ)F(τ,H(τ), I(τ), V (τ))‖D(Hβ)dτ,

≤
3r0
4

+

∫ t

0

C1
β,2(t − s)−β

∥∥∥F(τ,H(τ), I(τ), V (τ))
∥∥∥

2
dτ,

≤
3r0
4

+

∫ t

0

C1
β,2(t − s)−β

∥∥∥F(τ,H(τ), I(τ), V (τ))−F(τ,H0, I0, V0)

+F(τ,H0, I0, V0)
∥∥∥

2
dτ,

≤
3r0
4

+

∫ t

0

C1
β,2(t − s)−β

(∥∥∥F(τ,H(τ), I(τ), V (τ))−F(τ,H0, I0, V0)
∥∥∥

2

+
∥∥∥F(τ,H0, I0, V0)

∥∥∥
2

)
dτ,

‖PH(t)−H0‖D(Hβ) ≤
3r0
4

+

∫ t

0

C1
β,2(t − s)−β

(
K1

1 r0 +K1
2 r0 +K1

3 r0 +M1

)
dτ,

≤
3r0
4

+ C1
β,2

(
K1

1 r0 +K1
2 r0 +K1

3 r0 +M1

) ∫ t

0

(t − s)−βdτ,

≤
3r0
4

+ 0,

≤ r0 for 0 ≤ t ≤ T,

and similarly
‖QI(t)− I0‖D(Hβ) ≤ r0 and < ‖RV (t)− V0‖D(Hβ) ≤ r0.Showing that (P,Q,R) maps Br0 (H0, I0, V0) into itself. Furthermore, from property 1) in Corol-lary 3.5 and Lemma 3.6 we compute:

‖PH(t + h)− PH(t)‖D(Hβ)

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 31

=
∥∥∥G1(t + h)H0 − G1(t)H0 +

∫ t+h

0

G1(t + h − τ)F(τ,H(τ), I(τ), V (τ))dτ

−
∫ t

0

G1(t − τ)F(τ,H(τ), I(τ), V (τ))dτ
∥∥∥
D(Hβ)

,

=
∥∥∥G1(t)G1(h)H0 − G1(t)H0 +

∫ t

0

G1(h)G1(t − τ)F(τ,H(τ), I(τ), V (τ))dτ

+

∫ t+h

t

G1(h)G1(t − τ)F(τ,H(τ), I(τ), V (τ))dτ

−
∫ t

0

G1(t − τ)F(τ,H(τ), I(τ), V (τ))dτ
∥∥∥
D(Hβ)

,

=
∥∥∥G1(t)(G1(h)− Id)H0 +

∫ t

0

(G1(h)− Id)G1(t − τ)F(τ,H(τ), I(τ), V (τ))dτ

+

∫ t+h

t

G1(t + h − τ)F(τ,H(τ), I(τ), V (τ))dτ
∥∥∥
D(Hβ)

.

The strong continuity of the semigroup and Lemma 3.6 yields
‖PH(t + h)− PH(t)‖D(Hβ) → 0 as h → 0+.

Therefore P is continuous from [0,T] into D(Hβ).A similar Calculation shows that Q and R have the same properties and item (a) is proved.Presently let us show that (P,Q,R) is a strict contraction on Br0 (H0, I0, V0).Let (y1, y2, y3), (z1, z2, z3) ∈ Br0 (H0, I0, V0). Then
‖Py1(t)− Pz1(t)‖D(Hβ)

=
∥∥∥∫ t

0

G1(t − τ)
(
F(τ, y1(τ), y2(τ), y3(τ))−F(τ, z1(τ), z2(τ), z3(τ))

)
dτ
∥∥∥
D(Hβ)

,

≤
∫ t

0

∥∥∥G1(t − τ)
(
F(τ, y1(τ), y2(τ), y3(τ))−F(τ, z1(τ), z2(τ), z3(τ))

)∥∥∥
D(Hβ)

dτ,

≤
∫ t

0

C1
β,2(t − τ)−β

∥∥∥F(τ, y1(τ), y2(τ), y3(τ))−F(τ, z1(τ), z2(τ), z3(τ))
∥∥∥

2
dτ,

≤
∫ t

0

C1
β,2(t − τ)−β

(
K1

1‖y1 − z1‖D(Hβ) +K1
2‖y2 − z2‖D(Hβ) +K1

3‖y3 − z3‖D(Hβ)

)
dτ,

≤ C1
β,2

(
K1

1 sup
τ∈[0,t]

‖y1 − z1‖D(Hβ) +K1
2 sup
τ∈[0,t]

‖y2 − z2‖D(Hβ)

+ K1
3 sup
τ∈[0,t]

‖y3 − z3‖D(Hβ)

)∫ t

0

(t − τ)−βdτ.

That is :
‖Py1(t)− Pz1(t)‖D(Hβ)

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 32

≤ C1
β,2

(
K1

1 +K1
2 +K1

3 +
M1

r0

)(
sup
τ∈[0,t]

‖y1 − z1‖D(Hβ) + sup
τ∈[0,t]

‖y2 − z2‖D(Hβ)

+ sup
τ∈[0,t]

‖y3 − z3‖D(Hβ)

)∫ t

0

(t − τ)−βdτ,

≤
(

sup
τ∈[0,t]

‖y1 − z1‖D(Hβ) + sup
τ∈[0,t]

‖y2 − z2‖D(Hβ) + sup
τ∈[0,t]

‖y3 − z3‖D(Hβ)

)
×C1

β,2

(
K1

1 +K1
2 +K1

3 +
M1

r0

)∫ t

0

s−βds,

≤
(

sup
τ∈[0,t]

‖y1 − z1‖D(Hβ) + sup
τ∈[0,t]

‖y2 − z2‖D(Hβ) + sup
τ∈[0,t]

‖y3 − z3‖D(Hβ)

)
×

1

r0
C1
β,2

(
r0K

1
1 + r0K

1
2 + r0K

1
3 +M1

) ∫ T

0

s−βds,

≤
(

sup
τ∈[0,t]

‖y1 − z1‖D(Hβ) + sup
τ∈[0,t]

‖y2 − z2‖D(Hβ) + sup
τ∈[0,t]

‖y3 − z3‖D(Hβ)

) 1

r0
×
r0
4
,

≤
1

4

(
sup
τ∈[0,t]

‖y1 − z1‖D(Hβ) + sup
τ∈[0,t]

‖y2 − z2‖D(Hβ) + sup
τ∈[0,t]

‖y3 − z3‖D(Hβ)

)
,

for every t ∈ [0, T ]. Hence
sup
t∈[0,T ]

‖Py1(t)− Pz1(t)‖D(Hβ) ≤
1

4

(
sup
τ∈[0,t]

‖y1 − z1‖D(Hβ) + sup
τ∈[0,t]

‖y2 − z2‖D(Hβ) + sup
τ∈[0,t]

‖y3 − z3‖D(Hβ)

)
.

Similarly
sup
t∈[0,T ]

‖Qy1(t)−Qz1(t)‖D(Hβ) ≤
1

4

(
sup
τ∈[0,t]

‖y1 − z1‖D(Hβ) + sup
τ∈[0,t]

‖y2 − z2‖D(Hβ) + sup
τ∈[0,t]

‖y3 − z3‖D(Hβ)

)
and
sup
t∈[0,T ]

‖Ry1(t)− Rz1(t)‖D(Hβ) ≤
1

4

(
sup
τ∈[0,t]

‖y1 − z1‖D(Hβ) + sup
τ∈[0,t]

‖y2 − z2‖D(Hβ) + sup
τ∈[0,t]

‖y3 − z3‖D(Hβ)

)
.

Thus
sup

t∈[0,T ]
‖(P,Q,R)(y1(t), y2(t), y3(t))− (P,Q,R)(z1(t), z2(t), z3(t))‖D(Hβ)3

≤ sup
τ∈[0,t]

(
‖Py1 − Pz1‖D(Hβ) + ‖Qy2 −Qz2‖D(Hβ) + ‖Ry3 − Rz3‖D(Hβ)

)
,

≤ sup
τ∈[0,t]

‖Py1 − Pz1‖D(Hβ) + sup
τ∈[0,t]

‖Qy2 −Qz2‖D(Hβ) + sup
τ∈[0,t]

‖Ry3 − Rz3‖D(Hβ),

≤
3

4

(
sup
τ∈[0,t]

‖y1 − z1‖D(Hβ) + sup
τ∈[0,t]

‖y2 − z2‖D(Hβ) + sup
τ∈[0,t]

‖y3 − z3‖D(Hβ)

)
,

≤
3

4
sup

t∈[0,T ]
‖(y1(t), y2(t), y3(t))− (z1(t), z2(t), z3(t))‖D(Hβ)3 .

Hence (P,Q,R) is a strict contraction on Br0 (H0, I0, V0) and this proves Part b). According toBanach’s Fixed Point Theorem, (P,Q,R) has a unique fixed point in Br0 (H0, I0, V0). This is thesolution of (2.4) on [0,T] with initial value (H(0), I(0), V (0)) = (H0, I0, V0) in (D(Hβ))3. Thiscompletes the proof of Proposition 3.8.

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 33Appendix B. Proof of Theorem 3.11
We first prove the existence and positivity of the solution. From Theorem 3.10, there exist asequence (wmi ) and a function wi such that

wmi → wi in C0([0, T ],H),

with wi ≥ 0 and wi(0) = w0i .We check that
qi(w

m−1)→ qi(w) and qi(w
m−1)wmi → qi(w)wi in C0([0, T ],H).

We also have
fi(w

m−1)→ fi(w) in C0([0, T ],H)
⋂
L2((0, T ), E′).But, wmi is solution of〈∂wmi

∂t
, vi

〉
+
〈
Aiw

m
i , vi

〉
+
(
qi(w

m−1)wmi , vi
)

=
〈
fi(w

m−1), vi

〉
, ∀vi ∈ E. (B.1)

We take φ ∈ D((0, T )), such that φvi ∈ L2((0, T ), E),∫ T

0

〈∂wmi
∂t

, φvi

〉
dt+

∫ T

0

〈
Aiw

m
i , φvi

〉
dt+

∫ T

0

(
qi(w

m−1)wmi , φvi
)
dt =

∫ T

0

〈
fi(w

m−1), φvi

〉
dt.(B.2)The second term in the left side and the right side of the equality (B.2) converges due to the weakconvergence in L2((0, T ), E′). The third term in the left-hand side of (B.2) also converges, due tothe convergence in C([0, T ],H). We deduce that ∂wmi

∂t converges weakly in L2((0, T ), E′).But we have
wmi → wi in C0([0, T ],H).Then
∂wmi
∂t
→
∂wi
∂t

in D′((0, T ),H)Therefore, we obtain
∂wmi
∂t
→
∂wi
∂t

weakly in L2((0, T ), E′),and∫ T

0

〈∂wi
∂t
, φvi

〉
dt +

∫ T

0

〈
Aiwi , φvi

〉
dt +

∫ T

0

(
qi(w)wi , φvi

)
H
dt =

∫ T

0

〈
fi(w), φvi

〉
dt. (B.3)

This being true for all φ, one has〈∂wi
∂t
, vi

〉
+
〈
Aiwi , vi

〉
+
(
qi(w)wi , vi

)
H

=
〈
fi(w), vi

〉
, ∀vi ∈ E.

That is to say,
d

dt
(wi , vi)H + a(wi , vi) +

(
qi(w)wi , vi

)
H

=
〈
fi(w), vi

〉
, ∀vi ∈ E, (B.4)

∂wi
∂t

= fi(w)− Aiwi − qi(w)wi in L2((0, T ), E′). (B.5)

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 34According to (3.27) and (3.28) we have
wmi (t) = Gi(t)w0i +

∫ t

0

Gi(t − s)(−qi(wm−1)wmi + fi(w
m−1))(s)ds, (B.6)

and in addition, as qi(wm−1)wmi and fi(w
m−1) converge in C0([0, T ],H) and the operator Gi ,defined by the relation (3.29), is compact, using the limit in (B.6) one has,

wi(t) = Gi(t)w0i +

∫ t

0

Gi(t − s)(−qi(w)wi + fi(w))(s)ds. (B.7)
It remains to prove uniqueness.Let v be another solution of IBVP (3.21). Then

vi ∈ W (0, T, E, E′)⇒ vi ∈ C0([0, T ],H) and vi ≥ 0.

Consequently we obtain
qi(v)vi + fi(v) ∈ L2((0, T ), E′).

Thus, by Proposition 2.11 of [12], one has
vi(t) = Gi(t)w0i +

∫ t

0

Gi(t − s)(−qi(v)vi + fi(v))(s)ds.

Subtracting, we have
wi(t)− vi(t) =

∫ t

0

Gi(t − s)
(
− (qi(w)wi − qi(v)vi) + (fi(w)− fi(v))

)
(s)ds, (B.8)

with
qi(w)wi − qi(v)vi = qi(w)wi − qi(w)vi + qi(w)vi − qi(v)vi ,

= qi(w)(wi − vi) + (qi(w)− qi(v))vi .

Since wi is positive, one has∥∥∥∥∥ wj
α0 + α1wk + α2wj + α3wkwj

∥∥∥∥∥ ≤ 1

K
‖wj‖∞

where
‖wj‖∞ = ‖wj‖L∞((0,T ),H).

If we define
‖w‖∞ =

3∑
j=1

‖wj‖∞,

there is M1 > 0 such that
‖q(w)‖∞ ≤ M1‖wj‖∞.

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 35So, for r = 1, 2, 3, the numerator of qr (w) − qr (v) is the sum of terms of the form (wk − vk)vj or
(wj − vj)wk , and we can find M2 > 0 such that

∣∣qr (w)− qr (v)
∣∣
H(s) ≤ M2

( 3∑
j=1

|wj(s)− vj(s)|H
)
.

Also we can find M3 > 0 such that∣∣fr (w)− fr (v)
∣∣
H(s) ≤ M3

( 3∑
j=1

|wj(s)− vj(s)|H
)
.

Summing up |wj(s) − vj(s)|H and noting that ‖Gj(t − s)‖ ≤ NjeθjT with Nj , θj > 0, we can find
M > 0 such that

3∑
j=1

|wj(s)− vj(s)|H ≤ M‖w − v‖∞.

Replacing in (B.8), we obtain
3∑
j=1

|wj(s)− vj(s)|H ≤ M2‖w − v‖∞
∫ t

0

sds = M2 t
2

2
‖w − v‖∞.

By induction, we have
3∑
j=1

|wj(s)− vj(s)|H ≤
Mn

n!
T n‖w − v‖∞,

with
lim

n→+∞

Mn

n!
T n‖w − v‖∞ = 0.

Therefore w = v . This ends the proof of Theorem 3.11.
AcknowledgmentA.Nangue acknowledges support from the Faculty of Sciences of the University of Maroua, wherethis work was initiated.
Conflict of InterestsThe authors declare that they have no conflict of interests regarding the publication of this paper.
Authors’ ContributionA. Nangue provided the subject, wrote the introduction, the conclusion, checked the proofs andverified the calculation. He also managed the wellposedness of the initial and boundary valueproblem. B. Nde Tchiffo conceived the study and computed equilibria and their global stabilities.The two authors together undertook the numerical simulations on Mathematica. All the authorsread and approved the final manuscript.

https://doi.org/10.28924/ada/ma.3.1


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 36References
[1] R.A. Adams, Sobolev spaces, Academic Press, New York, (1975).[2] H. Amann, Dynamic theory of quasilinear parabolic equations—I. Abstract evolution equations, Nonlinear Anal.:Theory Methods Appl. 12 (1988) 895–919. https://doi.org/10.1016/0362-546x(88)90073-9.[3] H. Amann, Dynamic theory of quasilinear parabolic systems, Math Z. 202 (1989) 219–250. https://doi.org/10.

1007/bf01215256.[4] H. Amann, Dynamic theory of quasilinear parabolic equations-II. Reaction-diffusion systems, Differ. Integral Equ. 3(1990) 13-75.[5] J.R. Beddington, Mutual interference between parasites or predators and its effect on searching efficiency, J. AnimalEcol. 44 (1975) 331. https://doi.org/10.2307/3866.[6] A. Chatterjee, J. Guedj, A.S. Perelson, Mathematical modelling of HCV infection: what can it teach us in the era ofdirect-acting antiviral agents?, Antivir. Ther. 17 (2012) 1171–1182. https://doi.org/10.3851/imp2428.[7] M.S.F. Chong, M. Shahrill, L. Crossley, et al. The stability analyses of the mathematical models of hepatitis C virusinfection, Modern Appl. Sci. 9 (2015) 250-271. https://doi.org/10.5539/mas.v9n3p250.[8] S.M. Ciupe, R.M. Ribeiro, P.W. Nelson, A.S. Perelson, Modeling the mechanisms of acute hepatitis B virus infection,J. Theor. Biol. 247 (2007) 23–35. https://doi.org/10.1016/j.jtbi.2007.02.017.[9] P.H. Crowley, E.K. Martin, Functional responses and interference within and between year classes of a dragonflypopulation, J. North Amer. Benthol. Soc. 8 (1989) 211–221. https://doi.org/10.2307/1467324.[10] H. Dahari, R.M. Ribeiro, A.S. Perelson, Triphasic decline of hepatitis C virus RNA during antiviral therapy, Hepa-tology. 46 (2007) 16–21. https://doi.org/10.1002/hep.21657.[11] D.L. DeAngelis, R.A. Goldstein, R.V. O’Neill, A model for tropic interaction, Ecology. 56 (1975) 881–892. https:
//doi.org/10.2307/1936298.[12] C. Goudjo, B. Lèye, M. Sy, Weak solution to a parabolic nonlinear system arising in biological dynamic in the soil,Int. J. Differ. Equ. 2011 (2011) 831436. https://doi.org/10.1155/2011/831436.[13] J.W.H. So, S.A. Gourley, Dynamics of a food-limited population model incorporating nonlocal delays on a finitedomain, J. Math. Biol. 44 (2002) 49–78. https://doi.org/10.1007/s002850100109.[14] J. Guedj, A.U. Neumann, Understanding hepatitis C viral dynamics with direct-acting antiviral agents due to theinterplay between intracellular replication and cellular infection dynamics, J. Theor. Biol. 267 (2010) 330–340.
https://doi.org/10.1016/j.jtbi.2010.08.036.[15] K. Hattaf, N. Yousfi, Global stability for reaction–diffusion equations in biology, Computers Math. Appl. 66 (2013)1488–1497. https://doi.org/10.1016/j.camwa.2013.08.023.[16] K. Hattaf, N. Yousfi, Global stability of a virus dynamics model with cure rate and absorption, J. Egypt. Math. Soc.22 (2014) 386–389. https://doi.org/10.1016/j.joems.2013.12.010.[17] K. Hattaf, N. Yousfi, Global dynamics of a delay reaction–diffusion model for viral infection with specific functionalresponse, Comp. Appl. Math. 34 (2014) 807–818. https://doi.org/10.1007/s40314-014-0143-x.[18] K. Hattaf, N. Yousfi, A class of delayed viral infection models with general incidence rate and adaptive immuneresponse, Int. J. Dynam. Control. 4 (2015) 254–265. https://doi.org/10.1007/s40435-015-0158-1.[19] K. Hattaf, N. Yousfi, A. Tridane, Mathematical analysis of a virus dynamics model with general incidence rate andcure rate, Nonlinear Anal.: Real World Appl. 13 (2012) 1866–1872. https://doi.org/10.1016/j.nonrwa.2011.
12.015.[20] K. Hattaf, N. Yousfi, A. Tridane, A delay virus dynamics model with general incidence rate, Differ. Equ. Dyn. Syst.22 (2013) 181–190. https://doi.org/10.1007/s12591-013-0167-5.[21] K. Hattaf, N. Yousfi, A. Tridane, Stability analysis of a virus dynamics model with general incidence rate and twodelays, Appl. Math. Comput. 221 (2013) 514–521. https://doi.org/10.1016/j.amc.2013.07.005.

https://doi.org/10.28924/ada/ma.3.1
https://doi.org/10.1016/0362-546x(88)90073-9
https://doi.org/10.1007/bf01215256
https://doi.org/10.1007/bf01215256
https://doi.org/10.2307/3866
https://doi.org/10.3851/imp2428
https://doi.org/10.5539/mas.v9n3p250
https://doi.org/10.1016/j.jtbi.2007.02.017
https://doi.org/10.2307/1467324
https://doi.org/10.1002/hep.21657
https://doi.org/10.2307/1936298
https://doi.org/10.2307/1936298
https://doi.org/10.1155/2011/831436
https://doi.org/10.1007/s002850100109
https://doi.org/10.1016/j.jtbi.2010.08.036
https://doi.org/10.1016/j.camwa.2013.08.023
https://doi.org/10.1016/j.joems.2013.12.010
https://doi.org/10.1007/s40314-014-0143-x
https://doi.org/10.1007/s40435-015-0158-1
https://doi.org/10.1016/j.nonrwa.2011.12.015
https://doi.org/10.1016/j.nonrwa.2011.12.015
https://doi.org/10.1007/s12591-013-0167-5
https://doi.org/10.1016/j.amc.2013.07.005


Eur. J. Math. Anal. 10.28924/ada/ma.3.1 37

[22] D. Henry, Geometric theory of semilinear parabolic equations, 1st edition, Springer-Verlag, New York, (1981).[23] D. Henry, Geometric theory of semilinear parabolic equations, 840 of Lecture Notes in Mathematics, New York,USA (1993).[24] S.L. Hollis, R.H. Martin, Jr., M. Pierre, Global existence and boundedness in reaction-diffusion systems, SIAM J.Math. Anal. 18 (1987) 744–761. https://doi.org/10.1137/0518057.[25] G. Huang, W. Ma, Y. Takeuchi, Global properties for virus dynamics model with Beddington–DeAngelis functionalresponse, Appl. Math. Lett. 22 (2009) 1690–1693. https://doi.org/10.1016/j.aml.2009.06.004.[26] G. Huang, W. Ma, Y. Takeuchi, Global analysis for delay virus dynamics model with Beddington–DeAngelis functionalresponse, Appl. Math. Lett. 24 (2011) 1199–1203. https://doi.org/10.1016/j.aml.2011.02.007.[27] H. Khalil, Nonlinear systems, 3rd edn. Prentice Hall, New York, (2002).[28] D. Li, W. Ma, Asymptotic properties of a HIV-1 infection model with time delay, J. Math. Anal. Appl. 335 (2007)683–691. https://doi.org/10.1016/j.jmaa.2007.02.006.[29] L. Min, Y. Su, Y. Kuang, Mathematical analysis of a basic model of virus infection with application to HBV infection,Rocky Mount. J. Math. 38 (2008) 1573-1585. https://www.jstor.org/stable/44239519.[30] A.U. Neumann, N.P. Lam, H. Dahari, et al. Hepatitis C viral dynamics in vivo and the antiviral efficacy of interferon-αtherapy, Science. 282 (1998) 103–107. https://doi.org/10.1126/science.282.5386.103.[31] M.A. Nowak, C.R.M. Bangham, Population dynamics of immune responses to persistent viruses, Science. 272 (1996)74–79. https://doi.org/10.1126/science.272.5258.74.[32] M.A. Nowak, S. Bonhoeffer, A.M. Hill, et al. Viral dynamics in hepatitis B virus infection., Proc. Natl. Acad. Sci.U.S.A. 93 (1996) 4398–4402. https://doi.org/10.1073/pnas.93.9.4398.[33] M.H. Protter, H.P. Weinberger, Maximum principle in differential equations, Prentice-Hall, Englewood Cliffs, NewJersey, USA (1967).[34] D. Riad, K. Hattaf, N. Yousfi, Dynamics of capital-labour model with Hattaf-Yousfi functional response, Br. J. Math.Computer Sci. 18 (2016) 1–7. https://doi.org/10.9734/bjmcs/2016/28640.[35] L. Rong, A.S. Perelson, Mathematical analysis of multiscale models for hepatitis C virus dynamics under therapy withdirect-acting antiviral agents, Math. Biosc. 245 (2013) 22–30. https://doi.org/10.1016/j.mbs.2013.04.012.[36] X. Song, A.U. Neumann, Global stability and periodic solution of the viral dynamics, J. Math. Anal. Appl. 329 (2007)281–297. https://doi.org/10.1016/j.jmaa.2006.06.064.[37] H.R. Thieme, Spectral bound and reproduction number for infinite-dimensional population structure and time het-erogeneity, SIAM J. Appl. Math. 70 (2009) 188–211. https://doi.org/10.1137/080732870.[38] X. Wang, Y. Tao, X. Song, Global stability of a virus dynamics model with Beddington–DeAngelis incidence rate andCTL immune response, Nonlinear Dyn. 66 (2011) 825–830. https://doi.org/10.1007/s11071-011-9954-0.[39] K. Wang, W. Wang, Propagation of HBV with spatial dependence, Math. Biosci. 210 (2007) 78–95. https://doi.
org/10.1016/j.mbs.2007.05.004.[40] W. Wang, X.Q. Zhao, Basic reproduction numbers for reaction-diffusion epidemic models, SIAM J. Appl. Dyn. Syst.11 (2012) 1652–1673. https://doi.org/10.1137/120872942.[41] WHO, Global hepatitis report, World Health Organization, Geneva, 2017. http://apps.who.int/iris/

bitstream/10665/255016/1/9789241565455-eng.pdf, (accessed 29 January 2018).[42] Y. Zhang, Z. Xu, Dynamics of a diffusive HBV model with delayed Beddington–DeAngelis response, Nonlinear Anal.:Real World Appl. 15 (2014) 118–139. https://doi.org/10.1016/j.nonrwa.2013.06.005.[43] X. Zhou, J. Cui, Global stability of the viral dynamics with Crowley martin functional response, Bull. Korean Math.Soc. 48 (2011) 555-574. https://doi.org/10.4134/BKMS.2011.48.3.555.

https://doi.org/10.28924/ada/ma.3.1
https://doi.org/10.1137/0518057
https://doi.org/10.1016/j.aml.2009.06.004
https://doi.org/10.1016/j.aml.2011.02.007
https://doi.org/10.1016/j.jmaa.2007.02.006
https://www.jstor.org/stable/44239519
https://doi.org/10.1126/science.282.5386.103
https://doi.org/10.1126/science.272.5258.74
https://doi.org/10.1073/pnas.93.9.4398
https://doi.org/10.9734/bjmcs/2016/28640
https://doi.org/10.1016/j.mbs.2013.04.012
https://doi.org/10.1016/j.jmaa.2006.06.064
https://doi.org/10.1137/080732870
https://doi.org/10.1007/s11071-011-9954-0
https://doi.org/10.1016/j.mbs.2007.05.004
https://doi.org/10.1016/j.mbs.2007.05.004
https://doi.org/10.1137/120872942
http://apps.who.int/iris/bitstream/10665/255016/1/9789241565455- eng.pdf
http://apps.who.int/iris/bitstream/10665/255016/1/9789241565455- eng.pdf
https://doi.org/10.1016/j.nonrwa.2013.06.005
https://doi.org/10.4134/BKMS.2011.48.3.555

	1. Introduction
	2. Formulation of the PDE-cellular model
	2.1. Fluctuation of healthy hepatocytes
	2.2. Fluctuation of HCV infected cells
	2.3. Fluctuation of free HCV virions
	2.4. The initial boundary value problem associated to PDE-cellular model

	3. Qualitative and quantitative analysis and some properties of the solutions for IBVP (2.4) 
	3.1. Local existence and uniqueness of solutions for the IBVP (2.4)
	3.2. Boundedness of the solutions for IBVP (2.4)
	3.3. Global existence, uniqueness and positivity for the IBVP (2.4)

	4. Stability analysis of the spatially homogeneous equilibria
	4.1. HCV-spatial homogeneous uninfected equilibrium E0
	4.2. Basic reproduction number R0
	4.3. Existence and uniqueness of HCV-spatial homogeneous infected equilibrium E*
	4.4. Local stability of HCV-uninfected equilibrium
	4.5. Global stability of HCV-Uninfected equilibrium
	4.6. Local stability of HCV spatially homogeneous infected equilibrium
	4.7. Global stability of HCV-spatially homogeneous infected equilibrium

	5. Numerical simulations
	6. Conclusion
	Appendix A. Proof of Proposition 3.8
	Appendix B. Proof of Theorem 3.11
	References

