This is a preprint of: “Phase portraits of the Higgins-Selkov system”, Jaume Llibre, Marzieh Mousavi, Discrete Contin. Dyn. Syst. Ser. B, 2021. DOI: [10.3934/dcdsb.2021039] PHASE PORTRAITS OF THE HIGGINS–SELKOV SYSTEM JAUME LLIBRE Departament de Matemàtiques, Universitat Autònoma de Barcelona 08193 Bellaterra, Barcelona, Catalonia, Spain MARZIEH MOUSAVI Department of Mathematical Sciences, Isfahan University of Technology, Isfahan 84156-83111, Iran Abstract. In this paper we study the dynamics of the Higgins–Selkov system ẋ = 1 − xy γ , ẏ = αy(xy γ−1 − 1), where α is a real parameter and γ > 1 is an integer. We classify the phase por- traits of this system for γ = 3, 4, 5, 6, in the Poincaré disc for all the values of the parameter α. Moreover, we determine in function of the parameter α the regions of the phase space with biological meaning. 1. Introduction and statement of the main results The glycolysis is the first step in the transformation of glucose into energy for cellular metabolism. Higgins in [9] used a mathematical model for investigating sustained oscillations in the glycolysis process. But his model did not show any limit cycle corresponding to sustained oscillations observed in the experiments. For solving this problem Selkov in [12] introduced another mathematical model for studying the glycolysis now called the Higgins–Selkov system. This system is dx dy = 1 − xy γ , = αy(xy γ−1 − 1), (1) dt dt where x and y are state variables, t is the time variable, α is a real parameter and γ > 1 is an integer. The point (1, 1) is the unique singular point of system (1) for all values of parameter α 6= 0 and γ > 1. It is shown that a supercritical Hopf bifurcation occurs at α = 1/(γ − 1). In 2018, Artés et al. [2] characterized the dynamics of system (1) with γ = 2 and α ∈ R\(1, 3), and proposed a conjecture about existence and uniqueness of a limit cycle when α ∈ (1, 3). Chen and Tang in [6] have proved this conjecture which completes the phase portraits of system (1) with γ = 2. Brechmann and Rendall in [3] researched the uniqueness of limit cycles, and additionally they proved that there are no limit cycles when α ∈ (0, 1/(γ − 1)) for all γ > 1. Recently, it appears in arXiv the preprint of Brechmann and Rendall [4] where they show that for α > 0 the top right hand quarter of each phase portrait is topologically equivalent to that for γ = 2 for any integer γ > 1. 2010 Mathematics Subject Classification. 34C05, 34C07, 34C08. Key words and phrases. Higgins–Selkov system, limit cycle, Phase portrait, Poincaré compactifica- tion. 1 2 DYNAMICS OF HIGGINS–SELKOV SYSTEMS In this paper we provide a complete study of the dynamics of system (1) for γ = 3, 4, 5 and 6, and characterize the existence or non–existence of limit cycles in function of the parameter α. Moreover, we provide the qualitative behaviour of the system in a neighborhood of the infinity for all values of the integer γ > 1 using the Poincaré compactification. Finally we classify the phase portraits of system (1) in the Poincaré disc for γ = 3, 4, 5 and 6. For a definition of the Poincaré compactification and the Poincaré disc see Chapter 5 of [7], for the notions of separatrix, canonical regions and the concept of topologically equivalent differential systems see section 1.9 of [7]. Here, the number of separatrices and canonical region of a phase portrait will be denoted by S and R, respectively. Our main results are the following two theorems. Theorem 1. The phase portraits of the Higgins–Selkov system (1) in the Poincaré disc for γ = 3 and 5 are topologically equivalent to one of the phase portraits in Figure 1 as follows: (i) Figure 1(a) corresponds to α < 0 with S = 19 and R = 6, (ii) Figure 1(b) corresponds to α = 0 with S = ∞, (iii) Figure 1(c) corresponds to α ∈ (0, 1/(γ − 1)] with S = 17 and R = 4, (iv) Figure 1(d) corresponds to α ∈ (1/(γ − 1), αγ∗ ) where αγ∗ ∈ (1/(γ − 1), αγ ), √ √ α3 = (1261 + 57 57)/686 and α5 = 24(411 + 41 41)/3125 with S = 18 and R = 5, (v) Figure 1(e) corresponds to α = αγ∗ with S = 16 and R = 4, (vi) Figure 1(f ) corresponds to α > αγ∗ with S = 17 and R = 4. Moreover for all odd integer γ > 1 the phase portraits of Figure 1(a) − (c) hold for α < 0, α = 0 and α ∈ [0, 1/(γ − 1)), respectively; and the qualitative behaviour of system (1) at infinity is the same for any given α ∈ R. Theorem 2. The phase portraits of the Higgins–Selkov system (1) in the Poincaré disc for γ = 4 and 6 are topologically equivalent to one of the phase portraits in Figure 2 as follows: (i) Figure 2(a) corresponds to α < 0 with S = 19 and R = 6, (ii) Figure 2(b) corresponds to α = 0 with S = ∞, (iii) Figure 2(c) corresponds to α ∈ (0, 1/(γ − 1)] with S = 17 and R = 4 , (iv) Figure 2(d) corresponds to α ∈ (1/(γ − 1), αγ∗ ) where αγ∗ ∈ (1/(γ − 1), αγ ), √ α4 = 46875/14641 and α6 = (189384302 + 10121703 681)/47045881 with S = 18 and R = 5, (v) Figure 2(e) corresponds to α = αγ∗ with S = 16 and R = 4, (vi) Figure 2(f ) corresponds to α > αγ∗ with S = 17 and R = 4. Moreover for all even integer γ > 1 the phase portraits of Figure 1(a) − (c) hold for α < 0, α = 0 and α ∈ [0, 1/(γ − 1)), respectively; and the qualitative behaviour of system (1) at infinity is the same for any given α ∈ R. DYNAMICS OF HIGGINS–SELKOV SYSTEMS 3 Figure 1. The phase portraits of system (1) for γ = 3 and 5 in the Poincaré disc. The shaded areas correspond to the initial conditions of the orbits hav- ing a final finite evolution, so these are the initial conditions with biological meaning. In the phase portrait (c) the final behaviour is a stable singular point, and in the phase portrait (d) is a stable limit cycle. Figure 2. The phase portraits of system (1) for γ = 4 and 6 in the Poincaré disc. The shaded areas correspond to the initial conditions of the orbits hav- ing a final finite evolution, so these are the initial conditions with biological meaning. In the phase portrait (c) the final behaviour is a stable singular point, and in the phase portrait (d) is a stable limit cycle. 1 2. The proof of Theorems 1 and 2 2.1. The local phase portraits of the finite singular points. The Higgins– Selkov system (1) has a unique singular point (1, 1) and the Jacobian matrix of the 4 DYNAMICS OF HIGGINS–SELKOV SYSTEMS system at this point is   −1 −γ J= . α α(γ − 1) The determinant and trace of the matrix J are equal to α and α(γ−1)−1, respectively. The eigenvalues of the Jacobian matrix are 1h p i α(γ − 1) − 1 ± [α(γ − 1) − 1]2 − 4α . 2 Therefore the singular point (1, 1) is a hyperbolic saddle if α < 0;  √ 2 a stable hyperbolic node if α ∈ 0, (γ + 1 − 2 γ)/(γ − 1) ; √ 2 a stable focus if α ∈ (γ + 1 − 2 γ)/(γ − 1) , 1/(γ − 1) ; a stable weak focus if α = 1/(γ − 1), in this case a Hopf bifurcation occurs;  √ an unstable hyperbolic focus if α ∈ 1/(γ − 1), (γ + 1 + 2 γ)/(γ − 1)2 ; and √ an unstable hyperbolic node if α ≥ (γ + 1 + 2 γ)/(γ − 1)2 . For α = 1/(γ − 1) the Jacobian matrix J has a pair of imaginary eigenvalues. If it is considered the real part of the eigenvalues as a function of parameter α then the derivative of this function with respect to α at α = 1/(γ − 1) is nonzero. Therefore, according to Theorem 3.4.2 of [8] a Hopf bifurcation occurs. Then, by using the formula (3.4.11) of [8] we obtain the first Lyapunov coefficient which is equal to −1/16. Hence, we have a supercritical bifurcation for α slightly larger than 1/(γ −1), and the unique bifurcated limit cycle is stable. Moreover, if α = 0, then all the points on the curve x = y −γ are singular points of system (1) and all the parallel lines to the x-axis are invariant, therefore the phase portraits of system (1) are given in Figure 1(b) and 2(b) corresponding to the cases that the integer γ is odd and even, respectively. 2.2. The local phase portraits of the infinite singular points. In order to study the qualitative behaviour of the solutions of 2–dimensional polynomial differential systems at infinity, we can use the Poincaré compactification technique. The base of this technique is projecting a vector field in R2 onto a sphere which is a compact set. For more details about the Poincaré compactification, see chapter 5 of [7]. In the local chart U1 system (1) becomes u̇ = −uv γ+1 + uγ+1 + αuγ − αuv γ , v̇ = −v(v γ+1 − uγ ). (2) Let α 6= 0, then the infinite singular points of system (1) in local chart U1 are (0, 0) and (−α, 0). The Jacobian matrix of system (2) at point (0, 0) is identically zero. The blow up technique is used for studying the dynamics of the system in a neighborhood of the linearly zero singular points. For more details about the blow up, see [1]. Here we shall apply the blow up technique separately when the integer γ is odd and even. If γ is odd, for desingularizing the origin of system (2), we use two vertical blow ups. First we do a vertical blow up (u, v) 7→ (w1 v, v) where w1 is the new variable. Therefore system (2) in the new variables and after cancelling the common factor v γ−1 becomes ẇ1 = −αw1 (v − w1γ−1 ), v̇ = −v 2 (v − w1γ ), (3) where the origin is the unique singular point on the line v = 0, and the Jacobian matrix of the system at (0,0) is again identically zero. Hence, we use another vertical DYNAMICS OF HIGGINS–SELKOV SYSTEMS 5 blow up (w1 , v) 7→ (w2 v, v) where w2 is the new variable. Thus system (3) in the new variables and after eliminating the common factor v is ẇ2 = w2 (v − α − v γ w2γ + αv γ−2 w2γ−1 ), v̇ = v 2 (−1 + v γ−1 w2γ ). (4) The origin is the unique singular point of system (4) on the line v = 0, it is semi- hyperbolic and its local phase portrait is determined using Theorem 2.19 of [7]. By going back from (w2 , v) to (w1 , v) and from (w1 , v) to (u, v), we get that the point (0, 0) of the original system (2) is the union of one parabolic and four hyperbolic sectors when α > 0, see Figure 1(c)-(f); and the origin is a stable node when α < 0, see Figure 1(a). If γ is even we do a horizontal blow up (u, v) 7→ (u, uw) where w is the new variable. Then system (2) in the new variables writes u̇ = u(u − u2 wγ+1 + α − αuwγ ), ẇ = αw(−1 + uwγ ). (5) where by rescaling of time, the common factor uγ−1 has been eliminated. The point (0, 0) is the unique singular point of system (5) on u = 0, it is a hyperbolic saddle when α 6= 0. Now by performing a blowing down and returning from (u, w) to (u, v), we obtain that the origin of local chart U1 is formed by two hyperbolic and two parabolic sectors, see Figure 2(a) and 2(c)-(f). Now we consider the singular point (−α, 0) of system (2) in the local chart U1 . If γ is odd then the Jacobian matrix at this point has two positive eigenvalues when α < 0. Therefore (−α, 0) is a hyperbolic unstable node, see Figure 1(a). Also, when α > 0, the eigenvalues of the Jacobian matrix are negative and this implies that the point (−α, 0) is a hyperbolic stable node, see Figure 1(c)-(f). If γ is even, then the Jacobian matrix of system (2) at the point (−α, 0) has two positive eigenvalues when α 6= 0. It implies that (−α, 0) is a hyperbolic unstable node, see Figure 2(a) and 2(c)-(f). In the local chart U2 system (1) becomes u̇ = v γ+1 − u − αu2 + αuv γ , v̇ = −αv(u − v γ ). (6) The origin is a singular point of system (6) and the Jacobian matrix at this point has a unique zero eigenvalue. Therefore the origin is a semi-hyperbolic singular point and for determining the local behaviour of system (6) in the neighborhood of the origin we use Theorem 2.19 of [7] and we obtain that if γ is odd, then the origin is union of one parabolic and two hyperbolic sectors when α 6= 0, see Figure 1. In the case γ is even, the origin is a saddle when α > 0 and a stable node when α < 0, see Figure 2. 2.3. The periodic orbits. According to Corollary 2 of section 3.12 of [11], we know that the region limited by a periodic orbit of a C 1 planar differential system contains at least one singular point and if there are a finite number of singular points in this region, then the sum of the topological indices of these singular points is equal to one. Since system (1) for α < 0 has only one finite singular point (a saddle point) and the index of this point is −1 (see Theorem 7 of section 3.12 of [11]), then system (1) for α < 0 has no periodic orbits. Moreover if α = 0 then ẏ = 0, and in this case system (1) has no periodic solutions. 6 DYNAMICS OF HIGGINS–SELKOV SYSTEMS Lemma 3. The Higgins–Selkov system (1) for √ γ = 3, 4, 5 and 6, has no periodic solutions if √ α > αγ , where α3 = (1261 + 57 57)/686,√α4 = 46875/14641, α5 = 24(411 + 41 41)/3125 and α6 = (189384302 + 10121703 681)/47045881. Proof. For system (1) we have ẋ|x=0 = 1 > 0 and ẏ|y=0 = 0, and also if there exists a periodic orbit, it must surround the singular point (1, 1). Therefore if system (1) has a periodic solution, then it must be contained in the quadrant Q = {(x, y) ∈ R2 : x > 0, y > 0} for all γ > 1. The divergence of system (1) is f (x, y) = −y γ + αγxy γ−1 − α. The curve f (x, y) = 0 separates the quadrant Q into two open components f (x, y) > 0 and f (x, y) < 0. We prove that for γ = 3, 4, 5, 6, this curve {(x, √ y) : f (x, y) = 0} ∩ Q is transversal √ to the flow of system (1) if α > √ (1261 + 57 57)/686, 46875/14641, 24(411+41 41)/3125 and (189384302+10121703 681)/47045881, respectively. There- fore a periodic orbit cannot intersect the curve f (x, y) = 0, and it must be contained in a simple connected region either {f (x, y) > 0} ∩ Q, or {f (x, y) < 0} ∩ Q. Then the existence of a such periodic solution is impossible by the Bendixson criterion, see Theorem 7.10 of [7]. By solving f (x, y) = 0 we obtain x = (y γ + α)/αγy γ−1 . Consider the function   ∂f ∂f 1 Tα (x, y) = ẋ + ẏ y γ +α := Pα (y) ∂x ∂y x= αγy γ−1 γ −(1 + γ)y 2γ + α(γ − 2)y γ + αγ 2 y γ−1 − α2 (γ − 1)2 = . γ Hence for all positive integer γ, Pα is a polynomial of degree 2γ. In order to prove that Tα (x, y) in the first quadrant does not change sign, it is sufficient to verify that polynomial Pα has no real roots. Here, we use the notions of discrimination matrix, discriminant sequence and Theorem 8 in the Appendix to prove that Pα for γ = 3, 4, 5, 6 does not have real roots. First we√consider γ = 3 and prove that the polynomial Pα of degree 6 for α > (1261 + 57 57)/686 has no real roots. The elements of the discriminant sequence of polynomial Pα are D1 = d2 = 96, D2 = d4 = 0, D3 = d6 = −13824α2 , D4 = d8 = 186624α3 (−1536 + 7α), D5 = d10 = 30233088α5 (−1024 + 621α), D6 = d12 = −2176782336α8 (1024 − 1261α + 343α2 ). √ Therefore, according to the Appendix for α > (1261 + 57 57)/686 the sign list of sequence {D1 , D2 , ..., D6 } is [1, 0, −1, sign (D4 ), 1, −1], so we have the following revised sign list [1, −1, −1, 1 or − 1, 1, −1]. DYNAMICS OF HIGGINS–SELKOV SYSTEMS 7 Then for any sign of the fourth element in the above list, the number of the sign changes in the revised list is 3. Hence, according to Theorem 8, Pα (y) has 3 pairs of distinct conjugate imaginary roots. Hence it has no real roots and this implies that Tα (x, y) < 0 in Q. This concludes the proof of the lemma for γ = 3. Now we consider γ = 4. Again by using discrimination matrix and Theorem 8 we show that the polynomial Pα of degree 8 for α > 46875/14641 has no real roots. The elements of the discriminant sequence of polynomial Pα are D1 = d2 = 200, D2 = d4 = 0, D3 = d6 = 0, D4 = d8 = −12800000α3 , D5 = d10 = 409600000α4 (12500 + 11α), D6 = d12 = 39321600000α6 (62500 + 4387α), D7 = d14 = −11324620800000α9 (−15625 + 48257α), D8 = d16 = 9784472371200000α12 (1 + α)(−46875 + 14641α). For α > 46875/14641 the sign list of sequence {D1 , D2 , ..., D8 } is [1, 0, 0, −1, 1, 1, −1, 1], Therefore the corresponding revised sign list is [1, −1, −1, −1, 1, 1, −1, 1]. The number of the sign changes in this revised list is 4. Thus Pα (y) has 4 pairs of distinct conjugate imaginary roots and no real roots. It implies that Tα (x, y) < 0 in Q. This concludes the proof of the lemma for γ = 4. If γ = 5 in a similar way to √ the previous cases, we show that the polynomial Pα of degree 10 for α > 24(411 + 41 41)/3125 has only imaginary roots. The elements of the discriminant sequence of the polynomial Pα are D1 = d2 = 360, D2 = d4 = 0, D3 = d6 = 0, D4 = d8 = 0, D5 = d10 = 23619600000α4 , D6 = d12 = −22143375000000α5 (−9600 + α), D7 = d14 = −6643012500000000α7 (−48000 + 343α), D8 = d16 = −326517350400000000000000α10 , D9 = d18 = 26572050000000000000000α13 (−20352 + 8575α), D10 = d20 = −70858800000000000000000000α16 (1 + α)(18432 − 19728α + 3125α2 ). √ For α > 24(411 + 41 41)/3125 the sign list of sequence {D1 , D2 , ..., D10 } is [1, 0, 0, 0, 1, sign (D6 ), sign (D7 ), −1, 1, −1], so the corresponding revised sign list is [1, −1, −1, 1, 1, 1 or − 1, 1 or − 1, −1, 1, −1]. 8 DYNAMICS OF HIGGINS–SELKOV SYSTEMS The case that sign (D6 ) = −1 and sign (D7 ) = 1 is impossible, in the other cases the number of the sign changes in this revised list is 5. Therefore Pα (y) has 5 pairs of distinct conjugate imaginary roots and no real roots. It implies that Tα (x, y) < 0 in Q. This completes the proof of the lemma for γ = 5. For γ = 6 the elements of discriminant sequence of the polynomial Pα of degree 12 are D1 = d2 = 588, D2 = d4 = 0, D3 = d6 = 0, D4 = d8 = 0, D5 = d10 = 0, D6 = d12 = 78690759081984α5 , D7 = d14 = −265581311901696α6 (66706983 + 608α), D8 = d16 = −430241725280747520α8 (155649627 + 93272α), D9 = d18 = −15944515252329719468851200000α11 , D10 = d20 = 408976816222257304376033280000000α14 , D11 = d22 = −753691861204384855818240000000α17 (−1647086 + 898529α), D12 = d24 = 38108015454154290462720000000000α20 (−720600125 − 378768604α + 47045881α2 ). √ For α > (189384302+10121703 681)/47045881 the sign list of sequence {D1 , D2 , ..., D12 } is [1, 0, 0, 0, 0, 1, −1, −1, −1, 1, −1, 1], then the corresponding revised sign list is [1, −1, −1, 1, 1, 1, −1, −1, −1, 1, −1, 1]. The number of the sign changes in this revised list is 6. According to Theorem 8, Pα (y) has 6 pairs of distinct conjugate imaginary roots. It implies that Tα (x, y) < 0 in Q. This completes the proof of the lemma for γ = 6.  To simplify the investigation about the limit cycles of the Higgins–Selkov system (1), it is introduced the following change of variables α−X −Y dt 1 x= , y = 1 + X, = γ. α dτ y By performing this change of variables and again denoting X, Y and τ by x, y and t, respectively, system (1) transforms into the following Liénard differential system dx α dy αx = −y + α − x − , = . (7) dt (1 + x)γ−1 dt (1 + x)γ In [3] by considering Liénard system (7) for x > −1, it was proved that if the unique singular point of system (1) is unstable for a given value of parameter α, then there is at most one periodic orbit and if such a solution exists, it must be asymptotically stable. Moreover the authors of [3], by means of Theorem 2.1 of [10], have shown if α ∈ (0, 1/(γ − 1)) then system (7) has no limit cycles. This implies that when the unique singular point of the Higgins–Selkov system (1) is stable, then there is no limit cycles (see Theorem 7 of [3]). DYNAMICS OF HIGGINS–SELKOV SYSTEMS 9 In the following result we show that system (7) and consequently the Higgins– Selkov (1) in the case γ = 3, 4, 5, 6 has no limit cycle when α = 1/(γ − 1). Lemma 4. If α = 1/(γ − 1) then system (7) for γ = 3, 4, 5, 6 has no limit cycles. Proof. By doing the change of variable y → −y on system (7), we can use Corollary 10 of [5] such that α αx F̄ (x) = α − x − γ−1 , ḡ(x) = , (1 + x) (1 + x)γ and x ∈ (−1, ∞). If α = 1/(γ − 1), it can be easily checked that the conditions (i)-(iii) of Proposition 9 of [5] hold. To complete this proof, it suffices to show that there is no solution for the system F̄ (x1 ) = F̄ (x2 ), λ̄(x1 ) = λ̄(x2 ), (8) where F̄ 0 (x) = f¯(x), λ̄(x) = ḡ(x)/f¯(x), and −1 < x2 < 0 < x1 . From λ̄(x1 ) = λ̄(x2 ) and α = 1/(γ − 1), the following equation is obtained (x2 − x1 ) + x1 (1 + x2 )γ − x2 (1 + x1 )γ = 0, by using the binomial expansion and after simplifying, we can rewrite the above equation as        γ γ 2 2 γ γ−1 γ−1 x1 x2 (x2 − x1 ) + (x2 − x1 ) + ... + (x2 − x1 ) = 0, (9) 2 3 γ  where nγ = (γ−n)!n! γ! and n < γ is a positive integer. Now for γ = 3, 4, 5, 6, we show that there is no solution for equation (9) such that −1 < x2 < 0 < x1 , thus there is no solution for system (8). By substituting γ = 3, 4, 5, 6 into (9), we obtain • x1 x2 (x2 − x1 )(3 + x1 + x2 ) = 0 for γ = 3,   • x1 x2 (x2 − x1 ) (6 + 4x2 + x22 ) + x1 (4 + x2 ) + x21 = 0 for γ = 4,   • x1 x2 (x2 − x1 ) 10 + x31 + x2 (10 + 5x2 + x22 ) + x21 (5 + x2 ) + x1 (10 + 5x2 + x22 ) = 0 for γ = 5,  • x1 x2 (x2 − x1 ) (15 + 20x2 + 15x22 + 6x32 ) + x42 + x31 (6 + x2 ) + x21 (15 + 6x2 + x22 ) + x1 (20 + 15x2 + 6x22 + x32 ) + x41 = 0 for γ = 6. In all four equations listed above, the condition −1 < x2 < 0 < x1 implies that there is no solution for these equations. Therefore system (8) has no solutions and by using Corollary 10 of [5] we conclude that system (7) for γ = 3, 4, 5, 6 and α = 1/(γ − 1) has no limit cycles.  Now we study the existence of limit cycle of system (1) when α > √ 1/(γ − 1).√For this purpose we perform the linear change of variables (x, y, t) → (x, αy, −1/ αt) in system (7) when α > 0, and we obtain the following system √ dx √ 1 α dy −x =y+ α− √ x− γ−1 , = . (10) dt α (1 + x) dt (1 + x)γ Then we show that the vector field of system (10) is a generalized rotated vector field (see Definition 3.3 of [14]). 10 DYNAMICS OF HIGGINS–SELKOV SYSTEMS Lemma 5. If γ > 1 is a fixed integer, then the vector field of system (10) for all real α > 0 is a generalized rotated vector field. Proof. The unique singular point of system (10) at the origin for all values of the parameter α remains unchanged and also the value of the following determinant √ √ 1 α1 −x y + α1 − √ x − γ−1 α1 (1 + √x) (1 + x)γ α2 , √ 1 −x y + α2 − √ x − α2 (1 + x)γ−1 (1 + x)γ for all x > −1 and y ∈ R, is equal to   √ √ x((x + 1)γ−1 − 1) 1 x2 ( α2 − α1 ) +√ > 0, (1 + x)2γ−1 α1 α2 (1 + x)γ where α1 < α2 and α1 , α2 > 0. The value of this determinant can only be zero for x = 0. Therefore, according to Definition 3.3 of [14], the vector field of system (10) is a generalized rotated vector field.  Remark 6. Since system (10) has been obtained from system (7) by a linear change of variables, hence all facts about existence and number of limit cycles of the both systems are the same. Lemma 7. For each γ = 3, 4, 5 and 6, there exists a real number αγ∗ ∈ (1/(γ − 1), αγ ] such that system (10) has a unique limit cycle when α ∈ (1/(γ − 1), αγ∗ ) and no periodic solution when α ∈ (αγ∗ , αγ ] where αγ for γ = 3, 4, 5, 6, has been introduced in Lemma 3. Proof. According to Lemma 5 the vector field of system (10) is a generalized rotated vector field. Hence, if this system has a limit cycle for a given values of α, then by Theorem 3.5 of [6] this limit cycle can neither split nor disappear when parameter α changes monotonically and also the amplitude of the limit cycle will increase or decrease monotonically. As we know when α = 1/(γ − 1) a Hopf bifurcation occurs, therefore for α slightly larger than 1/(γ − 1) system (10) has a unique limit cycle. Theorem 3.5 of [6] implies that the amplitude of the limit cycle increases when parameter α increases. Moreover system (10) has no periodic orbit when α ∈ (−∞, 1/(γ − 1)] ∪ (αγ , +∞). Thus for each γ = 3, 4, 5, 6, there exists αγ∗ ∈ (1/(γ − 1), αγ ] such that when α tends to α∗ , the amplitude of the limit cycle tends to infinity. This completes the proof of the lemma.  2.4. Phase portraits in the Poincaré disc. By considering all information of system (1) which have been obtained until now, and taking into account that the straight line y = 0 is invariant, we get the phase portraits of Figures 1 and 2 described in Theorems 1 and 2, respectively. Therefore the proofs of Theorems 1 and 2 are completed. DYNAMICS OF HIGGINS–SELKOV SYSTEMS 11 3. Appendix In this appendix we state the definitions of the discrimination matrix and discrim- inant sequence which have been introduced in [13] and some notations and result related to them. For more details, see [13]. Definition 1. Consider a general polynomial of degree n, f (x) = a0 xn + a1 xn−1 + ... + an , the following (2n + 1) × (2n + 1) matrix that is formed by the coefficients of the polynomial,   a0 a1 a2 ... an  0 na0 (n − 1)a1 ... an−1     a0 a1 ... an−1 an     na0 ... 2an−2 an−1     ... ...   ... ...     a a ... a   0 1 n   0 na0 ... an−1  a0 a1 ... an where there are zeros in the blank spaces, is called the discrimination matrix of f (x), and denoted by Discr(f ). The determinant of the submatrix Discr(f ), formed by the first k rows and the first k columns, is denoted by dk , for k = 1, .., 2n + 1. Definition 2. The sequence {D1 , ..., Dk } is called the discriminant sequence of poly- nomial f , where Dk = d2k , for k = 1, ..., n. Definition 3. The list [sign(D1 ), ..., sign(Dn )] is called the sign list of a given se- quence {D1 , ..., Dk }, where   −1 if x < 0, sign (x) = 0 if x = 0,  +1 if x > 0. Definition 4. Let [s1 , s2 , ....sn ] be a sign list. The revised sign list [ε1 , ε1 , ..., εn ] is constructed according of the following rules: • If [si , si+1 , ..., si+j ] is a part of the given list, where si 6= 0, si+1 = si+2 = ... = si+j−1 = 0, si+j 6= 0; then, we replace the subsection [si+1 , si+2 , ..., si+j−1 ] with [−si , −si , si , si , −si , −si , si , si , ...], keeping the number of terms. • Otherwise, let εk = sk , i.e. no changes for other terms. Theorem 8. [13] Consider a polynomial f (x) which has real coefficients, f (x) = a0 xn + a1 xn−1 + ... + an . Let ν be the number of the sign changes of the revised sign list of {D1 , D2 , ..., Dk }. Then the number of the pairs of distinct conjugate imaginary roots of polynomial f (x) 12 DYNAMICS OF HIGGINS–SELKOV SYSTEMS is equal to ν. Moreover, if the number of non-zero members of the revised sign list is l, then, l − 2ν is the number of the distinct real roots of polynomial f (x). Acknowledgements We thank to the reviewers their comments and suggestions which help us to improve the presentation of this paper. The first author is partially supported by the Ministerio de Economı̀a, Indus- tria y competitividad, Agencia Estatal de Investigación grant MTM2016-77278-P (FEDER), the Agència de Gestió d’Ajusts Universitaris i de Recerca grant 2017 SGR 1617, and the European project Dynamics-H2020-MSCA-RISE-2017-777911. The second author is supported by Isfahan University of Technology (IUT). References [1] M. J. Álvarez, A. Ferragut and X. Jarque, A survey on the blow up technique, Int. J. Bifurcation and Chaos. 21 (2011), 3103–3118. [2] J. C. Arteś, J. Llibre and C. Valls, Dynamics of the Higgins–Selkov and Selkov system, Chaos. Sol. Frac. 114 (2018), 145–150. [3] P. Brechmann and A. D, Rendall, Dynamics of the Selkov oscillator, Mathematical Bio- sciences, 306 (2018), 152–159. [4] P. Brechmann and A. D, Rendall, Unbounded solutions of models for glycolysis, preprint arXiv (https://arxiv.org/abs/2003.07140), (2020). [5] H. Chen, J. Llibre and Y. Tang , Global dynamics of a SD oscillator, Nonlinear Dynamics, 91 (2018), 1755–1777. [6] H. Chen and Y. Tang, Proof of Artés-Llibre-Valls’s conjectures for the Higgins–Selkov and Selkov systems, J. Differential Equations 266 (2019), 7638–7657. [7] F. Dumortier, J. Llibre and J.C. Artés, Qualitative Theory of Planar Differential Sys- tems, Springer Verlag, New York, 2006. [8] J. Guckenheimer and P. Holmes, Nonlinear oscillations, dynamical systems, and bifurca- tions of vector fields, in: Applied Mathematical Sciences, vol. 42, Springer–Verlag, New York, 1986. [9] J. Higgins, A chemical mechanism for oscillation of glycolytic intermediates in yeast cells, Proc. Natl. Acad. Sci. (USA) 51 (1964), 989–994. [10] T. W. Hwang and H. J. Tsai, Uniqueness of limit cycles in theoretical models of certain oscillating chemical reactions, J. Phys. A 38 (2005), 8211–8223. [11] L. Perko, Differential Equations and Dynamical Systems, (3rd Ed., Springer, 2001). [12] E. E. Selkov, Self-oscillations in glycolysis. I. A simple kinetic model, Eur. J. Biochem.4 (1968), 79–86. [13] L. Yang, Recent advances on determining the number of real roots of parametric polynomials, J. Symb. Comput. 28 (1999), 225–242. [14] Z. Zhang, T. Ding and W. Huang, Qualitative Theory of Differential Equations, Transl. Math. Monogr. Amer. Soc., Providence, RI, 1992. Email address: jllibre@mat.uab.cat Email address: marzieh.mousavi@math.iut.ac.ir