This is a preprint of: “Global dynamics of a virus model with invariant algebraic surfaces”, Fabio Scalco Dias, Jaume Llibre, Clàudia Valls, Rend. Circ. Mat. Palermo (2), 2019. DOI: [10.1007/s12215-019-00417-0] GLOBAL DYNAMICS OF A VIRUS MODEL WITH INVARIANT ALGEBRAIC SURFACES FABIO SCALCO DIAS1 , JAUME LLIBRE2 AND CLAUDIA VALLS3 Abstract. In this paper by using the Poincaré compactification in R3 we make a global analysis for the virus system ẋ = λ − dx − βxz, ẏ = −ay + βxz ż = ky − µz with (x, y, z) ∈ R3 , β > 0, λ, a, d, k and µ are nonnegative pa- rameters due to their biological meaning. We give the complete description of its dynamics on the sphere at infinity. For two sets of the parameter values the system has invariant algebraic surfaces. For these two sets we provide the global phase portraits of the virus system in the Poincaré ball (i.e. in the compactification of R3 with the sphere S2 of the infinity). 1. Introduction and statement of the main results In this work we do a global analysis of the virus system ẋ = P (x, y, z) = λ − dx − βxz, (1) ẏ = Q(x, y, z) = −ay + βxz, ż = R(x, y, z) = ky − µz, where the state variables are (x, y, z) ∈ R3 , β > 0, λ, a, d, k and µ are nonnegative parameters due to their biological meaning. As usual, the dots denote derivative with respect to the time t. Note that β cannot be zero otherwise the differential system becomes linear, and consequently its study is trivial. System (1) was proposed in the book [7]. For more studies related with virus models, including the importance of our model (1), we refer the readers to [5, 8]. Let R[x, y, z] be the ring of real polynomials in the variables the variables x, y and z. We say that F = F (x, y, z) ∈ R[x, y, z] is a Dar- boux polynomial of system (1) if it satisfies (∇F ) · (P, Q, R) = KF , where K = K(x, y, z) is a real polynomial of degree at most 1, called the cofactor of F (x, y, z). If the cofactor is zero, then F (x, y, z) is a 2010 Mathematics Subject Classification. 34A34, 34C05, 34C14 . Key words and phrases. Invariant algebraic surfaces, Poincaré compactification, phase portrait, dynamics at infinity, virus model. 1 2 F.S. DIAS, J. LLIBRE AND C. VALLS polynomial first integral of system (1). If F (x, y, z) is a Darboux poly- nomial, then the algebraic surface F (x, y, z) = 0 is called an invariant algebraic surface; i.e. if an orbit of system (1) has a point on this surface, then the whole orbit is contained in it. We say that a C 1 function I(x, y, z, t) is an invariant of the differen- tial system (1) if dI/dt = 0 on the trajectories of the system. When an invariant function is independent of the time t, then it is a first integral. When a system has a Darboux polynomial F with a constant factor k = k0 ∈ R, then the function I(x, y, z, t) = F (x, y, z)e−k0 t is called a Darboux invariant of that system, see Chapter 8 of [3]. We say that the Darboux polynomials F1 , . . . , Fn generate all the Darboux polynomials of system (1) if any Darboux polynomial F of system (1) can be written as F = F1α1 · · · Fnαn with αk ∈ N for k = 1, . . . , n. The following theorem proved in [9] summarizes the results on the existence of invariant algebraic surfaces for system (1). Theorem 1. When βk(λ2 + d2 ) 6= 0, a set of generators for the set of all Darboux polynomials of system (1) are given in Table 1. Darboux Polynomial Cofactor Parameters F0 = x, K = −d − βz, λ = 0, F1 = x + y − λ/d, K = −d, a = d 6= 0, F2 = x + y + (a − d)z/k − λ/d, K = −d, µ = d 6= 0. Table 1. Invariant algebraic surfaces of system (1). As any polynomial differential system the virus system (1) can be extended to an analytic system on a closed ball of radius one, whose interior is diffeomorphic to R3 and its boundary, the 2–dimensional sphere S2 , plays the role of the infinity because in R3 one can go or come from infinity in as many directions as points have S2 . This closed ball is denoted by B and called the Poincaré ball. The technique for doing such an extension is precisely the Poincaré compactification for a polynomial differential system in R3 , which is described in details in [2], see also [1]. By using this compactification technique the dynamics of system (1) at infinity can be determined. Theorem 2. The phase portrait of the Poincaré compactification of system (1) at the infinity S2 is topologically equivalent to the one de- scribed in Figure 1. 3 q′ z q y x Figure 1. Phase portrait of the virus system (1) at the infinity S2 of the Poincaré ball. The circles of S2 con- tained in x = 0 and z = 0 are filled of singular points. Additionally there are two diammetrally opposite singu- lar points on z = 0 of nodal type. We go further and describe the global dynamics of system (1) for fixed values of the parameters having the invariant planes F1 = 0 or F2 = 0 described in Table 1. Theorem 3. The following statements hold for system (1) when βk(λ2 + d2 ) 6= 0. (a) For d = a 6= 0 system (1) has the Darboux polynomial F1 (x, y, z) = x + y − λ/d and the Darboux invariant Id (x, y, z, t) = (x + y − λ/d)edt . The singular points on F1 (x, y, z) = 0 are described in Table 2. The phase portraits of system (1) restricted to the compactified invariant plane F1 (x, y, z) = 0 are described in Figure 7. (b) For d = µ 6= 0 system (1) has the Darboux polynomial (a − d) λ F2 (x, y, z) = x + y + z− k d and the Darboux invariant   (a − d) λ dt Id (x, y, z, t) = x + y + z− e . k d The singular points on F2 (x, y, z) = 0 are described in Table 2. The phase portraits of system (1) restricted to the compactified invariant plane F2 (x, y, z) = 0 are described in Figure 7. 4 F.S. DIAS, J. LLIBRE AND C. VALLS We remark that system (1) restricted to the invariant plane F0 = x = 0 with λ = 0 is linear. We say that a set V ⊂ B is invariant by the flow of system (1) if for any p ∈ V the whole orbit passing through p is contained in V . The sphere of the infinity always is an invariant set. Let ϕ(t) = ϕ(t, p) be the solution of the compactified system (1) passing through the point p ∈ B when t = 0, defined on its maximal interval Ip = R (because B is compact). Then the α–limit set of ϕ is the invariant set α(ϕ) = {q ∈ B : ∃ {tn } such that tn → −∞ and ϕ(tn ) → q as n → ∞}. In a similar way the ω–limit set of ϕ is the invariant set ω(ϕ) = {q ∈ B : ∃ {tn } such that tn → ∞ and ϕ(tn ) → q as n → ∞}. The following theorem gives a complete description of the α– and ω– limit sets in the Poincaré ball, including its boundary S2 , for system (1) having the invariant algebraic surfaces of Table 1. Theorem 4. Let γ be an orbit of the compactified system (1) with βk 6= 0. Assume a = d 6= 0, or µ = d 6= 0. If γ is contained in the interior of B and outside the invariant plane F1 = 0 if a = d 6= 0, or F2 = 0 if µ = d 6= 0, then the following two statements hold. (a) The α–limit set of γ is some of the four infinite singular points q, q ′ , E2 , E2′ contained in the boundary of the plane Fi = 0 for i = 1, 2. (b) The ω–limit set of γ is contained in the closure of the invariant plane Fi = 0 in the Poincaré ball. More precisely, If µ = λ = 0, or λ = a = 0, then ω(γ) ⊂ {(0, 0, z), z ∈ R}, i.e. the straight line x = 0 of Figure 7. If µ = 0 and λ 6= 0, or a = 0 and λ 6= 0, then ω(γ) ⊂ {P0 , E2 , E2′ }. If µ 6= 0 and βkλ = d2 µ, or a 6= 0 and βkλ = ad2 , then ω(γ) = {P0 }. If µ 6= 0 and βkλ 6= d2 µ, or a 6= 0 and βkλ 6= ad2 , then ω(γ) ⊂ {P0 , P1 }. Here P0 , P1 are singular points of system (1), contained in the plane F1 = 0 (see Table 2), E2 is the origin of the local chart U2 , q = (−1, 0, 0) is a singular point in the local chart U1 , and E2′ and q ′ are the diametrally opposite points to E2 and q in the Poincaré sphere. 5 The paper is organized as follows. In section 2 we prove Theorem 2 by using the Poincaré compactification for a polynomial vector field in R3 . In section 3 we prove Theorem 3 by using the Poincaré compactification for a polynomial vector field in R2 . For precise definitions of all these notions see Chapter 5 of [3]. Finally the proof of Theorem 4 is given in section 4. 2. Proof of Theorem 2 In this section we analyze the flow of system (1) at infinity. Consider system (1) or equivalently its associated polynomial vector field X = (λ − dx − βxz, −ay + βxz, ky − µz). In the next three subsections we will study the Poincaré compactifi- cation, p(X), of system (1) in the local charts Ui and Vi , i = 1, 2, 3 in order to understand the global behavior of the solutions near infinity. See [1] for more details on these charts. 2.1. In the local chart U3 . Using the results of [1] we have that the expression of the Poincaré compactification Z3 = p(X) of system (1) in the local chart U3 is given by z˙1 = −βz1 + z3 ((µ − d)z1 + λz3 − kz1 z2 ), (2) z˙2 = βz1 + z3 ((µ − a)z2 − kz22 ), z˙3 = −z32 (kz2 − µ). In the points of the sphere S2 that correspond to the points at infinity we have z3 = 0, and so system (2) becomes z˙1 = −βz1 , z˙2 = βz1 , z˙3 = 0. After eliminating the common factor z1 (by rescaling of the time) we do not obtain any singular point in the local chart U3 which is dis- tinguished. From this system we see that system (2) has at infinity a continuum of singular points given by (0, z2 , 0). Moreover DZ3 (0, z2 , 0) has, for each z2 , the eigenvalues −β and 0 with multiplicity 2. By the normal hyperbolicity theory (see [4] for details) to each point (0, z2 , 0) it arrives two orbits, see Figure 2. Note that the equilibria (0, z2 , 0) with z2 ∈ R represent half of a great circle through the north pole of S2 . The orbits of the system in the local chart U3 at infinity have the phase portrait given in Figures 2 and 1. 6 F.S. DIAS, J. LLIBRE AND C. VALLS z2 z1 Figure 2. The phase portrait in the local chart U3 at infinity. 2.2. In the local chart U1 . In the same way the Poincaré compacti- fication Z1 = p(X) of system (1) in the local chart U1 is z˙1 = βz2 (z1 + 1) + z3 ((d − a)z1 − λz1 z3 ), (3) z˙2 = βz22 + z3 ((d − µ)z2 + kz1 − λz2 z3 ), z˙3 = z3 (dz3 + βz2 − λz32 ). System (3) restricted to infinity (that is with z3 = 0 ) becomes z˙1 = βz2 (z1 + 1), z˙2 = βz22 , z˙3 = 0. So there is a line of singular points which is z2 = 0. These singular points are such that DZ1 (z1 , 0, 0) has, for each z1 the eigenvalue 0 with multiplicity 3. After eliminating the common factor z2 (by a rescaling of the time) we get that the unique singular point is q = (−1, 0, 0) which is an unstable node. Note that the equilibria (z1 , 0, 0) for z1 ∈ R represent in the Poincaré sphere S2 half of the equator. The orbits of the system in the local chart U1 at infinity have the phase portrait given in Figures 3 and 1. 2.3. In the local chart U2 . The expression of the Poincaré compact- ification Z2 = p(X) in the local chart U2 is z˙1 = −βz1 z2 (z1 + 1) + z3 ((a − d)z1 + λz3 ), (4) z˙2 = −βz1 z22 + z3 (k + (a − µ)z2 ), z˙3 = z3 (az3 − βz1 z2 ). Now the unique point of the local chart U2 which is not covered by the local charts U1 , V1 , U3 and V3 is the origin of coordinates of U2 . System (4) restricted at the infinity (that is with z3 = 0 ) becomes z˙1 = −βz1 z2 (z1 + 1), z˙2 = −βz1 z22 , z˙3 = 0. 7 z2 z1 q Figure 3. The phase portrait in the local chart U1 at infinity. There are two lines of singular points z1 = 0 and z2 = 0. After removing these two lines of singular points, there is a singular point q = (−1, 0, 0) which is a stable node and its local phase portrait is the one given in Figure 4. z2 z1 q Figure 4. Local phase portrait at the origin of U2 at infinity. We observe that the flows in the Vi charts for i = 1, 2, 3 are the same as the ones in the respective Ui charts for i = 1, 2, 3 but with the time reversed because the compactified vector field p(X) in Vi coincides with the vector field in Ui multiplied by −1 for each i = 1, 2, 3. See [1]. From subsections 2.1, 2.2 and 2.3 it follows the proof of Theorem 2. 3. Proof of Theorem 3 In Theorem 1 it was proved that the virus system (1) has the invari- ant algebraic surfaces F1 (x, y, z) = x + y − λ/d = 0 when a = d 6= 0 8 F.S. DIAS, J. LLIBRE AND C. VALLS and F2 = x + y + (a − d)z/k − λ/d = 0 when µ = d 6= 0. It is easy to see that Id (x, y, z, t) = Fi (x, y, z)edt satisfies   dId ∂Id ∂Id ∂Id ∂Id dt dFi = + ẋ + ẏ + ż = e Fi (x, y, z)d + = 0, dt ∂t ∂x ∂y ∂z dt where i = {1, 2}. Therefore Id is a Darboux invariant of system (1). This completes the proof of the first part of statements (a) and (b) of Theorem 3. System (1) has two isolated singular points in R3 , namely     λ aµ λ dµ λk d (5) p0 = , 0, 0 , p1 = , − , − . d βk a βk aµ β Note that these two singular points are on the two invariant algebraic surfaces F1 (x, y, z) = 0 when a = d 6= 0 and F2 (x, y, z) = 0 when µ = d 6= 0. In what follows we shall analyze the flow of the system in these planes and relate their dynamics with the dynamics at infinity. 3.1. Phase portrait on {F1 (x, y, z) = 0}. We can consider in this subsection system (1) with a = d 6= 0 and βkd(λ2 + d2 ) 6= 0. System (1) restricted to the invariant plane F1 (x, y, z) = 0 is given by the differential system kλ (6) ẋ = λ − dx − βxz, ż = −kx − µz + . d Doing the change of variables (with a rescaling of time) d2 λ d (7) x= X+ , z= Z, t = dT, βk d β system (6) becomes (8) ẋ = −x − a1 z − xz, ż = −x − a2 z 3 where a1 = βkλ/d , a2 = µ/d ∈ [0, ∞), and we have renamed the new variables (X, Z) by (x, z), and the dot denotes derivative with respect to the new time T . We now study the infinite singular points of system (8). We will use the Poincaré compactification for a polynomial vector field in R2 which is described in Chapter 5 of [3]. 3.1.1. Infinite singular points. System (8) in the local chart U1 becomes ż1 = −z2 + z1 (z1 + (1 − a2 )z2 + a1 z1 z2 ), ż2 = z2 (z1 + z2 + a1 z1 z2 ). When z2 = 0 the only infinite singular point is the origin. Since the origin is nilpotent, using Theorem 3.5 of [3] we obtain that the phase 9 portrait of the origin consists of one hyperbolic, one elliptic and two parabolic sectors, see Figure 5. z1 Figure 5. Local phase portrait at the origin of U1 at infinity. Now it only remains to study if the point (0, 0) of the local chart U2 is a singular point. System (8) in the local chart U2 writes (9) ż1 = −z1 + z2 (−a1 + (a2 − 1)z1 + z12 ), ż2 = z22 (a2 + z1 ). By a linear change of coordinates this system can be written as ż1 = −a2 z12 + a1 z13 − z12 z2 , ż2 = z2 − a1 z12 + (1 − a2 )z1 z2 + a1 z12 z2 − z1 z22 . Therefore (0, 0) is a semi–hyperbolic point. By Theorem 2.19 of [3] if a2 6= 0 we have that (0, 0) is a saddle–node point and going back to system (9) its local phase portrait is the one described in Figure 6. Figure 6. Local phase portrait at the origin of U2 at infinity. On the other hand, if a2 = 0 but a1 6= 0, (0, 0) is a stable topological node. Finally, when a2 = a1 = 0, the system becomes ż1 = −z12 z2 , ż2 = z2 (1 + z1 − z1 z2 ), so there is a line of singular points which is z2 = 0. After eliminating the common factor z2 (by a rescaling of the time) we do not obtain any singular point in the local chart U2 which is distinguished. 10 F.S. DIAS, J. LLIBRE AND C. VALLS 3.1.2. Finite singular points. In what follows we shall study the local stability at the singular points of system (8) for the parameters a1 and a2 . Case 1: a2 = 0. In this case system (8) becomes ẋ = −x − a1 z − xz, ż = −x. When a1 = 0 the system has a continuum of singular points (0, z) on the straight line x = 0. When a1 6= 0, the unique singular point is (0, 0). The eigenvalues are √ √ −1 + 1 + 4a1 −1 − 1 + 4a1 λ1 = and λ2 = . 2 2 As a1 > 0, the origin is a saddle. Case 2: a2 6= 0. In this case system (8) has exactly the two finite singular points:   a1 − a2 P0 = (0, 0), and P1 = a2 − a1 , . a2 The eigenvalues of the linear part at these singular points are p −(a2 + 1) ± (a2 + 1)2 + 4(a1 − a2 ) for P0 , and 2 p −(a1 + a22 ) ± (a1 + a22 )2 + 4a22 (a2 − a1 ) for P1 . 2 When a2 < a1 , P0 is a saddle and P1 is an attractor node. On the other hand, when a1 < a2 , P0 is an attractor node and P1 is a saddle. Finally, when a1 = a2 , we have that P0 = P1 . It is semi-hyperbolic and using Theorem 2.19 of [3] we conclude that it is a saddle-node. In short we have Table 2. Parameters Finite singular points Figures a1 = a2 = 0 line x = 0 Figure 7 (a) a2 = 0, a1 > 0 P0 saddle Figure 7 (b) a2 > 0, a1 = a2 P0 = P1 saddle–node Figure 7 (c) a2 > 0, a1 < a2 P0 stable node, P1 saddle Figure 7 (d) a2 > 0, a1 > a2 P0 saddle , P1 stable node Figure 7 (d) Table 2. Finite singular points for the differential sys- tem (8) and for the parameters a1 = βkλ/d3 and a2 = µ/d if a = d 6= 0, and a1 = βkλ/d3 and a2 = a/d if µ = d 6= 0. 11 Taking into account the local information on the finite and infinite singular points of system (8) we get that the global phase portraits of systems (8) are topologically equivalent to those described in the Figure 7. E2 E2 (a) (b) q′ q q′ P0 q E2′ E2′ E2 E2 (c) (d) q′ P0 q q′ P0 q P1 E2′ E2′ Figure 7. Phase portrait in the Poincaré disc of the virus system (1) on the invariant planes F1 = x + y − λ/d = 0 if a = d 6= 0. 3.2. Phase portrait on F2 (x, y, z) = 0. In this subsection we consider system (1) with µ = d 6= 0 and βkd 6= 0. System (1) restricted to the invariant plane F2 (x, y, z) = 0 is given by the differential system kλ (10) ẋ = λ − dx − βxz, ż = −kx − az + . d Doing the change of variables in (7) system (10) can be written as (11) ẋ = −x − a1 z − xz, ż = −x − a2 z, where a1 = βkλ/d3 and a2 = a/d and again we have renamed (X, Z) as (x, z), and the dot means derivative with respect to the new time T . 12 F.S. DIAS, J. LLIBRE AND C. VALLS Systems (8) and (11) are the same and so their global phase portraits are exactly the same. In short, taking into account the results of this subsection we have proved Theorem 3. 4. Proof of Theorem 4 In this section we describe the α– and the ω–limit sets of all orbits of system (1) contained in B and outside the invariant planes F1 = 0 and F2 = 0. In section 3 we proved that Id (x, y, z, t) = Fi (x, y, z)edt is an invariant of system (1) when a = d 6= 0 if i = 1 and µ = d 6= 0 if i = 2. Let γ(t) = {(x(t), y(t), z(t)) : t ∈ R} be an orbit of the compactified system (1) outside the planes Fi = 0 for i = 1, 2. Therefore Id (t, x(t), y(t), z(t)) = h. with h ∈ R \ {0}. So (12) Fi (x(t), y(t), z(t)) edt = h, for all t ∈ R, i = {1, 2}. Taking the limit in (12) when t → −∞, we obtain that lim Fi (x(t), y(t), z(t)) = ±∞. t→∞ Therefore the α–limit is contained in the boundary of the plane Fi = 0 at infinity. Then from Figure 1 these α-limit sets can only be one of the four singular points contained in the boundary of Fi = 0 at infinity, i.e., q, q ′ , E2 and E2′ (see Figure 7). So, statement (a) is proved. Now taking the limit in (12) when t → ∞, we obtain that lim Fi (x(t), y(t), z(t)) = 0. t→∞ Therefore the ω–limit is contained in the closed plane Fi = 0 in the Poincaré ball. To complete the study of the ω-limit sets we first study the finite singular points in R3 whenever either d = a, or d = µ. Since the study is the same we will focus in the case in which d = a. We consider different cases. Case 1: λ = µ = 0. In this case the singular points are located at x = y = 0. Computing the eigenvalues of the Jacobian matrix at the points of this line, we get that they are 0, −a and −a − zβ. It follows from the normal hyperbolicity theory (see again [4]) that if −a−z0 β < 0 then the singular point (0, 0, z0 ) has a 2-dimensional stable manifold, whereas if −a − z0 β > 0 then the singular point (0, 0, z0 ) has 13 a 1-dimensional stable manifold and 1-dimensional unstable invariant manifold. Case 2: µ = 0 and λ 6= 0. In this case the unique singular point is (λ/a, 0, 0). Computing the eigenvalues of the Jacobian matrix at this point we get p −a2 ± a4 + 4akβλ −a, . 2a Since akβλ > 0 we get that the singular point (λ/a, 0, 0) has a 2- dimensional stable manifold and 1-dimensional unstable manifold. Case 3: µ = kβλ/a2 . In this case the unique singular point is again (λ/a, 0, 0). Computing the eigenvalues of the Jacobian matrix at this point we obtain a3 + kβλ 0, −a, − a2 So by the normal hyperbolicity theory the singular point (λ/a, 0, 0) has either a 2-dimensional stable manifold. Case 4: µ 6= kβλ/a2 . In this case there are two singular points which are given in (5) (p0 and p1 ). Computing the eigenvalues of the Jacobian matrix at these singular points we get p a(a + µ) ± (a2 + aµ)2 − 4a(a2 µ − kβλ) −a, − 2a for p0 , and p aµ2 + kβλ ± (aµ2 + kβλ)2 + 4au2 (a2 µ − kβλ) −a, − 2aµ for p1 . So, p0 has a 2-dimensional stable manifold and 1-dimensional unstable manifold when a2 µ − kβλ < 0, and a 3-dimensional stable manifold when a2 µ − kβλ > 0. On the other hand, p1 has a 3- dimensional stable manifold when a2 µ − kβλ < 0, and a 2-dimensional stable manifold and 1-dimensional unstable manifold when a2 µ−kβλ > 0. Now we can finish the proof of statement (b) of Theorem 4. Indeed, the first assertion in statement (b) of Theorem 4 follows directly from Case 1 and Figure 7 whereas the second assertion in statement (b) of Theorem 4 follows directly from Case 2 and Figure 7. Moreover, the third and four assertions in statement (b) of Theorem 4 follow from Cases 3 and 4 and Figure 7. 14 F.S. DIAS, J. LLIBRE AND C. VALLS Acknowledgements The first author is partially supported by FAPEMIG Grants APQ- 01086-15 and APQ-01158-17. The second author is partially supported by the Ministerio de Economı́a, Industria y Competitividad, Agen- cia Estatal de Investigación grant MTM2016-77278-P (FEDER), the Agència de Gestió d’Ajuts Universitaris i de Recerca grant 2017 SGR 1617, and the European project Dynamics-H2020-MSCA-RISE-2017- 777911. The third author is partially supported by FCT/Portugal through UID/MAT/04459/2013. References [1] J. Alavez-Ramı́rez, G. Blé, V. Castellanos and J. Llibre, On the global flow of a 3-dimensional Lotka–Volterra system, Nonlinear Anal. 75 (2012), 4114–4125. [2] A. Cima and J. Llibre, Bounded polynomial vector fields, Trans. Amer. Math. Soc. 318 (1990), 557–579. [3] F. Dumortier, J. Llibre and J.C. Artés, Qualitative Theory of Planar Differential Systems, Springer Verlag, New York, 2006. [4] H.W. Hirsch, C.C. Pugh and M. Shub, Invariant manfifolds, Lecture Notes in Math. 583, Springer Verlag, Berlin, 1977. [5] P. De Leenheer and H.L. Smith, Virus dynamics: a global analysis, SIAM J. Appl. Math. 63 (2003), 1313–1327. [6] J. Llibre, M. Messias and A.C. Reinol, Quadratic polynomial differential systems in R3 having invariant planes with total multiplicity 9, Rendiconti del Circolo Matematico di Palermo 67 (2018), 569–580. [7] R.M. May and M.A. Nowak, Virus Dynamics. Mathematical Principles of Immunology and Virology, Oxford University Press, New York, 2000. [8] L. Min, Y. Su and Y. Kuang, Mathematical analysis of a basic virus in- fection model with application to HBV infection, Rocky Mount. J. Math. 38 (2008), 1–13. [9] C. Valls, Invariant algebraic surfaces fo a virus dynamics, Z. Angew. Math. Phys. 66 (2015), 1315–1328. 1 Instituto de Matemática e Computação, Universidade Federal de Itajubá, Avenida BPS 1303, Pinheirinho, CEP 37.500–903, Itajubá, MG, Brazil. E-mail address: scalco@unifei.edu.br 2 Departament de Matemàtiques, Universitat Autònoma de Barcelona, 08193 Bellaterra, Barcelona, Spain. E-mail address: jllibre@mat.uab.cat 3 Departamento de Matemàtica, Instituto Superior Técnico, Univer- sidade de Lisboa, 1049-001 Lisboa, Portugal. E-mail address: cvalls@math.ist.utl.pt