This is a preprint of: “The local period function for Hamiltonian systems with applications”, Claudio A. Buzzi, Yagor Romano Carvalho, Armengol Gasull, J. Differential Equations, vol. 280, 590–617, 2021. DOI: [10.1016/j.jde.2021.01.033] THE LOCAL PERIOD FUNCTION FOR HAMILTONIAN SYSTEMS WITH APPLICATIONS. CLAUDIO A. BUZZI, YAGOR ROMANO CARVALHO, AND ARMENGOL GASULL Abstract. In the first part of the paper we develop a constructive procedure to obtain the Taylor expansion, in terms of the energy, of the period function for a non-degenerated center of any planar analytic Hamiltonian system. We apply it to several examples, including the whirling pendulum and a cubic Hamiltonian system. The knowledge of this Taylor expansion of the period function for this system is one of the key points to study the number of zeroes of an Abelian integral that controls the number of limit cycles bifurcating from the periodic orbits of a planar Hamiltonian system that is inspired by a physical model on capillarity. Several other classical tools, like for instance Chebyshev systems are applied to study this number of zeroes. The approach introduced can also be applied in other situations. 1. Introduction and main results Let γs , with s ∈ I ∈ R, be a parameterized continua of periodic orbits of a planar autonomous differential system. In general, I is either an open interval or an interval of the form [s0 , s1 ). The function that assigns to each s the minimal period of γs is called period function and it is denoted by T (s). Similarly, the function that assigns to each γs the area surrounded by this closed curve is denoted by A(s) and called area function. The period function is important to study theoretical properties of planar ordinary differential equations and their perturbations, see for instance [9, pp. 369-370]; to understand some mathematical models in physics or ecology, see [14, 17, 39, 45] and the references therein; in the description of the dynamics of some discrete dynamical systems, see [6, 11, 12]; or for counting the solutions of some boundary value problems, see [7, 8]. When the system is Hamiltonian, with Hamiltonian function H and γh ⊂ {H = h}, it is natural to consider s = h and write T = T (h). Given a planar analytic Hamiltonian system ẋ = −Hy (x, y), ẏ = Hx (x, y), (1) with a non-degenerated center at the origin (that without loss of generality we will as- sociated to h = 0 and then I = [0, h1 ) ⊂ R) it is known that T (h), in a neighborhood of h = 0, is an analytic function of the energy h and it is given by the derivative with respect h of the area function A(h), see [33]. There are several authors that compute the Taylor series of T at h = 0 for particular Hamiltonian systems but, to the best of our knowledge, most examples deal with Hamiltonian functions with separated variables H(x, y) = F (x) + G(y), see for instance [5, 16] and their references. Our first result provides a systematic constructive approach for finding this Taylor series up to any order for any Hamiltonian system. Theorem 1.1. Let H be an analytic function with H(0, 0) = 0 and assume that the Hamiltonian system (1) has a non-degenerate center at the origin. Then: 2010 Mathematics Subject Classification. Primary: 34C08. Secondary: 34C25, 37G15, 37J45. Key words and phrases. Period function; Limit cycles; Abelian integrals; Extended complete Cheby- shev systems; Picard-Fuchs differential equations. 1 2 C. A. BUZZI, Y. R. CARVALHO, AND A. GASULL (i) In a neighborhood of h = 0 the period function T (h) and the area function A(h) are analytic and satisfy A′ (h) = T (h). (ii) Let ⊤(ρ), with ρ > 0, be the period of the orbit starting at the point (x, y) = (0, ρ). Then: • The function ⊤ is analytic at ρ = 0 and its Taylor development at ρ = 0 can be obtained algorithmically from the expression of (1) in polar coordinates. • The equation H(0,√ ρ) = h has, in a positive neighborhood of h = 0, a positive solution ρ = S( h) with S analytic at zero. Moreover its Taylor development at zero can also be obtained √ algorithmically from the one of H(0, ρ). Rh • It holds that T (h) = ⊤(S( h )). Then A(h) = 0 T (s) ds. We want to stress that our contribution restricts to item (ii). As we have already said, item (i) is proved in [33]. We put both together for the sake of clarity. We remark that if the system is only of class C k , for some k ∈ N, our approach can be adapted to this setting providing several terms of the Taylor expansion of T (h) and A(h). It is also worthwhile to comment that sometimes it is also possible to obtain a closed integral expression for T (h) or to prove that it satisfies some differential equation. Hence other totally different approach consists in studying this integral and try to obtain its Taylor series at the origin or to study a particular solution of this differential equation. As far as we know, our two steps procedure is new, totally different and applicable to all integrable planar systems with a non-degenerated center. It is tedious, but so systematic that can be used with any computer algebra system. In this paper, these computations are done with Maple. In Section 2 we will apply it to several examples, including the whirling pendulum and a family of quadratic Hamiltonian systems. Our main application is given in the second part of this paper that studies, with this point of view, the number of limit cycles that appear in the study of a planar system motivated by a physical model on capillarity that we briefly describe. Capillary action is the physical property that fluids have to go down or up in extremely thin tubes. The capillary rise in a narrow vertical tube is a remarkable physical phenome- non that can also be observed in many other everyday situations, such as water transport in the soil or plants. The increased use of capillary flow as an application in the industry has made a substantial growth in the search for more appropriate mathematical models. For the description of the model and more details see [38]. In that paper we can see this model is  ′ u = v, √ (2) v ′ = 1 + Kv − 2u, and it is defined on u > 0, where the prime denotes a differentiation with respect to the  real time variable t and K ∈ R. Notice that (2) has a unique equilibrium point at 21 , 0 . √ Moreover when K = 0, it is a Hamiltonian system with H(u, v) = −u + 8u3 /3 + v 2 /2. Its critical point is a local minimum, so it is of center type. Consider the following perturbation √ of the Hamiltonian system (2), with K = 0, mo- tivated by the appearance of 2u in its expression,  ′ √ u = v + εPe( 2u), √ √ (3) v ′ = 1 − 2u + εQ(e 2u), √ √ √ e 2u) are polynomials in the variable 2u and defined on u > 0, where Pe( 2u) and Q( √ ε ∈ R is a small parameter. Now we perform the change of variable x = 2u − 1 and THE PERIOD FUNCTION WITH APPLICATIONS 3 y = −v, and a time rescaling. So, system (3) becomes  ′ x = −y + εP (x), (4) y ′ = x + x2 + εQ(x), defined on x + 1 > 0, with P and Q, again polynomials with respective degrees n + 1 and m. Clearly, when ε = 0 the above system is a Hamiltonian system, with a center at the origin, and x2 x3 y 2 H(x, y) = + + , (5) 2 3 2 the energy levels {H(x, y) = h} ∩ {x + 1 > 0} of the Hamiltonian (5) for h ∈ (0, 16 ) are ovals. By using the classical approach, see Subsection 3.1 for more details, the number of limit cycles that bifurcate from the periodic orbits γh , by a first order analysis in ε for system (4) is given by the maximum number of simple zeroes of the Abelian integral associated to (4), Z Z X n Z ′ I(h) = P (x) dy = − P (x)y dx = αj Ij (h) with Ij (h) = xj y dx, (6) γh γh j=0 γh where the parameters αj are arbitrary real parameters that depend on the ones of P. A first result relates the above integrals with the first part of our paper. It will be proved in Subsection 3.2. Proposition 1.2. Let A(h) and T (h) be the area and period functions associated to the Hamiltonian system with H(x, y) = x2 /2 + y 2 /2 + x3 /3. (i) The functions A and T verify the 2 × 2 Picard-Fuchs equations  ′     A (h) 0 6(6h − 1)h A(h) 6(6h − 1)h = . T ′ (h) −5 0 T (h) (ii) The function A verifies the Hill’s equation 5 A′′ (h) = A(h). (7) 6(1 − 6h)h (iii) The function p(h) = A(h)/T (h) satisfies the Riccati differential equation 5 p′ (h) = −p2 (h) + . 6(1 − 6h)h (iv) Consider the Abelian integrals defined in (6). Then I0 (h) = −A(h) and for all 0 ≤ n ∈ N, I3n+1 (h) =an+1 (h)A(h) + bn (h)(6h − 1)hT (h), I3n+2 (h) =cn+1 (h)A(h) + dn (h)(6h − 1)hT (h), I3n+3 (h) =en+1 (h)A(h) + fn (h)(6h − 1)hT (h), where g ∈ {a, b, c, d, e, f }, and gk ∈ Q(h) denotes a polynomial of degree k. Moreover, I3n+2 is a linear combination of several Ij with j < 3n + 2 and j 6≡ 2 (mod 3). We prove next theorem, see Subsection 3.4 to recall the definition of extended complete Chebyshev (ECT) system. Theorem 1.3. Set I(h) given in (6), where γh , for h ∈ L := (0, 1/6), are the ovals of {x2 /2 + x3 /3 + y 2 /2 = h}. Then, Xn n X I(h) = αj Ij (h) = aj Ij (h) j=0 j=0, j6≡2 (mod 3) 4 C. A. BUZZI, Y. R. CARVALHO, AND A. GASULL where aj can be taken free real parameters and depend on the coefficients of P. Moreover: (i) The functions Ij (h), j 6≡ 2 (mod 3), are linearly independent. In particular, there are values of the constants aj such that I(h) has N(n) simple zeroes in L where N(n) = Card{j : 0 ≤ j ≤ n and j 6≡ 2 (mod 3)} − 1 = n − [(n + 1)/3] and, as usual, [· ] denotes the integer part function. (ii) There exists h1 > 0 such that for each n ≤ 50 the functions Ij (h), 0 ≤ j ≤ n and j 6≡ 2 (mod 3), form an ECT system on (0, h1 ). In particular, for n ≤ 50 the maximum number of zeroes of I(h) in (0, h1 ), taking into account their multiplicities is N(n).   (iii) The functions I0 (h), I1 (h) and I3 (h), I1 (h), I0 (h) form ECT systems on L. In particular, for n = 0, 1, 2, 3, the maximum number of zeroes of I(h) in L, taking into account their multiplicities is N(n). A similar result could also be proved taking in (4) more general perturbations with P and Q depending also on y. In fact, while we were ending the above proof we realized that this more general problem was addresses many years ago by Petrov with an equivalent expression of the Hamiltonian (he took H(x, y) = y 2 + x3 − x), see [37], or the second part of the book [10] with the also equivalent expression H(x, y) = x − x3 /3 + y 2 /2. He got a more general result by using the complexification of the corresponding Ij (h), the Picard-Fuchs equations satisfied for them and the argument principle. In particular he proved that the maximum number of zeroes of I(h) in L, taking into account their multiplicities is always N(n). Our proof of item (i) is similar to the one of Petrov but the proof of item (ii) strongly uses that knowing the function T (h) near h = 0 suffices to get lower bounds of the number of zeroes of I(h) and gives a different computational approach to the problem, that is valid for a given n. We stop at n = 50, but it is easy to go further in our computations. In fact, although Petrov approach gives strong results in this case, as we will explain in Remark 4.1, our point of view can also be used for studying several perturbations of many Hamiltonian systems. In all these cases a lower bound of the number of limit cycles follows from the computation of the Wronskian of some polynomials on h at h = 0 and these polynomials can be obtained simply from the knowledge of the Taylor’s series T (h) or A(h) at h = 0. The key point for proving item (iii) is the method introduced in [26, 33]. It provides an alternative approach to the one of Petrov for small n. A straightforward corollary of Theorem 1.3 is: Corollary 1.4. For ε small enough, system (4) has at least n − [(n + 1)/3] limit cycles surrounding the origin, that bifurcate from their periodic orbits γh , h ∈ (0, 1/6). The paper is organized as follows: In Section 2 we prove Theorem 1.1 and apply it to several Hamiltonian systems. In Section 3 we include same preliminaries devoted to prove Theorem 1.3 and we also prove Proposition 1.2. More concretely, there is a subsection devoted to recall the relation between limit cycles and Abelian integrals; a second one dedicated to find the Picard-Fuchs, Hill and Riccati differential equations for I0 (h) and I1 (h), see Proposition 3.4, and to present other relations among all the involved Abelian integrals; a third one about the parameterization of genus 0 planar polynomials curves and its application to our problem; and the last one on Chebyshev systems and how to use them in our situation. Finally, in Section 4 we prove Theorem 1.3 and Corollary 1.4. 2. Proof of Theorem 1.1 and some applications Before proving Theorem 1.1 and for the sake of completeness we prove Q a preliminary result. Notice that the solution of an equation of the form F0 (w, z) = ni=1 (w − αi z) = 0, THE PERIOD FUNCTION WITH APPLICATIONS 5 with all αi different, is given by the n straight lines w = αi z, i = 1, 2, . . . , n. Next lemma asserts that when we consider a more general analytic equation of the form F (w, z) = bF0 (w, z) + O(n + 1) = 0, with b 6= 0, and the O(n + 1) part denotes terms with degree at least n + 1 in w and z, then its solutions near (0, 0) is given by n analytic branches w = Wi (z), i = 1, 2, . . . , n, that are tangent to these n lines. In fact, next result and more general ones can be obtained and proved by using the so called Newton’s polygon, see [3] for more details. Q Lemma 2.1. Consider an analytic function F (w, z) = b ni=1 (w − αi z) + O(n + 1), with b 6= 0 and all αi different. Then, in a neighborhood of (0, 0), the solutions of F (w, z) = 0 are given by n branches w = Wi (z) = αi z + O(z 2 ), i = 1, 2, . . . , n, where all the functions Wi are analytic at zero and their Taylor series at the origin can be obtained by implicit derivation. Proof. To solve F (w, z) = 0, we divide it by b and we make the blow up w = uz. Then we get the following equivalent equation n Y n z (u − αi ) + z n+1 h(u, z) = 0, i=1 with h(u, z) being an analytic function at (0, 0). Hence it suffices to consider Yn G(u, z) = (u − αi ) + zh(u, z) = 0. i=1 Since, for all j ∈ {1, . . . , n}, it holds that Yn ∂ G(αj , 0) = 0 and G (αj , 0) = (αj − αi ) 6= 0, ∂u i=1, i6=j by the Implicit Function Theorem it follows that for every j there is an analytic function Uj in variable z, that satisfies Uj (0) = αl and G(Uj (z), z) = 0 for all z in a neighborhood of 0. As w = uz in a neighborhood of (0, 0), the solutions of F (w, z) = 0 are w = Wj (z) = zUj (z) = αj z + O(z 2 ), and the Taylor series of each of them can be obtained simply by implicit derivation.  Proof of Theorem 1.1. (i) This result is proved in [33]. (ii) We will prove that the function ⊤(ρ) is analytic at ρ = 0 for any non-degenerated center, not necessarily Hamiltonian. We simply follow the approach developed in [3], see also [18]. It is not restrictive to write the differential system as  ẋ = −y + f (x, y), (8) ẏ = x + g(x, y), where f and g are analytic functions in a neighborhood of the origin starting with terms at least of degree two. Passing system (8) to polar coordinates (r, θ) we obtain ṙ = S(r, θ) and θ̇ = 1+T1 (r, θ). Now we leave θ as the new independent variable and so we get the following analytic differential equation dr S(r, θ) = = F2 (θ)r 2 + · · · + Fn (θ)r n + O(r n+1). (9) dθ 1 + T1 (r, θ) We consider the initial condition r(0) = ρ > 0. We can write the solution of (9) with this initial condition as r(θ, ρ) = rρ (θ) = ρ + u1(θ)ρ + · · · + un (θ)ρn + O(ρn+1 ), 6 C. A. BUZZI, Y. R. CARVALHO, AND A. GASULL which is also analytic. By plugging it in equation (9) we find each ui (θ), i ≥ 1, by solving simple differential equations with initial conditions ui (0) = 0. Notice also that (1+T1 (r, θ))−1 is analytic at r = 0 because T1 (0, θ) ≡ 0 and so, in a suitable neighborhood of r = 0 X∞ 1 = 1+ gk (θ)r k . 1 + T1 (r, θ) k=1 Next we consider the following differential equation X ∞ X∞ dt 1 k = =1+ gk (θ) r(θ, ρ) = 1 + Gk (θ)ρk−1 . dθ 1 + T1 (r(θ, ρ), θ) k=1 k=2 Hence, Z ⊤(ρ) Z 2π Z 2π  X∞  dθ ⊤(ρ) = dt = = 1+ Gk (θ)ρk dθ 0 0 1 + T1 (r(θ, ρ), θ) 0 k=1 X∞ Z 2π = 2π + tk ρk , where tk = Gk (θ)dθ, k=1 0 and we have used that r(ρ, θ) converges uniformly towards 0 when ρ tends to 0. To prove the second item of statement (ii) we will use Lemma 2.1 with n = 2. It is not restrictive to assume that H(x, y) = x2 /2 + y 2/2 + O(3), because otherwise a linear change plus a rescaling of the time can be done before starting the study. Hence there is a relation between ρ > 0 and h = z 2 > 0, given by ρ2 1 √  √  F (ρ, z) = H(ρ, 0)−h = H(ρ, 0)−z 2 = −z 2 +O(3) = ρ+ 2z ρ− 2z +O(3) = 0. 2 2 2 Moreover, when ρ > 0 and h = z > 0 are such that F (ρ, z) = 0 then T (h) = ⊤(ρ) because both values give the period of the same periodic orbit. By Lemma 2.1, near (0, 0),√the above equation has two analytic solutions ρ = Sj (z), j = 1, 2 where Sj (z) = (−1)j 2z + O(z 2 ) satisfy locally F (Sj (z), z) ≡ 0. We are interested √  in S := S2 because it sends positive values into positive ones. Hence ρ = S h , with S analytic at 0 and Taylor series computable simply by implicit derivation, as explained in the proof of Lemma 2.1. √ Finally T (h) = ⊤(ρ) = ⊤(S( h)), as we wanted √ to prove. Notice that although from this result it simply seems that T is analytic on h, from item (i) we already know that all the odd terms of the Taylor series of ⊤ ◦ S at zero must cancell.  Remark 2.2. Assume that the analytic system we are interested in writes as  ′ u = −αv + f1 (u, v), y ′ = βu + g1 (u, v), where the prime denotes the derivative respect some time, say s, and α and β are both e v) = positive. Moreover, assume that it is Hamiltonian, with Hamiltonian function H(u, βu2 αv2 2 + 2 + O(3). If we introduce the following change of variables and time, x y t u= √ , v= √ and s = √ , (10) β α αβ   it writes as system (8) and has the Hamiltonian function H(x, y) = H e √x √v , α . Then β the period function with respect to time t can be obtained in terms of the energy levels, h, by the method developed √ in Theorem 1.1 and finally the time in the variable s is the previous one divided by αβ. THE PERIOD FUNCTION WITH APPLICATIONS 7 Next subsections are dedicated to apply the above results to three examples. In all the examples, by the sake of shortness, we only present a few terms of the Taylor expansion of T (h), but it is not difficult to obtain much more terms. Obviously, the first terms of A(h) can be obtained from the ones of T (h). 2.1. The whirling pendulum. In this example we calculate the first terms in h of the period function for the whirling pendulum. The motion of a whirling pendulum is considered for instance in [32]. It writes as g u′′ = − sin(u) + ω 2 sin(u) cos(u), u ∈ S1 , ℓ where ℓ is the length of pendulum, u its angle deviation, g is the gravity constant and ω is a constant rotation rate. Introducing a new variable v = −u′ converts this second order equation into the planar analytic Hamiltonian system u′ = −v, (11) y ′ = sin(u) (a − b cos(u)) , where a = g/ℓ > 0, b = ω 2 ≥ 0, and Hamiltonian function 2 e v) = v − a cos(u) + b cos2 (u) + a − b H(u, 2 2  2  2 2 v u b a 4 a b = + (a − b) + − u + − u6 + O(7). 2 2 6 24 720 45 When a − b > 0 we have a non degenerate center at the origin and this is the case that we will consider. Following the notation of Remark 2.2 we apply the change of variables and time (10) with α = 1, β = a − b and we obtain     4 e x x2 y 2 4b − a x H(x, y) = H √ ,v = + + 2 + O(5). a−b 2 2 24(a − b) 4 Then, applying the first step of item (ii) of Theorem 1.1 we obtain that (a − 4b)π 2 (11a2 − 16ab + 176b2 )π 4 ⊤(ρ) = 2π + ρ + ρ 8(a − b)2 1536(a − b)4 (−11072b3 + 173a3 − 2568ab2 − 708a2 b)π 6 + ρ + O(ρ8 ). 368640(a − b)6 Doing the computations detailed in the second step of the same item we arrive to √ √ √ √ 1/2 2(a − 4b) 3/2 2(3a2 − 16ab + 48b2 ) 5/2 S( h) = 2h − h + h 12(a − b)2 160(a − b)4 √ 2(120ab2 + 5a3 − 320b3 − 36a2 b) 7/2 − h + O(h9/2 ). 896(a − b)6 Then √ (a − 4b) 3(3a2 − 16ab + 48b2 ) 2 T (h) = ⊤(S( h)) = 2π + πh + πh 4(a − b)2 128(a − b)4 5(5a3 + 120ab2 − 320b3 − 36a2 b) 3 + πh + O(h4 ), 1024(a − b)6 8 C. A. BUZZI, Y. R. CARVALHO, AND A. GASULL is the period function √ for the √ new time introduced in the change of variables. Finally, we must divide it by αβ = a − b to obtain the actual period function of system (11),    2 2π (a − 4b) h 3(3a2 − 16ab + 48b2 ) h T (h) = √ 1+ + a−b 4 2(a − b)2 64 2(a − b)2  3 ! 5(5a3 + 120ab2 − 320b3 − 36a2 b) h + + O(h4 ). 256 2(a − b)2 We can note if b = 0 then the equation correspond to the simple pendulum. Replacing a = g/ℓ gives s    2  3 ! ℓ 1 ℓh 32 ℓh 52 ℓh T (h) = 2π 1+ 2 + 6 + 8 + O(h4 ). g 2 2g 2 2g 2 2g Notice that the above terms coincide with the first ones of the well-known expression of the period function of the pendulum given by Lagrange s ∞  2  n ℓ X (2n)! ℓh T (h) = 2π n 2 , g n=0 (2 n!) 2g obtained from the expression of this period function in terms of an elliptic integral. 2.2. A quadratic system. As a second application we take the simplest family of Hamil- tonians that are not of the form H(x, y) = F (x) + G(y). More concretely we consider the family of quadratic systems with Hamiltonian function x2 y 2 x3 y3 H(x, y) = + − + axy 2 − b . 2 2 3 3 These systems are studied for instance in [4, 27]. Applying our two steps procedure we obtain that π π ⊤(ρ) = 2π + (9a2 + 5b2 − 6a + 5)ρ2 − (9a2 + 5b2 − 6a + 5)ρ3 6 9 5π  + 189 a4 + 378a2 b2 + 77b4 − 180a3 − 84ab2 + 126a2 + 10b2 − 84a + 77 ρ4 +O(ρ5 ), 288 √ √ √ √ 1/2 2 5 2 3/2 32 2 77 2 5/2 896 3 S( h) = 2h + h + h + h + h + h + O(h7/2 ), 3 9 27 54 243 and finally, √ π 5π T (h) = ⊤(S( h)) =2π + (9a2 + 5b2 − 6a + 5)h + 189a4 + 378a2b2 + 77b4 3 72  − 180a3 − 84ab2 + 126a2 + 10b2 − 84a + 77 h2 + O(h3 ). 2.3. System (4) with ε = 0. The third example deals with the Hamiltonian system with x2 y 2 x3 H(x, y) = + + . 2 2 3 In this case, to prove item (iii) of Theorem 1.3 we need to obtain more terms of the Taylor’s development of T (h). As in the previous subsections we obtain first that 5 5 385 4 385 5 103565 6 85085 7 6551545 8 ⊤(ρ) = 2π+ πρ2 + πρ3 + πρ + πρ + πρ + πρ + πρ +O(ρ8 ), 6 9 288 216 31104 15552 663552 THE PERIOD FUNCTION WITH APPLICATIONS 9 and √ √ 2 5 √ 3/2 32 2 77 √ 5/2 896 3 2431 √ 7/2 ρ2 ( h) = 2h1/2 − h + 2h − h + 2h − h + 2h + O(h8 ). 3 9 27 54 243 486 Then, √ 5 385 2 85085 3 37182145 4 T (h) =⊤(S( h)) = 2π + πh + πh + πh + πh 3 72 3888 373248 1078282205 5 1169936192425 6 + πh + πh + O(h7 ). (12) 2239488 483729408 It is very intriguing the appearance of the primorial function in the coefficients of the Taylor’s series of T (h). Recall that if pn denotes the nth prime number, then the primorial of pn is denoted by pn # and is pn # = p1 p2 · · · pn . For instance recall that the numbers pn # + 1 play a key role in the proof of Arquimedes of the existence of infinitely many prime numbers. Computing some more terms of the expression of T (h) get  5# 11# 17# 23# 29# 5 · 7 · 31# 6 T (h) = π 2 + 2 h + 4 3 h2 + 5 6 h3 + 10 7 h4 + 11 8 h5 + h 2·3 2 3 23 2 3 2 3 214 311  5 · 41# 7 5 · 47# 8 5 · 7 · 53# 9 7 · 11 · 59# 10 + 15 12 h + 22 14 h + h + h + O(h11 ). 2 3 7 2 3 7 223 318 226 319 In fact, it is known that T (h) is a monotonous increasing function, defined for h ∈ [0, 1/6) and tending to infinity when h goes to 1/6, see for instance [17]. As we will prove in Proposition 1.2, it can be seen that the area function A, where ′ A (h) = T (h), satisfies the Hill’s equation (7), 5 A′′ (h) = A(h), A(0) = 0, A′ (0) = 2π. 6(1 − 6h)h From this second order differential equation it is easier to obtain more terms of the Taylor’s series of T at zero. We notice that this result is computationally simpler that our general approach but only works for some special Hamiltonian systems. The above initial value problem can be solved in terms of some hypergeometric function 2 F1 . It holds that 1 5  A(h) = 2πh 2 F1 , ; 2; 6h , 6 6 where recall that for |z| < 1, ∞ X (a)n (b)n z n 2 F1 (a, b; c; z) = , n=0 (c)n n! (d)n = d(d + 1)(d + 2) · · · (d + n − 1) is the Pochhammer symbol and (d)0 = 1. This expression helps to understand the appearance of all prime numbers in the coefficients of T because all prime numbers are of the form 6k + 1 or 6k + 5 and all these factors appear in the numerators of (1/6)n (5/6)n . 3. Definitions and preliminary results. This section reviews some definitions and prove some results that we will be use to prove Proposition 1.2 and Theorem 1.3. 10 C. A. BUZZI, Y. R. CARVALHO, AND A. GASULL 3.1. Limit cycles and Abelian integrals. The second part of the Hilbert’s 16th prob- lem asks about the maximum number of limit cycles and their relative locations in planar polynomial vector fields. It is one of the most famous and difficult open problems in mathematics, see [28, 41]. At the end of the last century there has been a very significant advance when Ilyashenko and Écalle independently proved the Dulac problem which is the case of individual finitude, that is, the number of limit cycles of a given planar poly- nomial differential system is finite. We address for a very particular case of a weaker version of Hilbert’s 16th problem, the so called infinitesimal Hilbert’s problem, which asks about an upper bound for the number of zeros of a particular Abelian integral. Let XH = (−Hy , Hx ) be the planar Hamiltonian vector field associated to (1) and consider a perturbation given by Xε = XH + εY , where Y = (P, Q) with P and Q polynomials. The Poincaré-Pontryagin functions or Melnikov functions of order k ∈ N are obtained from the coefficients of the displacement function of the first return Poincaré map as a Taylor’s series in the ε variable near 0, that is, if Pε is the first return Poincaré map of the planar system Xε then its displacement function is given by ∆ε (h) = Pε (h)−h and it has a Taylor’s series in the ε variable near 0 given by ∆ε (h) = εM1 (h) + ε2 M2 (h) + · · · + εk Mk (h) + O(εk+1), which converges to small values of ε. Thus, when M1 (h) ≡ M2 (h) ≡ Mk−1 (h) ≡ 0, the Poincaré-Pontryagin functions or Melnikov functions of order k ∈ N is given by Mk (h), k ∈ N. The Poincaré-Pontryagin Theorem ensures that Z M1 (h) = Q(x, y)dx − P (x, y)dy, γh and that from each simple root of the M1 bifurcates a single hyperbolic limit cycle. (m−1) ∗ Moreover, if there is an h∗ such that M1 (h∗ ) = M1′ (h∗ ) = . . . = M1 (h ) = 0 and (m) ∗ M1 (h ) 6= 0 we have at most m limit cycles bifurcating from γh . So the total number ∗ of the limit cycles, counting the multiplicities, bifurcating from a bounded continua of periodic orbits is at most the number of isolated zeroes, taking into account their mul- tiplicities, of the Abelian integral M1 (h). This is the way how isolated roots of Abelian integrals are related with the number of limit cycles of perturbed Hamiltonian systems. It is costumary to consider I(h) = −M1 (h) as the first Poincaré-Pontryagin function, that is Z Z ZZ ∂P (x, y) ∂Q(x, y) I(h) = ω= P (x, y)dy − Q(x, y)dx = + dxdy, γh γh Int(γh ) ∂x ∂y where in the last equality we have used Green’s Theorem. Hence in a few words, the number of isolated zeros of I(h), counted with their mul- tiplicities, gives an upper bound for the number limit cycles of Xε generated from the ovals of H near ε = 0. Moreover, if all these zeroes are simple this number of zeroes gives rise to the same number of hyperbolic limit cycles for the perturbed system. For more details see for instance [10, 44]. Notice that again from Green’s Theorem we know that for any i, j ∈ N, Z ZZ Z i j i−1 j i x y dy = ix y dxdy = − xi−1 y j+1dx, γh Int(γh ) γh j + 1 so, considering P and Q polynomials, the function I(h) can always be written as the linear combination ℓ X Z I(h) = βk Jk (h), where Jk (h) = xik y jk dx, (13) k=0 γh THE PERIOD FUNCTION WITH APPLICATIONS 11 for some ℓ ∈ N, where ik , jk ∈ N and all βk depend on the coefficients of P and Q. In view of expression (13) it is natural to study the number of zeros of linear combina- tions of ℓ + 1 functions. If all these functions are linearly independent it is not difficult to find linear combinations with exactly ℓ simple zeros in any given interval. On the other hand, when the Hamiltonian H is polynomial, many times some of the involved func- tions have some linear or functional relations, like for instance, the so called Picard-Fuchs equations that include also the derivatives of Ik . Once all these relations are taken into account there appear some other functions, say Jˆk (h), also involving Abelian integrals, and maybe other elementary functions, such that ℓ p X X I(h) = βk Jk (h) = αk Jˆk (h), k=0 k=0 for some ℓ ≥ p ∈ N. Then, the simplest situation is when these p + 1 functions form a so called extended complete Chebyshev system, see Section 3 for more information. This will be the case for the system considered in this paper. In the literature there are many works dealing with zeros of Abelian integrals, see again [10, 44] and their references. Without the aim of being exhaustive and for completeness we list some other techniques used elsewhere to approach the problem. For example, in some works (see [15, 27, 36]) there is a study the geometrical properties of the so-called centroid curve using the fact that it verifies a Riccati equation (which is itself deduced from a Picard-Fuchs system). On the other hand in [21, 23, 24, 25, 37], the authors use complex analysis and algebraic topology (analytic continuation, argument principle, monodromy, Picard-Lefschetz formula, . . . ). Other times it is proved that the p + 1 functions are a Chebyshev system with some accuracy k, meaning that the maximum number of limit cycles provided by the Abelian integral is p + k, see for instance [19]. To end this subsection we state a simple but useful general result, proved in [13], that will be used in the proof of item (i) of Theorem 1.3. Lemma 3.1. Set L ⊂ R an open real interval and let Fj : L → R, j = 0, 1, . . . , N, be N +1 linearly independent analytic functions. Assume also that one of them, say Fk , 0 ≤ k ≤ N, has constant sign on L. Then, there exist real constants cj , j = 0, 1, . . . , N, such that the P linear combination N j=0 cj Fj has at least N simple zeroes in L. 3.2. Relations among Abelian integrals. This subsection is devoted to find relations among the integrals Z Ik (h) = xk y dx, k = 0, 1, 2, . . . , (14) γh where γh , for h ∈ (0, 1/6), are the ovals of {H(x, y) = x2 /2 + x3 /3 + y 2/2 = h} and their derivatives. In particular we obtain the Picard-Fuchs equations satisfied by I0 and I1 and the Hill’s equation satisfied by I0 . All our computations are rather standard and we do not give all the details, see for instance [10]. Lemma 3.2. Consider the Abelian integrals defined in (14). Then, (i) For all 1 ≤ k ∈ N, Xk   k j k−j 3 2 (3k − j)I3k−j−1(h) = 0. (15) j=0 j In particular, 3I3 + 5I4 9I5 + 21I6 + 16I7 I2 = −I1 , I5 = − , I8 = − . 2 4 12 C. A. BUZZI, Y. R. CARVALHO, AND A. GASULL (ii) For 3 ≤ k ∈ N, (2k + 5)Ik (h) + 3(k + 1)Ik−1 (h) − 6(k − 2)hIk−3 (h) = 0. (16) In particular, −12I2 + 6hI0 −15I3 + 12hI1 I3 = , I4 = . (17) 11 13 Proof. (i) Since 3x2 + 3y 2 + 2x3 = 6h, for any k ≥ 1, we have Z Z X k   Z  2 k 2  3 k k j k−j 0= 6h − 3y dy = 3x + 2x dy = 32 x3k−j dy γh γh j=0 j γh X k   Z X k   k j k−j 3k−j−1 k j k−j =− 3 2 (3k − j) x y dx = − 3 2 (3k − j)I3k−j−1(h), j=0 j γh j=0 j where we have used Green’s Theorem. (ii) From x2 /2 + x3 /3 + y 2 /2 = h we know that on γh , (x + x2 )dx + ydy = 0. Hence, multiplying this equality for xk−2 y, integrating and using again Green’s Theorem, we obtain that Z  Z Z k−1 k  k−2 2 0= x + x y dx + x y dy = Ik−1 (h) + Ik (h) + xk−2 y 2 dy γh γ γh Zh k−2 =Ik−1(h) + Ik (h) − xk−3 y 3 dx 3 γ Zh  k−2 2  =Ik−1(h) + Ik (h) − xk−3 2h − x2 − x3 y dx 3 γh 3 k−2  2  =Ik−1(h) + Ik (h) − 2hIk−3 (h) − Ik−1 (h) − Ik (h) 3 3 2k + 5 k+1 2(k − 2) = Ik (h) + Ik−1 (h) − hIk−3 (h). 9 3 3 From the above equality the result follows.  Next result gives a relation between the Abelian integrals and their derivatives. Lemma 3.3. Considering the Abelian integrals defined in (14). Then: (i) It holds that I0 (h) = −A(h) and I0′ (h) = −T (h), where T is the period function associated to H and A(h) is the areal surrounded by the oval γh . (ii) For 0 ≤ k ∈ N, ′ ′ 2Ik+3 (h) + 3Ik+2 (h) + 3Ik (h) − 6hIk′ (h) = 0. (18) Proof. Recall that the Gelfand-Leray formula, see for instance [29, Thm. 26.32], allows to compute easily the derivative of Abelian integrals under suitable regularity conditions. It asserts that Z Z d ω= η, dh γh γh provided that dω = dH ∧ η. In particular, for 0 ≤ k ∈ N, by taking ω = xk y dx and η = xk /y dx, since dH = (x + x2 )dx + ydy, it holds that dω = dH ∧ η. Hence, Z Z ′ d k xk Ik (h) = x y dx = dx. (19) dh γh γh y THE PERIOD FUNCTION WITH APPLICATIONS 13 (i) By item (i) of Theorem 1.1 we know that A′ (h) = T (h), where A(h) is the area surrounded by γh . Hence, by Green’s Theorem, ZZ Z A(h) = dx dy = − y dx = −I0 (h). Int(γh ) γh Therefore, Z ′ 1 T (h) = A (h) = −I0′ (h) =− dx, γh y where we have used (19) for k = 0. In fact, for this particular Hamiltonian, the above relation simply follows by using the first differential equation of the Hamiltonian system, dx/dt = −y. (ii) By using (19) and the expression of H(x, y) = h it holds that Z ′  xk ′ ′ 6hIk (h) = 2x3 + 3x2 + 3y 2 dx = 2Ik+3 (h) + 3Ik+2 (h) + 3Ik (h), γh y as desired.  Proposition 3.4. Considering the Abelian integrals defined in (14). Then (i) The functions I0 and I1 verify the 2 × 2 Picard-Fuchs equations  ′     I0 (h) 6(5h − 1) −7 I0 (h) 6(6h − 1)h = . (20) I1′ (h) 6h 42h I1 (h) (ii) The function I0 verifies the Hill’s equation 5 I0′′ (h) = I0 (h). 6(1 − 6h)h (iii) The function p(h) = I0 (h)/I1 (h) satisfies the Riccati differential equation 1  p′ (h) = 7p2 (h) + 6(2h + 1)p(h) + 6h . 6(6h − 1)h (iv) It holds that I2 (h) = −I1 (h) and, for all 1 ≤ n ∈ N, I3n (h) =hpn−1 (h)I0 (h) + qn−1 (h)I1 (h), I3n+1 (h) =hrn−1 (h)I0 (h) + sn (h)I1 (h), I3n+2 (h) =hun−1 (h)I0 (h) + vn (h)I1 (h), where w ∈ {p, q, r, s, u, v}, and wk ∈ Q(h) denotes a polynomial of degree k. More- over, I3n+2 is a linear combination of several Ij with j < 3n + 2 and j 6≡ 2 (mod 3). Proof. (i) By using that I2 = −I1 , together with (17), expression (18) for k = 0 and k = 1, write as 6hI0′ (h) + I1′ (h) − 5I0 (h) = 0, 6hI0′ (h) + 6(2 − 11h)I1′ (h) + 6I0 (h) + 77I1 (h) = 0. Solving this system with respect to I0′ and I1′ gives the result. (ii) From the first equation of (20) we get that 6(5h − 1)I0 (h) − 6(6h − 1)hI0′ (h) I1 (h) = . (21) 7 By plugging this expression in the second equation of (20) we arrive to the Hill’s equation. (iii) It is a direct consequence of item (i). 14 C. A. BUZZI, Y. R. CARVALHO, AND A. GASULL (iv) The result follows by using induction on n and equalities (16) taking in each step blocks of three integrals I3n , I3n+1 and I3n+2 . For instance, 6 12 I3 (h) = hI0 (h) + I1 (h), 11 11 90  12 180  I4 (h) = − hI0 (h) + h− I1 (h), 143 13 143 108  30 216  I5 (h) = hI0 (h) + − h + I1 (h). 143 13 143 The final property for I3n+2 is a simple consequence of relation (15).  Proof of Proposition 1.2. Notice that this proposition simply restates Proposition 3.4 but changing I0 (h) and I1 (h), by the two functions A(h) and T (h). In fact, it suffices to use the results of Lemma 3.3, I0 (h) = −A(h) and I0′ (h) = −T (h), and that the expression (21) reads as 6 6 I1 (h) = (1 − 5h)A(h) + (6h − 1)hT (h). (22) 7 7 Then all the results are simple computations.  3.3. Involutions and rational parameterizations. Let A be a smooth function with a minimum at x = 0. Then, it has associated an involution σ defined on some open interval K = (xl , xr ) ∋ 0 that satisfies A(x) = A(σ(x)). Recall that a map σ is called an involution if σ ◦ σ = Id and σ 6= Id . By the results of [26] (see Theorem 3.12) this involution plays an important role when studying some Abelian integrals associated to the Hamiltonian H(x, y) = A(x) + y 2/2. In our case A(x) = x2 /2 + x3 /3 and hence z = σ(x) is defined implicitly by x2 x3 z 2 z 3 x + z x2 + xz + z 2 + − + = (x − z)S(x, y) = 0, where S(x, y) = + . (23) 2 3 2 3 2 3 Solving S(x, z) = 0 we get 1 p  z = Z ± (x) = −3 − 2x ± 3(3 − 4x − 4x2 ) . (24) 4 Then σ = Z + and K = (−1, 1/2). As we will see, when one wants to apply Theorem 3.12 we need to control the sign of functions of the form R(x, σ(x)), where R ∈ R(x, y) is a polynomial. In this situation it is very useful to introduce the so called rational parameterizations of algebraic curves. Given a planar algebraic curve R(x, y) = 0, it is said that admits a rational parameteri- zation if there exist two non-constant rational functions u(s) and v(s), s ∈ R, such that R(u(s), v(s)) ≡ 0. Cayley-Riemann’s Theorem ([1, 2]) ensures that R can be rationally parameterized if and only if its genus is zero. Moreover, in such case there are effective methods to find a parameterization, see for instance [40, Chap. 4&5]. In particular, when S is an irreducible quadratic polynomial, it has genus 0 and it can be rationally param- eterized. Next lemma gives one of its parameterizations and other useful expressions for our forthcoming computations. Its proof is straightforward. Lemma 3.5. A rational parameterization of the algebraic curve S(x, z) = 0 with S given in (23) is −3s(s − 2) −3s x = u(s) = 2 , z = v(s) = 2 , 2(s − 2s + 4) s − 2s + 4 and u and v map bijectively [0, 1] into [0, 1/2] and [−1, 0], respectively. Moreover 3(s − 2) σ(u(s)) = Z + (u(s)) = v(s) and Z − (u(s)) = s2 − 2s + 4 THE PERIOD FUNCTION WITH APPLICATIONS 15 Let us illustrate how to use it and its advantages with respect other approaches in a simple example that will be used later. Assume that we want to prove that the function M(x) = 1 + x + σ(x), x ∈ (0, 1/2) (25) does not vanish. A first naive way consists in trying to find its solutions, plugging the expression of σ = Z + , given in (24). Then, isolating the square root term and squaring in √ both sides we obtain the polynomial equation 2x2 + 2x − 1 = 0 that has the root x0 = ( 3 − 1)/2 in (0, 1/2), that in fact is not a solution of (25). So, this approach fails unless we discard this spurious solution. A second powerful approach consists on using resultants, see [43]. An advantage is that it can be utilized for any involution z = σ(x) defined implicitly by a polynomial relation S(x, z) = 0. This is the method used systematically in [26, 33]. In this case it reduces to prove that the following resultant  U(x) = Resz M (x, z), S(x, z) , where M(x, z) = 1 + x + z, does not vanish for x ∈ (0, 1/2). Of course, for this simple M there is no need of doing the resultant because M (x, z) = 0 is equivalent to z = −1 − x, but for higher degree functions this general approach can always be used. In this case U(x) = (2x2 +2x−1)/6, and as in the previous approach U(x0 ) = 0 and, as a consequence, we cannot assure that M(x) does not vanish on (0, 1/2). Finally, with the approach that we propose, we can prove our goal. The only disad- vantage is that it only works when S(x, z) = 0 has genus 0, but fortunately, this is the situation in the case we are dealing with. Notice that to prove that M(x) does not vanish on (0, 1/2) it suffices to prove that   s2 + 4s − 8 M u(s) = M u(s), v(s) = 1 + u(s) + v(s) = − 2 2(s − 2s + 4) does not vanish for s ∈ (0, 1), result that trivially holds. In fact, the reason why this third approach works in this case, while the two previous ones do not, is simple. The first two approaches consider simultaneous the other branch z = Z − (x) defined by S(x, z) = 0 and this branch is not taken into account in the third one. In fact,  s2 − 8s + 4 M u(s), Z − (u(s) = − 2 , 2(s − 2s + 4) √ and this function vanishes at s = s0 = 4 − 2 3 ∈ (0, 1) and u(s0 ) = x0 . The reader interested to see more utilities of the rational parameterizations in dynam- ical systems can take a look to [20]. We will use either the second or the third methods when we study the sign of functions R(x, σ(x)). In fact, both approaches lead to a final polynomial in one variable in Q(x), x ∈ (0, 1/2) or Q(s), s ∈ (0, 1). It is well known that the control of the zeroes of these polynomials in the respective intervals can be done by computing their Sturm sequences, see for instance [42]. We briefly recall this method that we will systematically use without giving the details. Definition 3.6 (Sturm’s sequence). A sequence (f0 , . . . , fm ) of continuous real functions on [a, b] is called a Sturm’s sequence for f = f0 on [a, b] if the following holds: (a) f0 is differentiable on [a, b]. (b) fm does not vanish on [a, b]. (c) if f (x0 ) = 0, x0 ∈ [a, b], then f1 (x0 )f0′ (x0 ) > 0. 16 C. A. BUZZI, Y. R. CARVALHO, AND A. GASULL (d) if fi (x0 ) = 0, x0 ∈ [a, b], then fi+1 (x0 )fi−1 (x0 ) < 0, i ∈ {1, . . . , m}. Theorem 3.7 (Sturm’s Theorem). Let (f0 , . . . , fm ) be a Sturm’s sequence for f = f0 on [a, b] with f (a)f (b) 6= 0. Then the number of roots of f in (a, b) is equal to V (a) − V (b), where V (c) is the number of changes of sign in the ordered sequence (f0 (c), . . . , fm (c)), where zeroes are not taken into account. A Sturm’s sequence for any polynomial f with simple roots can be easily found by a small variation of Euclid’s algorithm for finding the greatest common divisor, see again the classical book [42]. 3.4. Chebyshev systems. For the characterization of Chebyshev systems in an open interval we will use the following results which can be found in [30] and [34]. Definition 3.8. Let u0 , . . . , un−1, un be functions defined in an open interval L of R. (a) The set of functions (ui )ni=0 form a Chebyshev system, or for short T -system, on L if any nontrivial linear combination a0 u0 + · · · + an un has at most n isolated roots in L. (b) The ordered set of functions (ui )ni=0 form a complete Chebyshev system, or for short a CT -system, on L if (ui )ki=0 form a T -system for all k = 0, 1, . . . , n. (c) The ordered set of functions (ui )ni=0 form an extended complet Chebyshev system, or for short an ECT -system, on L if any nontrivial linear combination a0 u0 + · · · + ak uk has at most k isolated roots in L counting multiplicity, for every k = 0, 1, . . . , n. Notice that an ECT -system on L is also a CT -system on L. Definition 3.9. Let u0 , . . . , un functions that have derivatives until order n on L. The Wronskian of such functions in x ∈ L is given by u0 (x) ··· un (x) u′0 (x) ··· u′n (x) W (u0, . . . , un )(x) = .. .. .. . . . . (n) (n) u0 (x) · · · un (x) The following result is the most common approach to prove that a set of functions forms an ECT -system. Lemma 3.10. The ordered set of functions (u0 , . . . , un ) forms an ECT -system on L if, and only if, for every k = 0, . . . , n, W (u0, . . . , uk )(x) 6= 0 for every x ∈ L. Pn Remark 3.11. If (J0 , J1 , . . . , Jn ) forms Pn an ECT -system on L then i=0 αi Ji = 0 has the i same roots bifurcation diagram that i=0 βi t = 0 for the P simple ECT -system (1, t, . . . , tn ). In particular, the coefficients αi can be chosen such that ni=0 αi Ji = 0 has n simple roots in L. Next result was developed by Grau, Mañosas and Villadelprat ([26, 33]) and is an extension of a previous work of Li and Zhang [31] where the authors provided a sufficient condition for the monotonicity of the ratio of two Abelian integrals. It allows to prove that a set of Abelian integrals, of some special shape and for a special type of Hamiltonian system, form a Chebyshev system, simply proving that a similar property is satisfied by the integrands. Next, we state a version of Theorem B in [26] adapted to our interests. Recall that when A has a minimum in x = 0 then the origin of the Hamiltonian systems has a center. Moreover, A has an associated involution σ such that A(x) = A(σ(x)) for x ∈ (xl , xr ) ∋ 0. THE PERIOD FUNCTION WITH APPLICATIONS 17 Theorem 3.12. ([26]) Let us consider the n Abelian integrals Z Jk (h) = fk (x)y 2s−1 dx, 0 < s ∈ N, k = 0, . . . , n − 1, γh where each fk (x) is an analytic function and, for each h ∈ (0, h0 ), γh is the oval sur- rounding the origin contained in the level set γh = {A(x) + y 2/2 = h}. Let σ be the involution associated to A, and define fk (x) fk (σ(x)) ℓk (x) = ′ − . A (x) A′ (σ(x)) If (ℓ0 , . . . , ℓn−1 ) is a CT -system on (0, xr ) and s > 2(n − 2) then (J0 , . . . , Jn−1 ) is an ECT -system on (0, h0 ). When condition s > 2(n − 2) is not fulfilled it is possible, in some situations, to obtain equivalent expressions of the Abelian integrals for which the corresponding new “s” is large enough to verify the inequality, see next lemma. The procedure for obtaining these new Abelian integrals follows from Lemma 4.1 of [26] and other tricks developed in that paper. For the sake of completeness we also include its proof. Lemma 3.13. Let γh be an oval inside the level set {A(x) + y 2 /2 = h}. (i) If F is a smooth function such that F/A′ is analytic at x = 0, then for s ∈ N ∪ {0}, Z Z  ′ s F (x) F (x)y dx = ′ (x) y s+2dx. γh γh (s + 2)A (ii) If F is a smooth function such that F ·A/A′ is analytic at x = 0, then for s ∈ N∪{0}, Z Z  ′  s F (x)A(x) F (x) h F (x)y dx = ′ + y s+2dx. γh γh (s + 2)A (x) 2 Proof. (i) Notice that Z   Z  ′ Z F (x) s+2 F (x) s+2 F (x) s+1 0= d y = y dx + + y dy γh (s + 2)A′ (x) γh (s + 2)A′ (x) γh A′ (x) Z  ′ Z F (x) s+2 = ′ (x) y dx − F (x)y s dx, γh (s + 2)A γh where in the last equality we have used that A′ (x) dx + y dy = 0 on γh . (ii) In this case, Z Z   Z Z s y2 s s F (x) s+2 h F (x)y dx = A(x) + F (x)y dx = A(x)F (x)y dx + y dx γh γh 2 γh γh 2 Z  ′  F (x)A(x) F (x) = ′ (x) + y s+2dx, γh (s + 2)A 2 where in the last step we have used item (i).  4. Proof of Theorem 1.3 and Corollary 1.4 Proof of Theorem 1.3. (i) We start listing the set of ordered Abelian integrals I0 (h), I1 (h), I3 (h), I4 (h), I6 (h), . . . , I3k−2 (h), I3k (h), I3k+1(h), I3k+3 (h), . . . , Im (h), where we have removed from the list of all functions Ij , 0 ≤ j ≤ n, the ones with j ≡ 2 (mod 3) and m = n unless n ≡ 2 (mod 3), case in which m = n − 1. For the sake of notation we will denote the above functions by F0 (h), F1 (h), F2 (h), . . . , FN (n) (h), 18 C. A. BUZZI, Y. R. CARVALHO, AND A. GASULL keeping the same order. Notice that then, N (n) X I(h) = cj Fj (h), (26) j=0 where cj can be taken as arbitrary constants. We claim that these N(n) + 1 functions are linearly independent. Moreover, since F0 (h) = I0 (h) = −A(h), the function F0 does not vanish for h ∈ (0, 1/6). Hence, we can apply Lemma 3.1 and item (i) follows. Let us prove the claim. First of all, notice that I(h) can be written as I(h) = v[N (n)/2] (h)I0 (h) + w[(N (n)−1)/2] (h)I1 (h), where vk and wk are arbitrary polynomials of degree k, because, by using item (iv) of Proposition 3.4, we know that each time that for a given m we consider some new terms F2m (h) and F2m+1 (h) in I(h) it is equivalent to the appearance of two new terms of the form hm I0 (h) and hm I1 (h) in the expression of I(h). Notice that for each n, I(h) is expressed as a linear combination N(n) + 1 functions of the form hj I0 (h) and hk I1 (h), for suitable 0 ≤ j ≤ [N(n)/2] and 0 ≤ k ≤ [(N(n) − 1)/2]. In order to prove the claim, assume that we consider a linear combination of them that gives identically zero. Then I1 (h) v[N (n)/2] (h) v(h) ≡− =: , I0 (h) w[(N (n)−1)/2] (h) w(h) or, in other words, we had that I1 /I0 would be the rational function v/w. On the other hand we know that at A(h) ∼ A0 and T (h) ∼ k ln(1 − 6h) when h ↑ 1/6, for some 0 < k ∈ R. This is so for the area function A(h) because A0 is the area surrounded by the homoclinic loop contained in {x2 /2 + x3 /3 + y 2 /2 = 1/6} and for the period function T (h), because we know that limh↑1/6 T (h) = ∞ and its dominant asymptotic term is given by the passage time near the hyperbolic saddle (−1, 0) of the system (3) with ε = 0, see for instance [22]. By using (22) and that I0 (h) = −A(h) we get that     T (h) 1 7 I1 (h) 1 7 v(h) = 5h − 1 − = 5h − 1 − A(h) (6h − 1)h 6 I0 (h) (6h − 1)h 6 w(h) would be a rational function. This is a contradiction unless v = w = 0 because when h ↑ 1/6 the left-hand side of the above equality goes to infinity as k ln(1 − 6h) and the right-hand only can go to infinity with speed c(1 − 6h)−m for some 0 < m ∈ N and c ∈ R. (ii) From Proposition 1.2 we know that all the Abelian integrals Ij can be expressed in terms of polynomials of h, T (h) and A(h). Similarly in Proposition 3.4 we get a similar property but changing A(h) and T (h) by I0 (h) = −A(h) and I1 (h). Moreover A′ (h) = T (h) and I1 (h) can be obtained from I0 (h) and I0′ (h), see equation (21). In any case, given the Taylor series at h = 0 of any of the following functions T (h), A(h) or I0 (h) the Taylor series at h = 0 of all the other function Ij (h), 0 ≤ j can be easily obtained by using the results of these propositions and, equivalently the Taylor series of all the Fj (h), also at h = 0. For convenience we introduce the function Gj (h) = Fj (h)/(πh) for THE PERIOD FUNCTION WITH APPLICATIONS 19 all j ≥ 0 and the the expression of (26) for h ∈ (0, 1/6) can be written as N (n) I(h) X = dj Gj (h), πh j=0 for arbitrary real dj . For instance, from the expression of T (h) given in (12) given in Subsection 2.3 or the one of I0 (h) we obtain that 5 385 2 85085 3 7436429 4 1078282205 5 G0 (h) = − 2 − h − h − h − h − h + O(h6 ), 6 216 15552 373248 13436928 35 5005 3 323323 4 185910725 5 4775249765 6 G1 (h) =h + h2 + h + h + h + h + O(h7), 18 864 15552 2239488 13436928 5 385 3 17017 4 7436429 5 770201575 6 G2 (h) = h2 + h + h + h + h + O(h7 ), 3 72 864 93312 2239488 35 5005 4 2263261 5 26558675 6 G3 (h) = − h2 − h3 − h − h − h + O(h7 ), 8 288 31104 82944 5 77 85085 5 7436429 6 770201575 7 G4 (h) = − h3 − h4 − h − h − h + O(h8 ), 4 8 1728 31104 663552 and similarly we have obtained all Gj (h) for 0 ≤ j ≤ N(50) = 33 until order 50. Of N (50) course, we do not explicite them. To prove that the function (Gj (h))j=0 are an ECT system in a neighborhood (0, h1 ) of h = 0, from Lemma 3.10 it suffices to prove that the following Wronskians G0 (h) ··· Gk (h) G′0 (h) ··· G′k (h) Wk (h) = Wk (G0 , . . . , Gk )(h) = .. .. .. , . . . (k) (k) G0 (h) · · · Gk (h) do not vanish at h = 0, for k = 0, 1, . . . , N(50). After some tedious calculations we get that Wk (0) 6= 0 for all these values of k. For instance, W0 (0) = −2, W1 (0) = −2, W2 (0) = −20/3, W3 (0) = 140/3, W4 (0) = −12320/3, 11211200 83859776000 2899871054080000 W5 (0) = − , W6 (0) = − , W7 (0) = , 9 9 9 and so on. The result on the G′j s implies the one for the Fj′ s and, as a consequence, the desired result for the Ij , j 6≡ 2 (mod 3). (iii) We will apply Theorem 3.12 when A(x) = x2 /2 + x3 /3 to the Abelian integrals  I0 (h) and I1 (h). By this result it suffices to prove that (ℓ0 , ℓ1 ) is an ECT -system in 0, 12 , where 1 1 1 1 ℓ0 (x) = − and ℓ1 (x) = − . x(1 + x) σ(x)(1 + σ(x)) 1 + x 1 + σ(x) Derivating implicitly S(x, σ(x)) = 0, where S is given in (23), we get that 4x + 2z + 3 σ ′ (x) = − , 2x + 4z + 3 where z = σ(x). Moreover their Wronskians are (1 + x + z)(z − x) W0 (ℓ0 )(x) = ℓ0 (x) = , x(1 + x)z(1 + z) (z − x)3 (4x2 + 6xz + 4z 2 + 7x + 7z + 3) W1 (ℓ0 , ℓ1 )(x) = . x2 (1 + x)2 z 2 (1 + z)2 (2x + 4z + 3) 20 C. A. BUZZI, Y. R. CARVALHO, AND A. GASULL We will prove that all factors do not vanish when x ∈ (0, 1/2). For the first function this holds trivially for all factors but one, 1 + x + z = 1 + x + σ(x) which is precisely the one that we have studied in detail in Section 3.3 with our approach using rational parameterizations, see equation (25). Let us study the remaining factors R2 (x, z) = 4x2 + 6xz + 4z 2 + 7x + 7z + 3 and R1 (x, y) = 2x + 4z + 3. For them it suffices to use the resultants approach, also explained in Section 3.3. It holds that  4 8 2 2 1 Resz R2 (x, z), S(x, z) = x4 + x3 − x2 − x + , 9 9 9 3 2 and  Resz R1 (x, z), S(x, z) = (2x + 3)(2x − 1), both not vanishing for x ∈ (0, 1/2), as desired. In fact, for the first one, Sturm’s approach proves that it has no real roots. (ii) We want to use the same approach that in item (i) to prove that the functions I0 , I1 and I3 form a Chebyshev system, but we must reorder them because, otherwise our approach fails. We will prove that (I3 (h), I1 (h), I0 (h)) form a Chebyshev system on (0, 1/6), or equivalently that the functions (hI3 (h), hI1 (h), hI0 (h)) form an ECT. By item (ii) of Lemma 3.13 we have that Z Z 11x2 + 22x + 12 hI0 (h) =h y dx = f0 (x)y 3 dx, with f0 (x) = , γh γh 18(1 + x)2 Z Z x(13x2 + 27x + 15) hI1 (h) =h xy dx = f1 (x)y 3 dx, with f1 (x) = , γh γh 18(1 + x)2 Z Z 3 x3 (17x2 + 37x + 21) hI3 (h) =h x y dx = f3 (x)y 3 dx, with f3 (x) = . γh γh 18(1 + x)2 To apply Theorem 3.12, first we consider the functions fi (x) fi (σ(x)) ℓi (x) = ′ − ′ , i = 0, 1, 3. A (x) A (σ(x)) and compute the following Wronskians, where z = σ(x), (x − z)R6 (x, z) W0 (ℓ3 )(x) = ℓ3 (x) = , 18(1 + x)3 (1 + z)3 (x − z)3 R7 (x, z) W1 (ℓ3 , ℓ1 )(x) = , 324(2x + 4z + 3)(1 + x)5 (1 + z)5 (x − z)6 R12 (x, z) W2 (ℓ3 , ℓ1 , ℓ0 )(x) = 2916(2x + 4z + 3)3 x3 z 3 (1 + x)7 (1 + z)7 where R6 (x, z) =17x3 z 3 + 51x3 z 2 + 51x2 z 3 + 51x3 z + 141x2 z 2 + 51xz 3 + 17x3 + 128x2 z + 128xz 2 + 17z 3 + 37x2 + 100xz + 37z 2 + 21x + 21z, and the polynomials Rk ∈ Z[x, y] have degree k, and we do not explicit them for the sake of shortness. In the computations of W2 we have used that 12(4x2 + 4xz + 4z 2 + 6x + 6z + 3) σ ′′ (x) = − , (2x + 4z + 3)3 THE PERIOD FUNCTION WITH APPLICATIONS 21 expression obtained once more, derivating implicitly S(x, σ(x)) = 0. Finally, 9s2 (s − 4)2 S8 (s) R6 (u(s), σ(u(s))) = R6 (u(s), v(s)) = − , 8(s2 − 2s + 4)6 27(s − 4)2 S12 (s) R7 (u(s), σ(u(s))) = R7 (u(s), v(s)) = , 8(s2 − 2s + 4)6 729(s − 4)2 S22 (s) R12 (u(s), σ(u(s))) = R12 (u(s), v(s)) = , 32(s2 − 2s + 4)12 where S8 (s) = 5s8 − 55s7 + 197s6 − 43s5 − 1162s4 + 1100s3 + 5240s2 − 10240s + 5120, and Sk ∈ Z[s], have degree k, and we do not explicit S12 and S22 . By using Sturm’s approach we prove that none of them vanish in (0, 1), as desired.  Proof of Corollary 1.4. From the results of Subsection 3.1 the simple zeroes in (0, 1/6) of Z ZZ Z X n ′ ′ I(h) = P (x) dy = P (x) dxdy = − P (x)y dx = αj Ij (h) γh Int(γh γh j=0 give rise to limit cycles of (4) for ε small enough. By using item (i) of Theorem 1.3 the result follows.  Remark 4.1. It is known that for many Hamiltonian systems, some k × k Picard-Fuchs system differential equations is satisfied by several Abelian integrals, see for instance R [35] and their references. When one of these integrals is the area function A(h) = − γh y dx, from its Taylor’s expansion at h = 0 it is not difficult to get the Taylor expansion at h = 0 of all the other Abelian integrals involved in the system. Then the approach used to prove item (ii) of Theorem 1.3, that recall reduces to compute a Wronskian at h = 0, can be applied to study lower bounds of the number of limit cycles bifurcating from the periodic orbits of the Hamiltonian. Since A′ (h) = T (h), the Taylor’s expansion of A(h) near the center can be obtained either by using Theorem 1.1 or by using the linear k-th order differential equation satisfied by A(h) obtained from the Picard-Fuchs system. Notice that to study system (6) we have used both approaches. Acknowlegements This work has received funding from the Ministerio de Economı́a, Industria y Com- petitividad - Agencia Estatal de Investigación (MTM2016-77278-P FEDER grant), the Agència de Gestió d’Ajuts Universitaris i de Recerca (2017 SGR 1617 grant), CAPES grant 88881.068462/2014-01, CNPq grant 304798/2019-3, and São Paulo Paulo Research Foundation (FAPESP) grants 2019/10269-3, 2018/05098-2 and 2016/00242-2. References [1] Abhyankar, S. S. What is the difference between a parabola and a hyperbola? Math. Intelligencer 10, 4 (1988), 36–43. [2] Abhyankar, S. S., and Bajaj, C. L. Automatic parameterization of rational curves and surfaces. III. Algebraic plane curves. Comput. Aided Geom. Design 5, 4 (1988), 309–321. [3] Andronov, A. A., Leontovich, E. A., Gordon, I. I., and Maı̆er, A. G. Theory of bifur- cations of dynamic systems on a plane. Halsted Press [A division of John Wiley & Sons], New York-Toronto, Ont.; Israel Program for Scientific Translations, Jerusalem-London, 1973. Translated from the Russian. [4] Artés, J. C., and Llibre, J. Quadratic Hamiltonian vector fields. J. Differential Equations 107, 1 (1994), 80–95. 22 C. A. BUZZI, Y. R. CARVALHO, AND A. GASULL [5] Beléndez, A., Arribas, E., Márquez, A., Ortuño, M., and Gallego, S. Approximate expressions for the period of a simple pendulum using Taylor series expansion. European J. Phys. 32, 5 (2011), 1303–1310. [6] Beukers, F., and Cushman, R. Zeeman’s monotonicity conjecture. J. Differential Equations 143, 1 (1998), 191–200. [7] Chicone, C. The monotonicity of the period function for planar Hamiltonian vector fields. J. Differential Equations 69, 3 (1987), 310–321. [8] Chicone, C. Geometric methods for two-point nonlinear boundary value problems. J. Differential Equations 72, 2 (1988), 360–407. [9] Chow, S. N., and Hale, J. K. Methods of bifurcation theory, vol. 251 of Grundlehren der Mathe- matischen Wissenschaften [Fundamental Principles of Mathematical Science]. Springer-Verlag, New York-Berlin, 1982. [10] Christopher, C., and Li, C. Limit cycles of differential equations. Advanced Courses in Mathe- matics. CRM Barcelona. Birkhäuser Verlag, Basel, 2007. [11] Cima, A., Gasull, A., and Mañosa, V. Studying discrete dynamical systems through differential equations. J. Differential Equations 244, 3 (2008), 630–648. [12] Cima, A., Gasull, A., and Mañosa, V. Non-integrability of measure preserving maps via Lie symmetries. J. Differential Equations 259, 10 (2015), 5115–5136. [13] Coll, B., Gasull, A., and Prohens, R. Bifurcation of limit cycles from two families of centers. Dyn. Contin. Discrete Impuls. Syst. Ser. A Math. Anal. 12, 2 (2005), 275–287. [14] Constantin, A., and Villari, G. Particle trajectories in linear water waves. J. Math. Fluid Mech. 10, 1 (2008), 1–18. [15] Dumortier, F., Li, C., and Zhang, Z. Unfolding of a quadratic integrable system with two centers and two unbounded heteroclinic loops. J. Differential Equations 139, 1 (1997), 146–193. [16] Foschi, S., Mingari Scarpello, G., and Ritelli, D. Higher order approximation of the period- energy function for single degree of freedom Hamiltonian systems. Meccanica 39, 4 (2004), 357–368. [17] Freire, E., Gasull, A., and Guillamon, A. First derivative of the period function with appli- cations. J. Differential Equations 204, 1 (2004), 139–162. [18] Gasull, A., Guillamon, A., and Mañosa, V. An explicit expression of the first Liapunov and period constants with applications. J. Math. Anal. Appl. 211, 1 (1997), 190–212. [19] Gasull, A., Lázaro, J. T., and Torregrosa, J. On the Chebyshev property for a new family of functions. J. Math. Anal. Appl. 387, 2 (2012), 631–644. [20] Gasull, A., Lázaro, J. T., and Torregrosa, J. Rational parameterizations approach for solving equations in some dynamical systems problems. Qual. Theory Dyn. Syst. 18, 2 (2019), 583– 602. [21] Gasull, A., Li, W., Llibre, J., and Zhang, Z. Chebyshev property of complete elliptic integrals and its application to Abelian integrals. Pacific J. Math. 202, 2 (2002), 341–361. [22] Gasull, A., Mañosa, V., and Villadelprat, J. On the period of the limit cycles appearing in one-parameter bifurcations. J. Differential Equations 213, 2 (2005), 255–288. [23] Gautier, S., Gavrilov, L., and Iliev, I. D. Perturbations of quadratic centers of genus one. Discrete Contin. Dyn. Syst. 25, 2 (2009), 511–535. [24] Gavrilov, L. The infinitesimal 16th Hilbert problem in the quadratic case. Invent. Math. 143, 3 (2001), 449–497. [25] Gavrilov, L., and Iliev, I. D. Bifurcations of limit cycles from infinity in quadratic systems. Canad. J. Math. 54, 5 (2002), 1038–1064. [26] Grau, M., Mañosas, F., and Villadelprat, J. A Chebyshev criterion for abelian integrals. Trans. Amer. Math. Soc. 363, 1 (2011), 109–129. [27] Horozov, E., and Iliev, I. D. On the number of limit cycles in perturbations of quadratic Hamiltonian systems. Proc. London Math. Soc. (3) 69, 1 (1994), 198–224. [28] Ilyashenko, Y. Centennial history of Hilbert’s 16th problem. Bull. Amer. Math. Soc. (N.S.) 39, 3 (2002), 301–354. [29] Ilyashenko, Y., and Yakovenko, S. Lectures on analytic differential equations, vol. 86 of Grad- uate Studies in Mathematics. American Mathematical Society, Providence, RI, 2008. [30] Karlin, S., and Studden, W. J. Tchebycheff systems: With applications in analysis and sta- tistics. Pure and Applied Mathematics, Vol. XV. Interscience Publishers John Wiley & Sons, New York-London-Sydney, 1966. [31] Li, C., and Zhang, Z.-F. A criterion for determining the monotonicity of the ratio of two abelian integrals. journal of differential equations 124, 2 (1996), 407–424. THE PERIOD FUNCTION WITH APPLICATIONS 23 [32] Lichardová, H. Limit cycles in the equation of whirling pendulum with autonomous perturbation. Appl. Math. 44, 4 (1999), 271–288. [33] Mañosas, F., and Villadelprat, J. Area-preserving normalizations for centers of planar Hamil- tonian systems. J. Differential Equations 179, 2 (2002), 625–646. [34] Marděsić, P. Chebyshev systems and the versal unfolding of the cusps of order n, vol. 57 of Travaux en Cours [Works in Progress]. Hermann, Paris, 1998. [35] Novikov, D., and Yakovenko, S. Redundant Picard-Fuchs system for abelian integrals. J. Dif- ferential Equations 177, 2 (2001), 267–306. [36] Peng, L. P. Unfolding of a quadratic integrable system with a homoclinic loop. Acta Math. Sin. (Engl. Ser.) 18, 4 (2002), 737–754. [37] Petrov, G. S. The Chebyshev property of elliptic integrals. Funktsional. Anal. i Prilozhen. 22, 1 (1988), 83–84. [38] Plociniczak, L., and Świtala, M. Monotonicity, oscillations and stability of a solution to a nonlinear equation modelling the capillary rise. Phys. D 362 (2018), 1–8. [39] Rothe, F. The periods of the Volterra-Lotka system. J. Reine Angew. Math. 355 (1985), 129–138. [40] Sendra, J. R., Winkler, F., and Pérez-Dı́az, S. Rational algebraic curves, vol. 22 of Algorithms and Computation in Mathematics. Springer, Berlin, 2008. A computer algebra approach. [41] Smale, S. Mathematical problems for the next century. Math. Intelligencer 20, 2 (1998), 7–15. [42] Stoer, J., and Bulirsch, R. Introduction to numerical analysis, third ed., vol. 12 of Texts in Applied Mathematics. Springer-Verlag, New York, 2002. Translated from the German by R. Bartels, W. Gautschi and C. Witzgall. [43] Sturmfels, B. Solving systems of polynomial equations, vol. 97 of CBMS Regional Conference Se- ries in Mathematics. Published for the Conference Board of the Mathematical Sciences, Washington, DC; by the American Mathematical Society, Providence, RI, 2002. [44] Uribe, M., and Movasati, H. Limit cycles, Abelian integral and Hilbert’s sixteenth problem. Publicações Matemáticas do IMPA. [IMPA Mathematical Publications]. Instituto Nacional de Matemática Pura e Aplicada (IMPA), Rio de Janeiro, 2017. 31o Colóquio Brasileiro de Matemática. [45] Waldvogel, J. The period in the Lotka-Volterra system is monotonic. J. Math. Anal. Appl. 114, 1 (1986), 178–184. Mathematics Department, Universidade Estadual Paulista Julio de Mesquita Filho, 15054-000 São José do Rio Preto, São Paulo, Brazil E-mail address: claudio.buzzi@unesp.br Mathematics Department, Universidade Estadual Paulista Julio de Mesquita Filho, 15054-000 São José do Rio Preto, São Paulo, Brazil E-mail address: yagor.carvalho@unesp.br Mathematics Department, Universitat Autònoma de Barcelona, 08193 Bellaterra, Barcelona, Catalonia, Spain E-mail address: gasull@mat.uab.cat