Xinyue Fan
email: fan.xinyue@163.com

Claude-Michel Brauner
email: cmbrauner@gmail.com

Linda Wittkop
email: linda.wittkop@isped.u-bordeaux2.fr

MATHEMATICAL ANALYSIS OF A HIV MODEL WITH QUADRATIC LOGISTIC GROWTH TERM

Keywords: 2000 Mathematics Subject Classification. Primary: 34A34; Secondary: 34D20, 92C50 HIV modeling, stability, Hopf bifurcation, orbits, snakes

. Numerical results are presented.

1. Introduction. CD4 + T lymphocytes (CD4 + cells) are the main target for the Human Immunodeficiency Virus (HIV) [START_REF] Douek | HIV preferentially infects HIV-specific CD4+ T cells[END_REF]. During acute infection an initial peak of HIV RNA (viral load) and a massive depletion in CD4 + cell count is observed and after 3-6 months of infection the viral load and CD4 + cell count reach a more stable level (steady state) [START_REF] Ho | Rapid turnover of plasma virions and CD4 lymphocytes in HIV-1 infection[END_REF]. Without antiretroviral treatment the asymptotic chronic state lasts for around 10 years [START_REF] Babiker | Time from HIV-1 seroconversion to AIDS and death before widespread use of highlyactive antiretroviral therapy: a collaborative re-analysis[END_REF] before an increase in viral load and an accelerated decrease in CD4 + cell count is observed. Ultimately, the progressive destruction of the immune system leads to the so-called Acquired Immunodeficiency Syndrome (AIDS) and death (see Fig. 1).

The immune response and the mechanisms for T-cell depletion during chronic infection are still not completely understood [START_REF] Mcmichael | The immune response during acute HIV-1 infection: clues for vaccine development[END_REF]. However, numerous mathematical models based on ordinary differential equations have been proposed to study the HIV-host interaction in primary HIV infection and helped to understand the HIV dynamics ( [START_REF] Ho | Rapid turnover of plasma virions and CD4 lymphocytes in HIV-1 infection[END_REF][START_REF] Perelson | Decay characteristics of HIV-1-infected compartments during combination therapy[END_REF][START_REF] Perelson | HIV-1 dynamics in vivo: virion clearance rate, infected cell life-span, and viral generation time[END_REF] ).

[31] for a review). Further, mathematical models helped to gain biological insights of viral load and CD4 + cell count evolution after treatment start [START_REF] Perelson | Modelling viral and immune system dynamics[END_REF][START_REF] Perelson | Decay characteristics of HIV-1-infected compartments during combination therapy[END_REF] and helped to understand the cellular viral reservoir, e.g. HIV infected memory CD4 + cells whose half life was estimated to be 44 months [START_REF] Finzi | Latent infection of CD4+ T cells provides a mechanism for lifelong persistence of HIV-1, even in patients on effective combination therapy[END_REF][START_REF] Jones | Transient viremia, plasma viral load, and reservoir replenishment in HIVinfected patients on antiretroviral therapy[END_REF][START_REF] Siliciano | Long-term follow-up studies confirm the stability of the latent reservoir for HIV-1 in resting CD4+ T cells[END_REF]. Other cell types, e.g. monocytes or dendritic cells and anatomical sites, e.g. the central nervous system, may also constitute a reservoir for HIV (for a review, see [START_REF] Blankson | The challenge of viral reservoirs in HIV-1 infection[END_REF]).

In this article we will study a basic model for HIV-host interaction (1)-(3) including three variables [START_REF] Perelson | Mathematical analysis of HIV-I: dynamics in vivo[END_REF]. Let T denote the concentration of uninfected CD4 + T cells, I the concentration of infected CD4 + cells and V the concentration of free virus particles. For simplicity, we will call V the virus throughout the paper. The system reads:

dT dt = α -µ T T + rT ( 1 - T T max ) -γVT, (1)
dI dt = γVT -µ I I, ( 2
)
dV dt = Nµ I I -µ V V. ( 3
)
This model is based on the hypotheses that uninfected CD4 + T cells are produced from precursors in the thymus by a constant rate α and that they die by a constant rate µ T . Free virus infects CD4 + cells by a rate proportional to their abundance γVT and these infected CD4 + cells will die at a rate µ I . One infected CD4 + cell will produce N free virions during its lifespan [START_REF] Blankson | The challenge of viral reservoirs in HIV-1 infection[END_REF] and the term µ V represents the clearance rate of the free virus. In healthy individuals the CD4 + T-cell count is around 1000 cells/mm 3 and by homeostatic processes the body tries to maintain the T-cell count [START_REF] Jameson | Maintaining the norm: T-cell homeostasis[END_REF][START_REF] Jameson | T cell homeostasis: keeping useful T cells alive and live T cells useful[END_REF]. We introduce a quadratic logistic term as previously proposed [START_REF] Perelson | Dynamics of HIV infection of CD4+ T cells[END_REF][START_REF] Perelson | Mathematical analysis of HIV-I: dynamics in vivo[END_REF] to consider the homeostatic process for the CD4 + T-cell count. We thus assume that the uninfected CD4 + cells proliferate logistically as it is very likely that the proliferation will shut off once the normal/maximal level of CD4 + T-cell count is reached. The maximal CD4 + T-cell proliferation rate r (see [START_REF] Antia | The role of models in understanding CD8+ T-cell memory[END_REF]) would be obtained by the absence of a population limit. T max represents the normal/maximal density of CD4 + T cells at which CD4 + T-cell proliferation shuts off. Mathematically it will play the role of a large perturbation parameter. We refer to Table 1 for a gamut of numerical data (see, e.g., [30][4]). Typical evolution profiles are shown in Fig. 2. In the absence of virus, the T-cell population has a steady state value:

T 0 = (r -µ T ) + √ (r -µ T ) 2 + 4rα T max 2r T max .
Note that the quantity r -µ T , the net T-cell proliferation rate, needs not to be positive (see [30, p. 86]). We set X(t) = (T (t), I(t), V(t)) and write the system (1)-( 3) as the smooth Cauchy problem:

dX dt (t) = F (X(t)), X(0) = X 0 , (4)
where X 0 corresponds to non-negative initial conditions (see App. B).

It is easily seen that the system (4) has two possible nonnegative equilibria: (i) the uninfected steady state X u = (T u , I u , V u ):

T u = T 0 = (r -µ T ) + √ (r -µ T ) 2 + 4rα T max 2r T max , I u = 0, V u = 0, (5)
which corresponds to a positive equilibrium in case of no infection; Here N = 1000, r = 0.2. The other quantities are as in Table 1.

(ii) the infected steady state X i = (T i , I i , V i ):

T i = µ V γN , I i = α µ I - µ T µ V Nγµ I + rµ V Nγµ I ( 1 - µ V NγT max ) , V i = αN µ V - µ T γ + r γ ( 1 - µ V NγT max ) . (6)
The infected steady state or "seropositivity steady state" corresponds to a positive equilibrium in case of infection.

The basic reproduction number R 0 for the above described model reads:

R 0 = γNT 0 µ V . (7)
It denotes the average number of infected T cells derived from one infected T cell and is based on the definition of R 0 by epidemiologists [START_REF] Diekmann | Mathematical epidemiology of infectious diseases[END_REF]. As it can easily be seen, R 0 depends on the logistic parameter r, whereas its critical value is equal to unity for all r. R 0 was often considered as a bifurcation parameter to study stability properties of the two above described equilibria [START_REF] Leenheer | Virus dynamics: a global analysis[END_REF][START_REF] Nowak | Virus dynamics[END_REF]. De Leenheer and Smith [START_REF] Leenheer | Virus dynamics: a global analysis[END_REF] found sustained oscillations and that the steady state may be unstable when considering a system as described by ( 1)-( 3) for R 0 > 1. Perelson, Kirschner and De Boer also observed sustained oscillations in a model including a compartment of latently infected T cells [START_REF] Perelson | Dynamics of HIV infection of CD4+ T cells[END_REF]. However, their biological meaning remains unclear thus far. Perelson et al. found that N, the number of free virions produced by one actively infected CD4 + cell during its lifespan, needs to be above some critical level for successful HIV infection [START_REF] Perelson | Dynamics of HIV infection of CD4+ T cells[END_REF]. Therefore, in Section 2, we consider N and r as two independent parameters to study the properties of the "uninfected" and "infected" steady states. Simultaneously, T max will play the role of a large perturbation parameter.

First, we determine two domains in the parameter space, respectively (U) referring to "uninfected" (R 0 < 1) and (I) referring to "infected" (R 0 > 1), see Fig. 3. The uninfected equilibrium is asymptotically stable in (U) and unstable in (I). For determining these two domains we use the classical Routh-Hurwitz criterion. This criterion ensures that the eigenvalues of the Jacobian matrix have negative real part if and only if the corresponding Hurwitz determinants D j , j = 1, 2, 3, are positive. Note that the positive infected equilibrium exists in (I) only.

Second, there is an unbounded subdomain (P) in (I), the symbol (P) referring to "perturbation" (see below). In this subdomain (P): (i) the positive infected equilibrium becomes unstable whereas it is asymptotically stable in the rest of (I), and (ii) there are asymptotically stable periodic orbits with asymptotic phase. The boundary of (P) is precisely the locus of points (N, r) such that D 2 = 0, explicitly parameterized by the mappings r 1 (N) and r 2 (N). When the logistic term rT (1-T/T max ) is replaced by rT (1-(T + I)/T max ), numerical evidence shows that (P) is bounded (see App. A).

The subdomain (P) may be biologically interpreted as a perturbation of the infection by a specific or unspecific immune response against HIV. An infection with other pathogens could lead to a non-specific immune response against HIV or existing specific memory T cells could be stimulated by cross-reactive stimulation [START_REF] Antia | The role of models in understanding CD8+ T-cell memory[END_REF]. Further, interventions to boost the immune response such as Interleukin-7 or Interleukin-2 treatment could lead to containment of HIV infection [START_REF] Levy | Immunological and virological efficacy of a therapeutic immunization combined with interleukin-2 in chronically HIV-1 infected patients[END_REF][START_REF] Levy | Enhanced T cell recovery in HIV-1-infected adults through IL-7 treatment[END_REF].

Starting from Section 3, we fix N and take r as bifurcation parameter, which is logic as it is one of the crucial parameters in the host-pathogen interaction. We prove via a linear analysis that the infected equilibria X i is unstable in the interior of the region (P) and Hopf bifurcations occur on the boundary, namely at the points r 1 and r 2 .

The nonlinear analysis is fully developed in Section 4, which constitutes the core of the paper. We follow [START_REF] Hassard | Theory and applications of Hopf bifurcation[END_REF][START_REF] Kuznetsov | Elements of applied bifurcation theory[END_REF] to compute the first Lyapunov coefficient and determine the direction and stability of the bifurcated branch of periodic orbits at the Hopf points (X i (r 1 ), r 1 ) and (X i (r 2 ), r 2 ) . However, a relevant analysis is based on an asymptotic approach when T max is large, especially the observation that the system (1)-(3) has a limit when T max → +∞, in which a new bifurcation parameter is the net T-cell proliferation rate. In other words, the perturbation by the small parameter 1/T max is regular at r 1 and we take advantage of it.

The last issue we address is the existence of a "bridge" between the family of periodic orbits at the two Hopf points (X i (r 1 ), r 1 ) and (X i (r 2 ), r 2 ). We are able to prove the existence of a "snake", a maximal path of orbits, thanks to the theory of Mallet-Paret and Yorke [START_REF] Mallet-Paret | Snakes: Oriented families of periodic orbits, their sources, sinks, and continuation[END_REF].

Finally, we illustrate the results of this paper with a series of numerical computations.

Stability.

2.1. Equilibria. The notations are those of Section 1. More precisely, we state: (i) N > 0 and r ≥ 0 are parameters;

(ii) the quantities (with associated dimension) µ T , µ I , µ V , α and γ are fixed positive numbers throughout the paper;

(iii) T max is a large perturbation parameter, larger than any finite combination of µ T , µ I , µ V , α and γ of the same dimension (mm -3 ). In other words, we will make the hypothesis that

T max ≥ T 0 max (µ T , µ I , µ V , α, γ), (8)
where T 0 max is fixed and large, and depends only of µ T , µ I , µ V , α and γ. In particular, this hypothesis contains the condition T max > α/µ T of [30, p. 85].

We consider the possible nonnegative equilibria of the system (1)- [START_REF] Blankson | The challenge of viral reservoirs in HIV-1 infection[END_REF] given by the formulas ( 5) and [START_REF] Leenheer | Virus dynamics: a global analysis[END_REF]. At fixed r ≥ 0, it is not difficult to compute the critical value at which

V u = V i = 0: N crit (r) = µ V ( -(r -µ T ) + √ (r -µ T ) 2 + 4αr T max ) 2αγ = µ V γT 0 . (9)
Obviously, N = N crit (r) is equivalent to R 0 (N, r) = 1, where (see [START_REF] Dellwo | On the birth of isolas[END_REF]):

R 0 (N, r) = γNT 0 (r) µ V . ( 10
) At r = 0, N crit (0) = µ T µ V αγ
, and, as r → +∞,

N crit (r) → N crit (+∞) = µ V γT max .
It will follow from the proof of next theorem that the roles of parameters N and r can be inverted. A more convenient critical value is then defined by:

r crit (N) = (µ T µ V -γαN) + µ V ( 1 -µ V NγT max ) , N > N crit (+∞) = µ V γT max . ( 11
)
Theorem 2.1. We assume [START_REF] Diekmann | Mathematical epidemiology of infectious diseases[END_REF]. In the first quarter plane of the (N, r) parameter space, there are two regions, (U) for uninfected (R 0 (N, r) < 1), and (I) for infected (R 0 (N, r) > 1) such that: (i) in (U), the uninfected steady state X u = (T u , I u , V u ) is the only nonnegative equilibrium. In (I), there are two nonnegative equilibria, respectively X u and the infected steady state

X i = (T i , I i , V i ); (ii) it holds 0 < T u < T max and, in (I), 0 < T i < T max ; (iii) in (U)
, the uninfected steady state is asymptotically stable, in (I), the uninfected steady state is unstable.

Proof. (i) It is an exercise to see that the mapping r → N crit (r) is decreasing and convex for r > 0, since N ′′ crit (r) > 0, N ′ crit (0) < 0, and N ′ crit (r) → 0 as r → +∞. Therefore we can split the quarter plane N > 0, r ≥ 0 in two domains, first (I) which lies strictly above the graph of mapping r → N crit (r), second (U) which lies below the graph. (ii) Thanks to [START_REF] Diekmann | Mathematical epidemiology of infectious diseases[END_REF], on the one hand

T u = (r -µ T ) + √ (r -µ T ) 2 + 4rα T max 2r T max < (r -µ T ) + √ (r -µ T ) 2 + 4µ T r 2r T max = T max .
On the other hand,

T i = µ V /(γN) < T max whenever N > N crit (r) = µ V /(γT max ). (iii)
The study of the stability of the uninfected solution is easy, however we include it for the convenience of the reader. Let r ≥ 0 be fixed. We compute explicitly the eigenvalues of the Jacobian matrix

J(X u ) =           r -µ T -2rT u T max 0 -γT u 0 -µ I γT u 0 Nµ I -µ V           .
They are all real:

λ 1 = - √ (-r + µ T ) 2 + 4 rα T max < 0, λ 2 = - 1 2
( (µ I + µ V ) + √ (µ I + µ V ) 2 -4(µ I µ V -Nγµ I T u ) ) < 0, λ 3 = 1 2          -(µ I + µ V ) + √ (µ I + µ V ) 2 -4µ I µ V ( 1 - N N crit (r) )          .
It is clear that λ 3 < 0 for N < N crit (r), therefore the three eigenvalues are negative. As long there is no eigenvalue on the imaginary axis, it is well-known that linear stability yields nonlinear stability for a smooth system of the form (4), since, by the Hartman-Grobman theorem (see, e.g. [START_REF] Hartman | Ordinary Differential Equations[END_REF]IX]), near a hyperbolic equilibrium the system is locally topologically equivalent to its linearization.

N r N crit (r) µ V /(γ T max ) µ T µ V /(αγ) (U) (I) R 0 > 1 R 0 < 1 Figure 3. Profile of the curve r → N crit (r) (i.e. R 0 = 1)
which defines the two domains (U) and (I). With the values of Table 1, N crit decreases from µ T µ V /αγ = 666.67 to µ V /γT max = 6.67. See [START_REF] Fiedler | Global bifurcation of periodic solutions with symmetry[END_REF] for the reverse mapping N → r crit (N).

The issue of the stability of the infected solution is obviously more complicated. For (N, r) ∈ (I), we consider the Jacobian matrix:

J(X i ) =           r -µ T -2rT i T max -γV i 0 -γT i γV i -µ I γT i 0 Nµ I -µ V           . ( 12
)
The spectrum of J(X i ) is given by the roots (three real roots or one real and two complex conjugates) of the following characteristic equation:

λ 3 + d 1 λ 2 + d 2 λ + d 3 = 0, (13)
where the coefficients (functions of N and r) read: [START_REF] Gantmakher | The theory of matrices[END_REF]. We recall that the Routh-Hurwitz criterion gives necessary and sufficient conditions for all the roots of a real polynomial to have negative real parts (see, e.g., [START_REF] Gantmakher | The theory of matrices[END_REF]Chap. XV]). Therefore it provides a stability analysis only by examining the coefficients of the characteristic polynomial, without computing the eigenvalues. In the case of ( 13), the roots have negative real part if and only if the corresponding Hurwitz determinants are positive:

d 1 = µ I + µ V + γαN µ V + rµ V γT max N , d 2 = γNαµ I µ V + γαN + rµ I µ V γT max N + rµ 2 V γT max N , ( 14
)
d 3 =       µ I µ V - µ 2 V µ I γT max N       r + γNαµ I -µ V µ I µ T .

Study of the characteristic equation

D 1 = d 1 > 0, D 2 = d 1 d 2 -d 3 > 0, D 3 = d 3 D 2 > 0.
In other words, the roots of ( 13) have negative real part if and only if

d 1 > 0, d 3 > 0, D 2 = d 1 d 2 -d 3 > 0. (15)
Obviously, d 1 > 0. We compute d 3 = 0, which is equivalent to

αγN 2 + (rµ V -µ V µ T )N - rµ 2 V γT max = 0,
whose positive solution is N crit (r). Therefore d 3 > 0 for N > N crit (r), i.e., whenever (N, r) ∈ (I). Next, we compute the leading Hurwitz determinant

D 2 = d 1 d 2 -d 3 = 1 µ 2 V (T max ) 2 γ 2 N 2 (Ar 2 + Br + C) (16)
as a function of the parameters N and r. The coefficients read:

A =µ 5 V + µ 4 V µ I , B = -NγT max µ 2 V (Nγµ I µ V T max -2γαNµ V -2Nγµ I α -µ 3 V -3µ 2 V µ I -µ 2 I µ V ), C =N 2 γ 2 (T max ) 2 (Nγµ 3 V α + Nγµ I αµ 2 V + Nγµ 2 I µ V α + γ 2 N 2 α 2 µ V + γ 2 N 2 α 2 µ I + µ 3 V µ I µ T ).
The issue is to determine the subdomains of (I) where D 2 (N, r) is respectively positive and negative and, as far as possible, the locus of points where D 2 (N, r) = 0.

Obviously, D 2 has the sign of Ar 2 + Br + C. Note that A > 0 is fixed, only B and C depend on N. Whereas C > 0, B(N) vanishes at

N 0 = µ V (3µ I µ V + µ 2 I + µ 2 V ) γ(µ I µ V T max -2µ V α -2µ I α) . ( 17
)
Since T max is large with respect to the other quantities, N 0 is positive but small:

N crit (+∞) = µ V γT max < N 0 < N crit (0).
Therefore, for N crit (+∞) ≤ N < N 0 , we have B > 0, and B < 0 for N > N 0 .

To determine the sign of the polynomial Ar 2 + Br + C, we compute the discriminant

∆(N) = B 2 -4AC = N 2 γ 2 (T max ) 2 µ 5 V (aN 2 + bN + c), which in turn is of the sign of aN 2 + bN + c.
The coefficients a, b and c read:

a = γ 2 T max µ I (-4µ I α + µ I T max µ V -4µ V α) , b = -2γµ I µ V ( T max µ 2 I + 3µ I T max µ V -4µ I α + µ 2 V T max -4µ V α ) , c = µ V ( µ 4 V + 6µ 3 V µ I + 11µ 2 V µ 2 I -4µ 2 V µ I µ T + 6µ V µ 3 I -4µ 2 I µ V µ T + µ 4 I ) .
For T max large enough, a > 0 and b < 0. However, the sign of c, which does not contain T max , depends on the other quantities. Next we compute:

δ = b 2 -4ac = 16γ 2 µ I (µ I + µ V ) µ V (T max αµ 4 I + 4T max αµ V µ 3 I + 5µ 2 I µ 2 V T max α + µ 2 I µ 2 V (T max ) 2 µ T + 4µ 2 I µ V α 2 -4µ 2 I µ V T max αµ T + 4µ I µ 3 V T max α + 4µ I µ 2 V α 2 -4µ I T max αµ 2 V µ T + µ 4 V αT max ).
For large T max , δ has the sign of the coefficient of (T max ) 2 , hence δ > 0.

The roots of ∆ = 0, namely those of aN 2 + bN + c = 0, are N 1 < N 2 :

N 1 = -b - √ δ 2a , N 2 = -b + √ δ 2a . ( 18
)
It is clear that N 2 > 0, while N 1 has the sign of c. However, N 1 plays a role only if larger than N crit (+∞). We also remark that N 2 is small for large T max , therefore N 2 < N crit (0). It is easy to determine the position of N 0 , see [START_REF] Jameson | Maintaining the norm: T-cell homeostasis[END_REF], with respect to N 1 and

N 2 . Actually, since B(N 0 ) = 0, ∆(N 0 ) = -4AC < 0, hence aN 2 0 + bN 0 + c < 0. Therefore it holds: Lemma 2.2. N 1 < N 0 < N 2 < N crit (0).
Now we are in position to begin our discussion following the position of N.

Case 1. It corresponds to N 1 < N < N 2 . In such a case, aN 2 + bN + c < 0, hence ∆(N) < 0. Obviously, Ar 2 + Br + C has the sign of A which is positive. Case 2. It corresponds to N < N 1 or N > N 2 . (i) If N < N 1 , then, by Lemma 2.2, B > 0. Hence, Ar 2 + Br + C > 0 for any r ≥ 0 since A, B, C > 0. It thus follows that D 2 > 0. (ii) If N > N 2 ,
then, by Lemma 2.2, B < 0. Since ∆ > 0, the equation Ar 2 + Br + C = 0 admits the two real and positive roots

r 1 (N) = -B - √ ∆ 2A , r 2 (N) = -B + √ ∆ 2A . (19)
Consequently, the Hurwitz determinant D 2 is positive for r crit (N) ≤ r < r 1 (N) and r > r 2 (N), it vanishes at r = r 1 (N) and r = r 2 (N), it is negative for r 1 (N) < r < r 2 (N). Remark 1. It is important to see that r crit (N) < r 1 (N) whenever N ≥ N 2 . For N 2 ≤ N ≤ N crit (0), if r crit (N) = r 1 (N), then d 3 and D 2 vanish simultaneously which is impossible.

When N > N crit (0), r crit (N) = 0 and r 1 (N) > 0.

We are now in position to define the subdomain (P) of (I) by

(P) = {(N, r), N ≥ N 2 , r 1 (N) ≤ r ≤ r 2 (N)} , (20)
with

N 2 = -b + √ δ 2a , r 1 (N) = -B - √ ∆ 2A , r 2 (N) = -B + √ ∆ 2A ,
and at N = N 2 , the endpoint of (P), r 1 = r 2 = -B/2A.

Theorem 2.3. Under the hypothesis (8), the infected equilibrium X i = (T i , I i , V i ) is asymptotically stable for (N, r) ∈ (I\P).

Proof. Under the hypothesis (8), the leading Hurwitz determinant D 2 (N, r) is respectively negative in the interior of the subdomain (P) of (I), positive in (I\P), and it vanishes on the boundary of (P). Then, the theorem follows from the Ruth-Hurwitz criterion.

However, we can not infer from Theorem 2.3 that the infected equilibrium X i is unstable in (P) or in its interior. This issue will be addressed in Section 3. 3. Hopf bifurcation and instability. From now on, we take the logistic parameter r as a bifurcation parameter, whereas N > 0 is fixed. In Subsection 2.2, we already used the fact that N 2 , the endpoint of the subdomain (P), is small for large T max . More precisely, it is easy to see from [START_REF] Jameson | T cell homeostasis: keeping useful T cells alive and live T cells useful[END_REF] that

N 2 = µ 2 I + 3µ I µ V + µ 2 V + 2 √ µ I µ V µ T (µ I + µ V ) γµ I 1 T max + o( 1 T max
).

We assume hereafter that

T max ≥ T 1 max (µ T , µ I , µ V , α, γ, N), ( 21
)
where [START_REF] Diekmann | Mathematical epidemiology of infectious diseases[END_REF]) is chosen large enough such that 0 < N 2 < N.

T 1 max (µ T , µ I , µ V , α, γ, N) ≥ T 0 max (µ T , µ I , µ V , α, γ) (see
3.1. Parametrization by the logistic parameter. Since N > N 2 it is clear that the vertical line trough the point (N, 0) intersects the subdomain (P). We consider the family of infected equilibria X i (r), smoothly parameterized by the logistic parameter r ≥ r crit = r crit (N). The other critical values of r are r 1 (N) and r 2 (N). Below we simply denote them by r 1 and r 2 .

We rewrite the system (4) as

dX dt (t) = F (X(t), r), ( 22
)
where F is obviously a smooth function. In summary, the bifurcation parameter is r ∈ (r crit , +∞), the equilibria are X i (r), and the two critical values of the bifurcation parameter are r 1 and r 2 .

The Jacobian matrix [START_REF] Finzi | Latent infection of CD4+ T cells provides a mechanism for lifelong persistence of HIV-1, even in patients on effective combination therapy[END_REF] reads

J(X i (r)) =           r -µ T -2rT i T max -γV i (r) 0 -γT i γV i (r) -µ I γT i 0 Nµ I -µ V           , (23)
where

V i (r) = αN µ V - µ T γ + r γ ( 1 - µ V NγT max ) ,
and the associated characteristic equation ( 13) is

λ 3 + d 1 (r)λ 2 + d 2 (r)λ + d 3 (r) = 0, (24)
see [START_REF] Hartman | Ordinary Differential Equations[END_REF] for the expressions of the coefficients. We recall that d 3 > 0 in (I), hence d 3 (r) > 0 for r > r crit . Let us denote by λ i (r), i = 1, 2, 3, the three eigenvalues. The mappings r → λ i (r) are continuous, and smooth except at exceptional points where two of the λ i 's coincide as in Figure 5 (see [20, Chap. II]). At fixed r, one of the eigenvalues is real while the others are either real or complex conjugate. Since d 1 > 0, d 2 > 0 and d 3 > 0 in (I), we immediately observe:

Lemma 3.1. The polynomial [START_REF] Mallet-Paret | Snakes: Oriented families of periodic orbits, their sources, sinks, and continuation[END_REF] has no non-negative real root. At least one of the eigenvalues is negative whenever r > r crit .

We also have the notations, according to Subsection 13:

D 1 (r) = d 1 (r), D 2 (r) = d 1 (r)d 2 (r) -d 3 (r). ( 25
)
The relation between the Hurwitz determinant D 2 (r) and the sums of pairs of λ i 's is given by Orlando's formula (see, e.g., [13, Chap. XV]):

D 2 (r) = -(λ 1 (r) + λ 2 (r))(λ 2 (r) + λ 3 (r))(λ 1 (r) + λ 3 (r)). ( 26
)
From the results of Section 2 and Theorem 2.3, we know that the infected equilibrium X i (r) is asymptotically stable whenever r crit < r < r 1 or r > r 2 . Stability may be lost when either d 3 (r) = 0 or D 2 (r) = 0. The former case corresponds to r = r crit where, say, λ 3 = 0 and

λ 1,2 = 1 2 ( -d 1 ± √ d 2 1 -4d 2
) ,

(see Figure 5). The latter case occurs at r = r 1 where a pair of purely imaginary complex roots exists, as we are going to see in the next subsection. As a consequence of Lemma 3.1, only a pair of complex conjugate eigenvalues may cross the imaginary axis, which leads to a Hopf bifurcation. Increasing r, stability is regained when r > r 2 .

3.2.

Existence of two Hopf bifurcations. We are going to examine the critical values r = r 1 and r = r 2 and prove that they correspond to Hopf points. A detailed study of the Hopf bifurcations will be the subject of Section 4. Lemma 3.2. (i) At r = r k , k = 1, 2, the eigenvalues λ 1 (r k ) and λ 2 (r k ) are purely imaginary and conjugate:

λ 1 (r k ) = iω k , λ 2 (r k ) = -iω k , ω k = √ d 2 (r k ) = √ γNαµ I µ V + γαN + r k µ I µ V γT max N + r k µ 2 V γT max N ; ( 27
)
(ii) the third eigenvalue λ 3 (r k ) is real and negative:

λ 3 (r k ) = -d 1 (r k ) = - ( µ I + µ V + γαN µ V + r k µ V γT max N ) . (28)
Proof. At r = r k , D(r k ) = 0, therefore Orlando's formula (26) yields:

(λ 1 (r k ) + λ 2 (r k ))(λ 2 (r k ) + λ 3 (r k ))(λ 1 (r k ) + λ 3 (r k )) = 0. ( 29
)
As a consequence of Lemma 3.1, one of the eigenvalues, say λ 3 (r k ), is negative and λ 1 and λ 2 can not be real non-negative. Therefore (λ 2 (r k ) + λ 3 (r k ))(λ 1 (r k ) + λ 3 (r k )) can not vanish. It holds that λ 1 (r k )+λ 2 (r k ) = 0, hence λ 1 (r k ) and λ 2 (r k ) are purely imaginary and conjugate.

Let

λ 1 (r k ) = iω k , λ 2 (r k ) = -iω k , ω k > 0 (say). -d 1 (r k ) = λ 1 (r k ) + λ 2 (r k ) + λ 3 (r k ) = λ 3 (r k ), d 2 (r k ) = λ 1 (r k )λ 2 (r k ) + λ 2 (r k )λ 3 (r k ) + λ 3 (r k )λ 1 (r k ) = ω 2 k , we can infer that λ 3 (r k ) = -d 1 (r k ) and ω k = √ d 2 (r k ).
Theorem 3.3. Under the hypothesis [START_REF] Kuznetsov | Elements of applied bifurcation theory[END_REF], there is a Hopf bifurcation at the critical points r = r k , k = 1, 2 (the (X i (r k ), r k ) are called Hopf points).

Proof. We only need to prove that the derivative of r → λ 1 (r) does not vanish at r = r k (transversality condition). Let us consider the critical point r = r 1 . From Lemma 3.2 we know that λ 1,2 (r 1 ) = ±i ω 1 , ω 1 = √ d 2 (r 1 ). In a neighborhood of r 1 , the mapping r → λ 1 (r) is smooth. We compute, differentiating [START_REF] Mallet-Paret | Snakes: Oriented families of periodic orbits, their sources, sinks, and continuation[END_REF]:

dλ 1 dr (r) = - λ 2 1 (r)d ′ 1 (r) + λ 1 (r)d ′ 2 (r) + d ′ 3 (r) 3λ 2 1 (r) + 2d 1 (r)λ 1 (r) + d 2 (r)
.

At r = r 1 , it comes

dλ 1 dr (r 1 ) = - -d 2 (r 1 )d ′ 1 (r 1 ) + i √ d 2 (r 1 )d ′ 2 (r 1 ) + d ′ 3 (r 1 ) -2d 2 (r 1 ) + 2d 1 (r 1 )i √ d 2 (r 1 )
, hence:

d dr (Reλ 1 )(r 1 ) = -d 2 (r 1 )d ′ 1 (r 1 ) + d ′ 3 (r 1 ) -d 1 (r 1 )d ′ 2 (r 1 ) 2d 2 (r 1 ) + 2d 2 1 (r 1 ) = 1 γ 2 N 2 T 2 max (2d 2 + 2d 2 1 ) { µ I µ V γ 2 N 2 T 2 max -NγT max (2γNαµ V + 2γNαµ I + µ 3 V + 3µ 2 V µ I + µ V µ 2 I ) -2r 1 µ I µ 2 V -2r 1 µ 3 V } . ( 30
)
Under the hypothesis ( 21), T max is large, the other quantities including N being fixed. Then the sign of (d/dr)(Reλ 1 )(r 1 ) is clearly positive. On the other hand, Imλ 1 (r 1 ) = ω 1 = √ d 2 > 0, and Reλ 3 (r 1 ) = -Red 1 (r 1 ) < 0. According to the classical Hopf bifurcation theorem, there exists a Hopf bifurcation at the critical point r = r 1 .

A similar argument holds at the other critical point r = r 2 at which (d/dr)(Reλ 1 )(r 2 ) < 0.

3.3.

Instability. While the Routh-Hurwitz criterion ensures asymptotic stability, instability does not hold automatically when the criterion is not satisfied. Additional information is needed. We are now in position to complete the analysis. Theorem 3.4. Under the hypothesis (21), the infected equilibrium X i (r) is unstable whenever r 1 < r < r 2 .

Proof. Let N > N 2 and r ∈ (r 1 (N), r 2 (N)) be fixed. We have D 2 (r) < 0. According to Lemma 3.1, λ 3 (r) < 0, say. Hence, by the reciprocal of Ruth-Hurwitz criterion, one of the eigenvalues λ j (r), j = 1, 2 has a non-negative real part. Since no one can be real non-negative, it follows that λ 1 (r) and λ 2 (r) are complex conjugate. If Reλ 1 is positive, instability is proved. Therefore, assume Reλ 1 = 0 which implies λ 1 (r) + λ 2 (r) = 0. Then, by Orlando's formula [START_REF] Nowak | Virus dynamics[END_REF], D 2 (r) = 0 and we reach a contradiction. More precisely, we have that λ 3 (r) < 0 and λ 1 (r) and λ 2 (r) are complex conjugate with positive real part.

4. Nonlinear analysis. The framework is the same as in Section 3: N > 0 is fixed, T max is large as in [START_REF] Kuznetsov | Elements of applied bifurcation theory[END_REF]. Sharper hypotheses will be made below.

Our aim is to study the system [START_REF] Levy | Immunological and virological efficacy of a therapeutic immunization combined with interleukin-2 in chronically HIV-1 infected patients[END_REF] in some neighborhood of the critical points r 1 and r 2 , at which Hopf bifurcations occur.

4.1.

A local study near the critical values. For r in some neighborhood of r k , k = 1, 2, we subtract the infected equilibria and set X = X -X i (r), namely

X1 = T -T i , X2 = I -I i , X3 = V -V i .
The system [START_REF] Levy | Immunological and virological efficacy of a therapeutic immunization combined with interleukin-2 in chronically HIV-1 infected patients[END_REF] now reads for X:

d X dt = J(X i (r)) X + F ( X), (31)
where J(X i (r)) is given by [START_REF] Levy | Enhanced T cell recovery in HIV-1-infected adults through IL-7 treatment[END_REF] and

F ( X) = (- r T max X2 1 -γ X1 X3 , γ X1 X3 , 0).
Let us recall that the Jacobian matrix J(X i (r)) has the eigenvalues λ 1,2 (r k ) = ±iω k and λ 3 (r k ) = -d 1 (r k ), at the critical point r = r k . The corresponding eigenvectors are denoted by ν(r k ), ν(r k ) and u 3 (r k ). We set u 1 (k) = Re ν(r k ), u 2 (r k ) = Imν(r k ). It comes:

u 1 =               µ V s k N(s 2 k +ω 2 k ) µ V Nµ I 1               , u 2 =               µ V ω k N(s 2 k +ω 2 k ) ω k Nµ I 0               , u 3 =               µ V N(µ I +µ V ) -1 N -αγ µ V µ I -r k µ V µ I γN 2 T max 1              
,

where s k = - r k µ V γT max N - γNα µ V .
We define the matrix

P k = (u 1 (r k ), -u 2 (r k ), u 3 (r k )). As P k is clearly invertible, we set Q k = (P k ) -1 .
Omitting the subscript k, the coefficients of P k and Q k will be respectively denoted by P i, j and Q i, j . Under the linear transformation X = PY, the system (31) has the canonical form:

d dt           Y 1 Y 2 Y 3           =           0 -ω k 0 ω k 0 0 0 0 -d 1 (r k )                     Y 1 Y 2 Y 3           +           F 1 F 2 F 3           , (32)
where

F = Q k F (PY).
We also set

r = r k + µ k .
According to the Hopf bifurcation theorem (see [15, Thms. I and II]), it follows from Theorem 3.3 that, for k = 1, 2, the system (32) has, in some neighborhood of Y = 0, µ k = 0, a family of smooth periodic solutions (Y k # (t), µ k ) of period T k , smoothly parameterized by ε:

(Y k # (ε), T k (ε), µ k (ε)), ε ∈ (0, ε k 0 ), (33)
Obviously, ε k 0 depends on all the parameters, especially N and T max . The sign of µ k (ε) determines the direction of the bifurcation, in other words if it is supercritical or subcritical. Another important information is given by the Floquet exponents β k (ε): the periodic solution is orbitally asymptotically stable with asymptotic phase if β k (ε) < 0, unstable otherwise.

Since the nonlinearities are analytic, µ k , T k , β k are also analytic and have expansions, for k = 1, 2, 0 < ε < ε k 0 :

µ k (ε) = µ k 2 ε 2 + . . . , T k (ε) = 2π ω k (1 + τ k 2 ε 2 + . . .), β k (ε) = β k 2 ε 2 + . . . .
If the leading terms µ k 2 and β k 2 do not vanish, they provide us with the direction of bifurcation and the stability, respectively. 4.2. Computation of the first Lypunov coefficient. In order to analyze the bifurcation, we compute the first Lyapunov coefficient l 1 (r k ) of the restricted system on the center manifold at the critical parameter values r k . Loosely speaking, the latter two-dimensional system in the complex variable z reads (see [15, p. 59]),

dz dt = λ 1 (µ k )z + g(z, z; µ k ), (34)
where g(z, z; µ k ) is smooth and has a Taylor expansion with coefficients g i, j (µ k ). The main issue is to control the cubic term in the expansion of [START_REF] Wang | Mathematical analysis of the global dynamics of a model for HIV infection of CD4+ T cells[END_REF], namely the sign of the first Lyapunov coefficient (see [START_REF] Kuznetsov | Elements of applied bifurcation theory[END_REF]Sec. 3.4 and 5.4]), which should not vanish (nondegeneracy condition).

To compute the g i, j (0)'s, we follow the Recipe-Summary of [15, p. 88-90]. At Y = 0, µ k = 0:

g 11 = 1 4       ∂ 2 F 1 ∂Y 2 1 + ∂ 2 F 1 ∂Y 2 2 + i       ∂ 2 F 2 ∂Y 2 1 + ∂ 2 F 2 ∂Y 2 2             = - r k 2T max Q 11 (P 2 11 + P 2 12 ) + 1 2 γP 11 (Q 12 -Q 11 ) + i { - r k 2T max Q 21 (P 2 11 + P 2 12 ) + 1 2 γP 11 (Q 22 -Q 21 )
} ,

g 02 = 1 4       ∂ 2 F 1 ∂Y 2 1 - ∂ 2 F 1 ∂Y 2 2 -2 ∂ 2 F 2 ∂Y 1 ∂Y 2 + i       ∂ 2 F 2 ∂Y 2 1 - ∂ 2 F 2 ∂Y 2 2 + 2 ∂ 2 F 1 ∂Y 1 ∂Y 2             = 1 2 { r k T max (Q 11 (P
+ γP 11 (Q 22 -Q 21 ) + γP 12 (Q 12 -Q 11 )] } , g 20 = 1 4       ∂ 2 F 1 ∂Y 2 1 - ∂ 2 F 1 ∂Y 2 2 + 2 ∂ 2 F 2 ∂Y 1 ∂Y 2 + i       ∂ 2 F 2 ∂Y 2 1 - ∂ 2 F 2 ∂Y 2 2 -2 ∂ 2 F 1 ∂Y 1 ∂Y 2             = 1 2 { r k T max (Q 11 (P
+ γP 11 (Q 22 -Q 21 ) -γP 12 (Q 12 -Q 11 )] } , G 21 = 1 8 ( ∂ 3 F 1 ∂Y 3 1 + ∂ 3 F 1 ∂Y 1 ∂Y 2 2 + ∂ 3 F 2 ∂Y 2 1 ∂Y 2 + ∂ 3 F 2 ∂Y 3 2 + i       ∂ 3 F 2 ∂Y 3 1 + ∂ 3 F 2 ∂Y 1 ∂Y 2 2 - ∂ 3 F 1 ∂Y 2 1 ∂Y 2 - ∂ 3 F 1 ∂Y 3 2       ) = 0.
Since the system (32) is three-dimensional, additional computations are needed:

h 11 = 1 4       ∂ 2 F 3 ∂Y 2 1 + ∂ 2 F 3 ∂Y 2 2       = 1 2 [ - r k T max Q 31 (P 2 11 + P 2 12 ) + γP 11 (Q 32 -Q 31 )
] , } .

h 20 = 1 4       ∂ 2 F 3 ∂Y 2 1 - ∂ 2 F 3 ∂Y 2 2 -2i ∂ 2 F 3 ∂Y 1 ∂Y 2       = 1 2 { r k T max Q 31 (P
Solving the linear system

-d 1 w 11 = -h 11 , (-d 1 -2iω 1 )w 20 = -h 20 ,
we get:

w 11 = h 11 d 1 , w 20 = h 20 d 1 + 2i √ d 2 .
Furthermore we compute: ]} , ]} .

G 110 = 1 2 ( ∂ 2 F 1 ∂Y 1 ∂Y 3 + ∂ 2 F 2 ∂Y 2 ∂Y 3 + i ( ∂ 2 F 2 ∂Y 1 ∂Y 3 - ∂ 2 F 1 ∂Y 2 ∂Y 3 )) = 1 2 { -
G 101 = 1 2 ( ∂ 2 F 1 ∂Y 1 ∂Y 3 - ∂ 2 F 2 ∂Y 2 ∂Y 3 + i ( ∂ 2 F 2 ∂Y 1 ∂Y 3 + ∂ 2 F 1 ∂Y 2 ∂Y 3 )) = 1
Eventually we get:

g 21 = G 21 + 2G 110 w 11 + G 101 w 20 .
Let us denote (dλ 1 /dr)(r k ) by λ ′ 1 (r k ) which do not vanish, see Thm. Let N be fixed, 0 < N < N * (hence H < 0). Then, for k = 1, 2, for all T max ≥ T 2 max (given by Lemma 4.4), there exists an ε k 0 > 0 such that the system (32) has a smooth family of nontrivial periodic solutions (Y k # (ε), µ k (ε)) with period T k (ε), parameterized by ε ∈ (0, ε k 0 ). Moreover, (i) as ε → 0, T k (ε) → 2π/ω k (see [START_REF] Pantaleo | New concepts in the immunopathogenesis of human immunodeficiency virus infection[END_REF]); (ii) the bifurcation at µ 1 = 0 is supercritical and the periodic solution Y 1 # (ε) is asymptotically stable with asymptotic phase; (iii) the bifurcation at µ 2 = 0 is subcritical and the periodic solution Y 2 # (ε) is orbitally asymptotically stable with asymptotic phase.

It follows from Theorem 4.5:

Corollary 1. Let N be fixed, 0 < N < N * . Then, for k = 1, 2, for all T max ≥ T 2 max , there exists an ε k 0 > 0 such that the system (4) has a smooth family of nonconstant periodic solutions (X k # (ε), r k (ε)) associated with the Hopf points (X i (r k ), r k ), with period T k (ε), and parameterized by ε ∈ (0, ε k 0 ):

X k # (ε) = X i (r k ) + P k Y k # (ε), r k (ε) = r k + µ k (ε). ( 42
)
Moreover, (i) as ε → 0, T k (ε) → 2π/ω k (see [START_REF] Pantaleo | New concepts in the immunopathogenesis of human immunodeficiency virus infection[END_REF]);

(ii) the bifurcation at r = r 1 is supercritical and the periodic solution X 1 # (ε) is asymptotically stable with asymptotic phase; (iii) the bifurcation at r = r 2 is subcritical and the periodic solution X 2 # (ε) is orbitally asymptotically stable with asymptotic phase.

According to [24, p. 426-427], the number 2π/ω k is called the "virtual period" of the Hopf point (X i (r k ), r k ). 4.4. "Snake" of orbits. Hereafter we call a nonconstant periodic solution (X # (r), r) of ( 1)-(3) an orbit in the (X, r)-space. We have already proved the existence of two families of such orbits (X k # (ε), r k (ε)), locally near the Hopf point (X i (r k ), r k ). The next issue is to investigate whether there exists a "bridge" or "global connection" between the two Hopf points, by a family of orbits. The answer is positive thanks to the theory of "snakes" developed by Mallet-Paret and York [START_REF] Mallet-Paret | Snakes: Oriented families of periodic orbits, their sources, sinks, and continuation[END_REF] (see also [11, p. 27-34]).

According to [START_REF] Mallet-Paret | Snakes: Oriented families of periodic orbits, their sources, sinks, and continuation[END_REF], the Hopf points (X i (r k ), r k ) are isolated centers. At such an isolated center, the crossing number χ is the net number of pairs of eigenvalues crossing the imaginary axis (see [24, p. 427]). Let E(r) denote the sum of the multiplicities of the eigenvalues of the Jacobian matrix J(X i (r)) having strictly positive real part, then

χ(r k ) = 1 2 (E(r k +) -E(r k -)).
Clearly χ(r k ) is equal to the sign of λ ′ 1 (r k ), see Theorem 3.3, hence χ(r 1 ) = +1, χ(r 2 ) = -1. Next, the center index is denoted by the Chinese character "中" (zhong):

中(X i (r k ), r k ) = χ(r k ) (-1) E(r k ) .
An isolated center is called a source if 中= +1, a sink is 中= -1 (see [24, p. 428]). Here we have: Lemma 4.6. 中(X i (r 1 ), r 1 ) = +1, 中(X i (r 2 ), r 2 ) = -1, hence the Hopf points (X i (r 1 ), r 1 ) and (X i (r 2 ), r 2 ) are respectively a source and a sink.

In the terminology of [24, p. 436], a "snake" is defined in the following way: let γ(ε) be a path (see [24, p. 432]) with domain J and let Orb(γ) be the set of orbits (X k # (ε), r(ε)) for which there is an ε ∈ J such that γ(ε) is on that orbit. Then, a set of orbits S is called a "snake" if there is maximal oriented path γ such that Orb(γ) = S . Theorem 4.7. There is a unique "snake" S 1 with source (X 1 i , r 1 ) and sink (X 2 i , r 2 ). Proof. According to the terminology of [START_REF] Mallet-Paret | Snakes: Oriented families of periodic orbits, their sources, sinks, and continuation[END_REF], the Hopf theorem yields that there is a unique "snake" emanating from a Hopf point. Let us consider the snake S 1 emanating from the source (X 1 i , r 1 ), and let us apply the "Snake Termination Principle" (see [24, p. 438-439]). Since S 1 has at least one end, it is an open snake and its sink is either ∞ or a center. From Lemma B.1(see App. B), it can not be ∞. The only center which is not (X 1 i , r 1 ) is (X 2 i , r 2 ), which is a sink. Finally, (X 1 i , r 1 ) is called the source of S 1 , (X 2 i , r 2 ) its sink. snake Figure 8. Numerical representation of the snake between the source (X i (r 1 ), r 1 ) and the sink (X i (r 2 ), r 2 ). The dotted line represents T , the solid line V and the dashed line 10 I. The variable on the horizontal axis is rr 1 , the constant infected equilibrium being normalized to 0. Remark 2. A further issue is to know whether the elements of S 1 are orbitally asymptotically stable for all values of r, r 1 < r < r 2 . De Lenheer and Smith [START_REF] Leenheer | Virus dynamics: a global analysis[END_REF] observed that the system (1)-( 3) gives rise to a three-dimensional competitive dynamical system, in such a case a generalized Poincaré-Bendixson theorem holds in 3D. It follows from [6, Lemma 3.6] and [35, Thm. 1.2] that, for each r ∈ (r 1 , r 2 ), there is at most a finite number of orbits and at least one is orbitally stable. Unfortunately, we can not infer that the latter belongs to S 1 , although there is a numerical evidence.

5. Numerical results. In view of the numerical computations, we take N = 300 to fix ideas. Other parameters are as in Table 1. Initial data are T 0 = T u , I 0 = 0.0, V 0 = 0.0185.

We start from the value r crit = 0.05625; then, we increase the logistical parameter r monotonically, until a critical condition is reached at which any further change would result in instability.

We present the graphs of numerical solution of the system (1)-( 3) and the trajectory of the solutions in the T -V-I space, see Figures 9101112. In Figures 10 and11 corresponding to subdomain (P), the solution approaches periodic orbits. 6. Conclusion. In summary, the infected equilibria are asymptotically stable in a region (I), except in an explicit subdomain (P) which may be biologically interpreted as a perturbation of the infection by a specific or unspecific immune response against HIV. The existence of Hopf bifurcations and periodic orbits in (P) lead to some very interesting mathematical issues which, as far as we know, have not been fully addressed in previous articles on the logistic HIV model. In a forthcoming paper [START_REF] Fan | Two dimensional stability in a HIV model with logistic growth term[END_REF], we will take the viral diffusion into account (see [START_REF] Brauner | Heterogeneous viral environment in a HIV spatial model[END_REF]).

Studying the subdomain (P) in models taking antiretroviral treatment into account and/or studying the impact of immunotherapy may allow generating new hypotheses on HIV elimination.

L(t) ≤ M 1 m + (L(0) - M 1 m )e -mt ,
therefore L(t) is bounded, and so is I(t) = L(t) -T (t). Finally, it follows from (3) that V(t) is also bounded.

  and

2 XINYUEFigure 1 .

 21 Figure 1. Typical course of HIV infection. Patterns of CD4 + T-cell decline and virus load increase vary greatly from one patient to another, as do the actual values of viral RNA load. (Modified from Pantaleo et al. [27]).

Table 1 .

 1 Variables and Parameters Parameters and variables Values Dependent variables T Uninfected CD4 + T-cell population mm -3 I Infected CD4 + T-cell density mm -3 V HIV population size mm -3 Parameters and constants r Proliferation rate of the CD4 + T-cell population 0.2 day -1 N Number of virus produced by infected cells 1000 α Production rate for uninfected CD4 + T cells 1.5 day -1 mm -3 γ Infection rate of uninfected CD4 + T cells 0.001 mm 3 day -1 T max Maximal population level of CD4 + T cells at 1500 mm -3 which the CD4 + T-cell proliferation shuts off µ T Death rate of uninfected CD4 + T-cell population 0.1 day -1 µ I Death rate of infected CD4 + T-cell population 0.5 day -1 µ V Clearance rate of free virus 10 day -1 Derived variables T 0 CD4 + T-cell population for HIV negative persons mm -3

4 XINYUE 3 Figure 2 .

 432 Figure 2. Typical profiles of virus vs. target cells in the phase plane (front), virus (left) and target cells (right) in the case of stable infection.Here N = 1000, r = 0.2. The other quantities are as in Table1.

10 XINYUE 3 .

 103 FAN, CLAUDE-MICHEL BRAUNER AND LINDA WITTKOPCase It corresponds to N = N 1 or N = N 2 . In such a case, ∆ = 0 and the polynomial Ar 2 + Br + C has the double root r = -B/2A. However, this solution makes sense only if B < 0, namely N > N 0 . Then by Lemma 2.2, we discard the case N = N 1 . Therefore at N = N 2 , D 2 is positive except at r(N 2 ) = -B/2A where it vanishes.

Figure 4 .

 4 Figure 4. Full portrait of the domains (I), (U) and (P) in the (N, r) parameter space, in logarithmic scales. Numerical values are as in Table 1. Here N 0 = 154.35, N 2 = 169.81, N crit (0) = 666.67, r crit (N 2 ) = 0.08, r 1 (N 2 ) = r 2 (N 2 ) = 14.03.

2 Figure 5 .

 25 Figure 5. Real part (left) and imaginary part (right) of eigenvalues λ 1 and λ 2 as functions of the bifurcation parameter r at N = 300 (to fix ideas). The graphic starts at r crit = 0.05625, λ 1 = λ 2 at r = 0.05734. The real eigenvalue λ 3 < -10 is not pictured.

2

 2

2 Figure 7 .

 27 Figure 7. For larger r, real part (left) and imaginary part (right) of eigenvalues λ 1 and λ 2 as functions of the bifurcation parameter r at fixed N = 300. They cross the imaginary axis at r 2 = 464.1225.

  3.3. According to [15, p. 90][21, p. 98-99 and 178], it holds: Theorem 4.5.

22 XINYUEFigure 9 .mm 3 Figure 10 .

 229310 Figure 9. Numerical solution of the system (1)-(3) the domain (I). Parameter values are N = 300 and r crit = 0.05625 < r = 1.0 < r 1 = 2.1846. The infected equilibrium is stable.

3 Figure 11 .mm 3 Figure 12 .

 311312 Figure 11. Numerical solution of the system (1)-(3) in (P). Parameter values are N = 300 and r = 200.0 < r 2 = 464.1225. The infected equilibrium is unstable, a periodic orbit is observed.

Figure 13 . 2 -

 132 Figure[START_REF] Gantmakher | The theory of matrices[END_REF]. Comparison of the "perturbated" subdomains (P). The closed curve is the "isola" in the "full" logistic model.

  2 12 -P 2 11 ) + 2Q 21 P 11 P 12 ) + γP 11 (Q 12 -Q 11 )

	-γP 12 (Q 22 -Q 21 ) + i[	r k T max	(Q 21 (P 2 12 -P 2 11 ) -2Q 11 P 11 P 12 )

  2 12 -P 2 11 ) -2Q 21 P 11 P 12 ) + γP 11 (Q 12 -Q 11 )

	+ γP 12 (Q 22 -Q 21 ) + i[	r k T max	(Q 21 (P 2 12 -P 2 11 ) + 2Q 11 P 11 P 12 )

  31 P 11 P 12 -γP 12 (Q 32 -Q 31 )]

			2 12 -P 2 11 ) + γP 11 (Q 32 -Q 31 )
	+ i[	2r T max	Q

  2r k T max P 13 (Q 11 P 11 + Q 21 P 12 ) + γ(Q 12 -Q 11 )(P 11 + P 13 ) Q 11 P 12 -Q 21 P 11 ) + γ(Q 22 -Q 21 )(P 11 + P 13 ) -γP 12 (Q 12 -Q 11 )

	+ γP 12 (Q 22 -Q 21 ) + i	[	2r T max	P 13 (

  Q 21 P 12 -Q 11 P 11 ) + γ(Q 12 -Q 11 )(P 11 + P 13 ) 22 -Q 21 )(P 11 + P 13 ) + γP 12 (Q 12 -Q 11 )

	2 P 13 (-γP 12 (Q 22 -Q 21 ) + i { 2r k T max	[	2r T max	P 13 (Q 11 P 12 + Q 21 P 11 )
	+ γ(Q

Acknowledgments. We gratefully thank Odo Diekmann and Luca Lorenzi for their comments.

then the first Lyapunov coefficient at r k reads:

Re (ig 20 g 11 + ω k g 21 ) ,

and

An explicit computation of the l 1 (r k )'s is uneasy. We are going to evaluate their sign via an asymptotic method, as well as the signs of µ k 2 and β k 2 .

Asymptotic analysis.

In this part we will take advantage of the limit T max → +∞.

Passing formally to the limit in ( 1)-( 3), we get to the system

(2)(3) for (T ∞ , I ∞ , V ∞ ), with bifurcation parameter r -µ T . The latter is nothing but the net T-cell proliferation rate, which changes sign at r = µ T .

The system (38)(2)(3) is interesting didactically since quite simpler than the original system (1)-( 3). However, it has clearly no biological meaning since it would implicate a dysfunction of the homeostatic system. With the same notations and arguments as above, it is not difficult to compute the equilibria X ∞ u and X ∞ i and check that

Beside, a Hopf bifurcation occurs at the curve

The infected equilibrium X ∞ i is unstable for r > r ∞ 1 (N), the latter region corresponding to (P ∞ ).

Taking the limit T max → +∞ in [START_REF] Jones | Transient viremia, plasma viral load, and reservoir replenishment in HIVinfected patients on antiretroviral therapy[END_REF], elementary calculations provide the asymptotic behaviors of r k and ω k :

Then it follows immediately from (30) that:

A natural idea is to plug these asymptotic formulas in Subsection 4.2. Since the computations are technical but straightforward, we are going to skip them.

In Lemma 4.3, D and H are respectively given by

). Whereas D is always positive, the sign of H varies upon N.

Since H(0) < 0, lim N→+∞ H(N) = +∞, let us define N * as the (first) positive zero of H. Therefore H < 0 for 0 ≤ N < N * . With the values of Table 1, N * ≃ 3.9 × 10 5 which is far beyond the biological size of N, therefore this hypothesis is reasonable.

It easily follows from (37) and Lemmas 4.2-4.3 that

Clearly, if T max is large enough, the signs of all the above quantities are given by their respective leading terms. More precisely, it holds: Lemma 4.4. Let N be fixed, 0 < N < N * (hence H < 0). There exists a T 2 max ≥ T 1 max (see [START_REF] Kuznetsov | Elements of applied bifurcation theory[END_REF]), which depends of µ T , µ I , µ V , α, γ and N, such that, whenever T max ≥ T 2 max , we have:

We summarize our results in the following Theorem and Corollary which follow from Lemma 4.4 and [15, Thms. I and II]):

Appendix A. The "full" logistic model. It is possible to replace (1) by

(see [START_REF] Culshaw | A delay-differential equation model of HIV infection of CD4 + T-cells[END_REF][START_REF] Perelson | Mathematical analysis of HIV-I: dynamics in vivo[END_REF]). However, the proportion of infected cells is very small, on the order of 10 -4 to 10 -5 of T cells, therefore this correction in the logistic proliferation term may be ignored (see [31, p. 10]). Wang and Li [START_REF] Wang | Mathematical analysis of the global dynamics of a model for HIV infection of CD4+ T cells[END_REF] studied the global dynamics of the system (43)(2)(3). Again, the boundary of (P) is the locus of points (N, r) such that D 2 (N, r) = 0, where

) .

Numerical evidence (see Fig. 13) shows that the subdomain (P) is bounded, in contrast to the situation of the model with [START_REF] Antia | The role of models in understanding CD8+ T-cell memory[END_REF]. Moreover, the shape of the boundary of (P) is reminiscent of an "isola" in bifurcation theory (see [START_REF] Dellwo | On the birth of isolas[END_REF]). Finally, as other authors, we omitted the term -γVT in (3), but its omission is as an approximation (see e.g., [START_REF] Perelson | Dynamics of HIV infection of CD4+ T cells[END_REF]).

Appendix B. Boundness lemma.

Lemma B.1. The solution of the system (1)-( 3) remains bounded for all times t > 0.

Proof. Let (T (t), I(t), V(t)) the solution of the system (1)-( 3) with non-negative initial data (T (0), I(0), V(0)). It is not difficult to check that the non-negative octant R 3 is positively invariant with respect to the system.

From (1), we observe that dT dt ≤ αµ T T + rT (1 -T T max

).

Let ∆T (t) = T (t) -