This is a preprint of: “Rational parameterizations approach for solving equations in some dynamical systems problems”, Armengol Gasull, José Tomás Lázaro, Joan Torregrosa, Qual. Theory Dyn. Syst., vol. 18, 583–602, 2019. DOI: [10.1007/s12346-018-0300-5] RATIONAL PARAMETERIZATIONS APPROACH FOR SOLVING EQUATIONS IN SOME DYNAMICAL SYSTEMS PROBLEMS ARMENGOL GASULL, J. TOMÁS LÁZARO, AND JOAN TORREGROSA Abstract. We show how the use of rational parameterizations facilitates the study of the number of solutions of many systems of equations involving polynomials and square roots of polynomials. We illustrate the effectiveness of this approach, applying it to sev- eral problems appearing in the study of some dynamical systems. Our examples include Abelian integrals, Melnikov functions and a couple of questions in Celestial Mechanics: the computation of some relative equilibria and the study of some central configurations. 1. Introduction A very common way to solve equations involving polynomials and square roots of polynomials consists in squaring them several times until they become polynomials and look for their zeros using a numerical method. Besides the fact that the degree of such polynomials increases with any squaring, there exists a second inconvenience (harder than the first one): there appear many spurious solutions that must be discarded. This is usually done numerically, that is, finding (with a given precision) all the real roots of the final polynomial equation, replacing them into the original equation and imposing the final value to be smaller than a given tolerance. Of course, this is not a rigorous proof and to get a correct solution, one has to follow these two steps using interval arithmetic (see for instance the first example of this introduction). The problem becomes more complicated when, instead of working with an equation, we deal with systems of equations of similar features. Is this the best convenient procedure to follow?: The answer is “no” and we will illustrate how rational parameterizations can be, in many cases, used to approach this kind of problems. As it will be seen, this method eliminates the spurious solutions and determines a well-posed domain containing the actual ones. This assertion will be supported on several examples. The use of rational parameterizations to study the solutions of equations containing square roots is, in fact, an already known technique (see for example [1, 10, 20] and references therein). Although a few authors have introduced its use in dynamical systems (see, for instance [18]). In our opinion it has not received an special attention and its use in this setting is still infrequent. The aim of this work is to show, by means of some examples, that the use of rational parameterizations is more than an alternative method to find the number of solutions in problems appearing in dynamical systems. Indeed, it permits, in some cases, to reduce the original question to a well-posed polynomial problem. Once this is moved into a polynomial setting, there are standard tools to approach its solutions, like the use of resultants or Gröbner basis (see for instance [23]). The fundamental idea in the method is to impose that any expression appearing inside the square root symbol is equal to some new squared variable. This condition can be rewritten as an algebraic expression of the form Fλ (x, y) = 0, where λ ∈ Rk is a set of 2010 Mathematics Subject Classification. Primary: 34C23. Secondary: 13P15; 14E05, 34C08, 37N05. Key words and phrases. bifurcation, resultant, rational parameterization, Abelian integral, Poincaré– Melnikov–Pontryagin function, relative equilibria, central configuration. 1 2 ARMENGOL GASULL, J. TOMÁS LÁZARO, AND JOAN TORREGROSA parameters. Cayley–Riemann’s Theorem ([2]) ensures that this algebraic curve 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 [21, Chapters 4&5]). In fact, many programs of symbolic calculus have implemented some of these methods and algorithms. The situation where Fλ (x, y) = 0 has genus zero is the one that will appear henceforth in the examples and results presented in this work. Let us begin with a simple well-known toy problem, that of the two leaning ladders (see for instance [5, 13]). It reads as follows: consider two ladders, of length A and B, leaning in opposite directions between two parallel walls. Assume that both of them contact the base of the opposite wall and that they cross in a point at distance h from the ground floor. Can we determine the distance w between both walls? And the distances a and b between the ground floor and the point where each ladder contacts its corresponding wall? (see Figure 1). Applying Thales’ Theorem it follows that the relation 1/h = 1/a + 1/b A B a b h w Figure 1. The two leaning ladders problem must be satisfied or, using Pythagoras, that 1 1 1 =√ +√ , (1) h 2 A −W 2 B −W where W = w2 . Squaring this equation, recurrently, one gets the final polynomial equation of degree 4,   P (W ) := W 4 + −2A2 − 2B 2 + 4h2 W 3 + A4 + 4A2 B 2 − 6A2 h2 + B 4 − 6B 2 h2 W 2  + −2A4 B 2 + 2A4 h2 − 2A2 B 4 + 8A2 B 2 h2 + 2B 4 h2 W + B 4 A4 − 2B 2 h2 A4 + h4 A4 − 2B 4 h2 A2 − 2B 2 h4 A2 + B 4 h4 = 0. Depending on the values of the parameters A, B, and h, one should look for positive zeros of this polynomial. It can be seen that it has more than one of these positive zeros. Afterwards, one must substitute them into equation (1) to decide whether or not they are in fact solutions. An alternative way to solve (1) is to seek for a new variable t, such that W is a rational function of t and in such a way that the equivalent problem has no spurious solutions. To find such change of variables, we notice that from a2 = A2 −W , b2 = B 2 −W the variables (a, b) move on the hyperbola b2 − a2 = B 2 − A2 . This curve, in the (a, b)-plane, has genus zero and so, by Cayley–Riemann’s Theorem, it can be rationally parameterized. Indeed, let us find a rational parameterization for it. Having in mind that b2 − a2 = (b + a)(b − a), RATIONAL PARAMETERIZATIONS FOR SOLVING EQUATIONS 3 we introduce a real parameter t such that B 2 − A2 b + a = t, b−a= . t Isolating a and b in terms of t we get t2 + A2 − B 2 t2 + B 2 − A2 a= , b= . (2) 2t 2t p Since a, b > 0, it follows that t > |A2 − B 2 | and we can parametrise W in terms of t as  2 2 2 2 2 t + A2 − B 2 (A2 − B 2 )2 A2 + B 2 t2 W =A −a =A − =− + − =: ϕ(t). (3) 2t 4t2 2 4 Therefore to study P (W ) = 0 is equivalent to consider P (ϕ(t)) = 0, where 2 t4 − 4ht3 − (A2 − B 2 ) Q(t) P (ϕ(t)) = =: . h(t2 + A2 − B 2 )(t2 − A2 + B 2 ) h(t2 + A2 − B 2 )(t2 − A2 + B 2 ) Hence, instead of solving P (W ) = 0, we can solve Q(t) = 0. Each of the positive roots of this second equation produces, by using (3), a unique W which is a solution of (1). In fact, taking into account Bolzano’s Theorem and Descartes’ rule (it exhibits one change of sign when moving over its√coefficients), p it follows the existence of a unique positive root ∗ ∗ of Q, t = t such that w = W = ϕ(t∗ ) is the unique solution of our initial problem. ∗ The values of a and b are also unique and they can be determined similarly from (2). To illustrate the advantages of this approach, we provide the details in a particular case. So, take for instance A = 7, B = 6, and h = 3. Then, P (W ) = W 4 − 134W 3 + 6163W 2 − 106326W + 426465 = 0 has two positive solutions, W1 ' 5.635917 and W2 ' 33.066465. The actual solution is W1 , since one can discard analytically the spurious solution W2 using simple inequalities. Indeed, by Bolzano’s Theorem, we know that 33 < W2 < 331/10. Hence, from equation (1), it follows that 1 1 1 1 1 1 −√ 2 −√ 2 = −√ −√ h A − W2 B − W2 3 49 − W2 36 − W2 √ 1 1 1 1 3 < −√ −√ = − < 0, 3 49 − 33 36 − 33 12 3 and so W2 is not a solution of (1). It has artificially appeared when squaring the formulas. Notice that in general discarding spurious solutions using interval arithmetic computation can be quite complicated (see for instance the last example in Section 3). On the other hand our approach provides a direct answer: in the new variable t, we get that Q(t) = t4 − 12t3 − 169 = 0 has a unique positive root, t∗1 ' 12.095502, and then   ∗ 169 85 t2 ϕ(t1 ) = − 2 + − = W1 . 4t 2 4 t=t∗ 1 √ Hence, the solution of our problem is w = W1 ' 2.374008. Another easy application of the rational parameterization is to prove identities involving square roots. For instance the one appearing in [17], which says that for x > y > 0 it 4 ARMENGOL GASULL, J. TOMÁS LÁZARO, AND JOAN TORREGROSA holds that s p s p √ x 2 x −y 2 x x2 − y 2 x±y = + ± − , 2 2 2 2 follows straightforwardly from the parameterization t 2 + s2 t2 − s2 x= , y= , 2t 2t with t > s > 0. Indeed, s p √ t √ s x x2 − y 2 t±s x+y = √ , x−y = √ and ± = √ . t t 2 2 2 t This useful parameterization is obtained again from a genus zero planar curve. As before, it is natural to impose that x2 − y 2 = s2 , for some unknown value s. Hence the planar curve Fs (x, y) = x2 −y 2 −s2 , for s 6= 0, is irreducible and has genus zero. Therefore it admits a rational parameterization. The one given above is the simplest one. These two simple examples show the advantages of using rational parameterizations to solve equations involving square roots. Along the paper some more applications of this use are presented. We briefly describe its contents. In Section 2 we recall with a simple example how to deal with 1-parameter families of polynomials systems, a common final situation after applying rational parameterizations. We will focus our attention on three particular types of applications. The first one, in Section 3, deals with polynomial expressions with two square roots of polynomials. These equations arise in the study of the number of zeros for some Abelian integrals (see [15]). The second one, in Section 4, handles with central configuration problems in Celestial Mechanics. Here, square roots appear due to the Euclidean distances among the bodies. Two examples are presented: one concerning critical points of a certain potential (see [6]) and a second one dealing with central configurations (see [7, 22]). And the third type, in Section 5, where we seek for estimates for the number of limit cycles of planar non-autonomous differential equations by studying the zeros of some Melnikov functions. The problem considered was proposed in [9]. 2. Preliminary results: systems of polynomial equations The well-known Bezout’s Theorem asserts that, generically and taking into account the solutions at infinity, the number Qk of complex solutions x ∈ Ck that a polynomial system {Fj (x) = 0}j=1,...,k has is j=1 deg(Fj ). Moreover, this number only changes when some of these degrees vary. This situation is drastically different when one is interested only in its real solutions. In this sense, we focus our attention on 1-parameter families of polynomial systems and we are concerned about knowing when the number of their simple real solutions vary. To this end we introduce the following definition: Definition 2.1. Let Fi (x, λ) be real polynomials in x ∈ Rk , i = 1, . . . , k, and λ ∈ R. An effective bifurcation value is a real value λ∗ such that the number of simple real solutions of the system {Fi (x, λ∗ + ε) = 0}i is different for ε > 0 than for ε < 0, for any small enough ε. Next example illustrates how these effective bifurcation values can be obtained by applying the techniques developed in [23]. In this particular example, both polynomials have degree 2 and effective bifurcations occur. RATIONAL PARAMETERIZATIONS FOR SOLVING EQUATIONS 5 Proposition 2.2. The effective bifurcation values for the number of simple real solutions of the system of polynomial equations F (x, y, λ) = (−λ + 1)x2 + 3λy 2 − x − λy = 0, (4) G(x, y, λ) = λx2 + (2λ + 1)y 2 − 2y + λx = 0, are the four real roots of the polynomial P8 (λ) = 107λ8 − 1190λ7 + 2525λ6 − 3222λ5 + 3307λ4 + 568λ3 − 868λ2 + 160λ − 4, namely, λ1 ' −0.5081, λ2 ' 0.02967, λ3 ' 1.7826, λ4 ' 8.7774. Moreover, the number of simple real solutions in the intervals (−∞, λ1 ), (λ1 , λ2 ), (λ2 , λ3 ), (λ3 , λ4 ), (λ4 , +∞) is 2, 4, 2, 4, 2, respectively (see Figure 2). N (−0.65) = 2 N (−0.48) = 4 N (−0.2) = 4 N (0.01) = 4 N (0.1) = 2 N (1) = 2 N (1.9) = 4 N (2.2) = 4 N (3) = 4 N (1000) = 2 e denotes the number of real solutions of system (4) for λ = λ. Figure 2. N (λ) e Proof. The first step consists on computing the resultants of F and G with respect to x and y. This is based on the fact that Res(p(x), q(x), x) = 0 if and only p(x) and q(x) have a common root (see [23]). Thus, Rx (y) = Res(F, G, x) =y((5λ2 − λ − 1)2 y 3 − 2(5λ2 − λ − 1)(λ2 + 2λ − 2)y 2 + (−2λ4 + 8λ3 + 3λ2 − 6λ + 4)y + λ(λ − 2)(λ2 + 2)), Ry (x) = Res(F, G, y) =x((5λ2 − λ − 1)2 x3 + 2(3λ2 + 2λ + 1)(5λ2 − λ − 1)x2 + (11λ4 + 11λ3 − 4λ2 + 14λ + 1)x + λ(2λ − 5)(λ2 + 2)). bx (y) and R Afterwards, we define R by (x) as the polynomials containing the same factors as Rx (y) and Ry (x) but with simple multiplicity. In our case, there are no multiple factors 6 ARMENGOL GASULL, J. TOMÁS LÁZARO, AND JOAN TORREGROSA in the polynomials Rx (y) and Ry (x), so we have R bx (y) = Rx (y) and R by (x) = Ry (x). Next, following the approach proposed in [11, 12], we consider double resultants (see [16] for a general introduction to this setting). Computing them we get  Rxy (λ) = Res Rbx (y), R bx0 (y), y =λ4 (λ2 + 2)2 (5λ2 − λ − 1)4 (λ − 2)4 P8 (λ),  Ryx (λ) = Res Rby (x), R by0 (x), x =λ4 (λ2 + 2)2 (5λ2 − λ − 1)4 (2λ − 5)4 P8 (λ), where P8 is given in the statement. The real roots of each of these polynomials are the only candidates to provide the effective bifurcation values we are looking for. So, all these values are contained in the set of real roots of the polynomial RF G (λ) = λ(5λ2 − λ − 1)(λ − 2)(2λ − 5)(λ2 + 2)P8 (λ). The real roots of P8 (λ) are λ1 ' −0.5081, λ2 ' 0.02967, λ3 ' 1.7826, and λ4 ' 8.7774, and the ones of 5λ2 − λ − 1 are λ5 ' −0.3582 and λ6 ' 0.5582. Hence, the (ordered) real roots of RF G are λ1 < λ5 < 0 < λ2 < λ6 < λ3 < 2 < 5/2 < λ4 . To determine whether they are or not effective bifurcation values, we study the number of simple real solutions of the system on each one of the ten intervals defined by these nine real roots of RF G . To do this, it suffices to fix a rational value of λ in each interval and compute the number of real solutions for that specific value (see Figure 2). It follows that only four of them are effective bifurcation values, which correspond to the ones given in the statement.  3. Polynomial expressions with two square roots Let us consider the following system of ordinary differential equations ẋ = − y(x − c)(y − d) + εP (x, y), ẏ = x(x − c)(y − d) + εQ(x, y), with c > 0, d > 0, and P and Q polynomials. For ε = 0, this system has a linear center at the origin and a couple of straight lines of equilibrium points (precisely, x = c and y = d). Then, all the circles centered at the origin of radius less than min(c, d) define a continuum of periodic orbits. It is well-known that the limit cycles that bifurcate from them come from the simple zeros of the Poincaré–Melnikov–Pontryagin function Z Q(x, y) dx − P (x, y) dy I(r) = x2 +y 2 =r2 (x − c)(y − d) (see for instance [8]). Using polar coordinates, it admits the equivalent expression Z 2π r(Q(r cos θ, r sin θ) sin θ − P (r cos θ, r sin θ) cos θ) I(r) = dθ. 0 (r cos θ − c)(r sin θ − d) As a direct consequence of Lemma 2.2 in [15], this function can be written as   2 1 V (r2 ) W (r2 ) I(r) = U (r ) + 2 √ +√ , c + d2 − r 2 c2 − r 2 d2 − r2 for some suitable polynomials U, V, and W. Writing z = −r2 , the form of the function (c2 − r2 )(d2 − r2 )I(r) motivates the study of the number of zeros of functions of type √ √ Pn0 (z) + Pn1 (z) z + a + Pn2 (z) z + b, (5) with Pnj (z), arbitrary polynomials of degree nj , for j = 0, 1, 2, and a, b ∈ R. The problem of estimating the maximal number of zeros of these functions was studied in [14]. They were considered as elements of a family which is, in fact, an Extended RATIONAL PARAMETERIZATIONS FOR SOLVING EQUATIONS 7 Complete Chebyshev system. A rather different question is, given a particular function of this type (5), to compute effectively its zeros. The following theorem and example show how rational parameterizations are an efficient way to do it. However, the proof of the theorem does not work when the number of square roots involved is larger than 2, or the degrees of the polynomials inside the square roots are higher than 1. Theorem 3.1. Let us consider polynomials Pnj of degree nj , j = 0, 1, 2, with no common factors. Then, the study of the number of zeros of the function √ √ F (x) := Pn0 (x) + Pn1 (x) x + a + Pn2 (x) x + b in x ≥ −a > −b, is equivalent to the study of the number of positive zeros of the polynomial Q(t), of degree 2m, given by the numerator of F (ϕ(t)), where m := max(2n0 , 2n1 +1, 2n2 + 1) and 2 bt2 (2 + t)2 − a (2 + 2t + t2 ) R(t) ϕ(t) = 2 = . 4 (1 + t) (1 + t)2 Proof. We look for (u, v), with u, v ≥ 0, satisfying x + a = v 2 , x + b = u2 . These equalities imply that they must move on the hyperbola u2 − v 2 = b − a =: c > 0. To find a rational parameterization for this genus zero curve, first we take τ such that u√− v = τ , u + v = c/τ , which√leads to u = (c + τ 2 )/2τ , v = (c − τ 2 )/2τ for τ ∈ (0, c]. The transformation τ = c/(t + 1) provides the rational parameterization for t ∈ [0, +∞): 2 + 2t + t2 √ t(2 + t) √ u= c, v = c. (6) 2(1 + t) 2(1 + t) Using this transformation we obtain that x = v 2 − a = u2 − b = ϕ(t), with the function ϕ given in the statement. Hence p p F (ϕ(t)) = Pn0 (ϕ(t)) + Pn1 (ϕ(t)) ϕ(t) + a + Pn2 (ϕ(t)) ϕ(t) + b     R(t) R(t) t(2 + t) √ = Pn0 + P n1 c (1 + t)2 (1 + t)2 2(1 + t)   R(t) 2 + 2t + t2 √ Q(t) + Pn2 2 c= , (1 + t) 2(1 + t) (1 + t)m where R(t) is the polynomial of degree 4 provided by the statement and Q(t) is a suitable polynomial of degree 2m. Notice that there is a one-to-one correspondence between zeros of F (x) and zeros of Q(t), that is, no false zeros are added. Furthermore, the final domain [0, +∞) of Q(t) is very convenient to apply polynomial standard techniques.  Let us illustrate how our approach works with a particular example of function F (x). Example 3.2. Consider, for x ≥ 0, the following function: √ √ Fα (x) = A(x) + B(x) x + C(x) x + 1, with A, B, C the polynomials A(x) = 38340 + α + 5760 x, B(x) = −69993 − 144298 x, C(x) = 23400 − 143702 x, and α being a real parameter. Then, Fα (x) has no real solutions for α < −61740 and it has exactly one for α ≥ −61740. Proof. In this case a = 0 and b = 1 in the statement of Theorem 3.1. Taking the rational parameterization (6),  2 t(2 + t) x = ϕ(t) = , t ≥ 0, (7) 2(1 + t) 8 ARMENGOL GASULL, J. TOMÁS LÁZARO, AND JOAN TORREGROSA we get that the numerator of Fα (ϕ(t)) is Qα (t) = − 72000 t6 − 429120 t5 − 968044 t4 + (2α − 950056) t3 + (6α − 232009) t2 + (6α + 230454) t + (2α + 123480), defined for t ∈ [0, +∞). To apply Descartes’ rule, we consider the values of α that produce changes in the sign of the coefficients of Qα (t). Namely, they are 115227 232009 −61740, − , , 475028. 3 6 It is straightforward to verify that there is no change in the sign of the coefficients of Qα (t) if α < −61740 and one change otherwise. Thus we can assert that Qα (t) has no real roots for t ∈ (0, +∞) in the first case and has at most one in the second case. Moreover, if α ≥ −61740 then Qα (0) ≥ 0 and since limt→+∞ Qα (t) = −∞ the existence of this zero tα > 0 follows. From expression (7) we obtain the value xα = ϕ(tα ) such that Fα (xα ) = 0.  Finally, let us show with a fixed value of α the advantage of using our approach instead of the classical recurrent squaring process to solve Fα (x) = 0. With this classical approach, the corresponding polynomial equation is 2 Gα (x) = A2 (x) + (x + 1)C 2 (x) − xB 2 (x) − 4(x + 1)A2 (x)C 2 (x) = 0. Consider, for instance, α = 2. Then G2 (x) = 29463035904000000x6 − 597819133749657600x5 + 4274549116862945536x4 − 13150960109956884832x3 + 16873855567372135585x2 − 7555772293001326152x + 851096590977473296. Its roots, z1 , z2 , . . . , z6 , are located approximately at: 0.1676387, 0.5656784, 1.731563, 3.713337, 5.502619, 8.609643. They all are candidates to be zeros of the function F2 (x), but only the first one, z1 , is an actual root. The others can be discarded, for example, using interval arithmetics. Let us do it, for instance, for the second one, z2 ' 0.5656784. Indeed, it satisfies 1/2 < z2 < 3/5 and it holds that √ √ F2 (z) = 38342 + 5760z − (69993 + 144298z) z + (23400 − 143702z) z + 1  r  r 3 1 1 1 1 < 38342 + 5760 − 69993 + 144298 + 23400 − 143702 +1 5 2 2 2 2 √ 48451 √ = 41798 − 71071 2 − 6 < 0. 2 Similar computations can be done for the rest of solutions. However, if we apply the rational parameterization provided by (7), we obtain Q2 (t) = −72000t6 − 429120t5 − 968044t4 − 950052t3 − 231997t2 + 230466t + 123484. This polynomial has exactly one zero in [0, +∞), at the point t2 ' 0.4900103. Substituting this value into (7) one obtains x2 = ϕ(t2 ) ' 0.1676387, an approximation of the unique zero of F2 (x). RATIONAL PARAMETERIZATIONS FOR SOLVING EQUATIONS 9 4. Applications to Celestial Mechanics Distances among bodies and, therefore, square roots appear often in equations modeling many problems in Celestial Mechanics. In this section we show two applications of rational parameterization and polynomial tools in this field: a first one concerning critical points of a certain potential (see [6]) and a second one dealing with central configurations in a system of four charged particles (see [7, 22]). 4.1. Relative equilibria in the Schwarzschild isosceles three body problem. The so-called Schwarzschild isosceles three body problem with masses m1 = m2 = M > 0 and m3 = 1 is introduce in detail in [6]. In particular, in cylindrical coordinates, it has the following reduced Hamiltonian 1 2 2M + 1 2 Hred = P + Pz + U, (8) M R 4M with potential C2 A B 4A1 16B1 U (R, z) = 2 − − 3−√ − 2 , MR R R R2 + 4z 2 (R + 4z 2 )3/2 where A, A1 , B, B1 , and C are positive parameters. This model extends to a three mass point system the classical Schwarzschild model of interaction between two mass points, with potential U (r) = −A/r − B/r3 , where r is the relative distance between them. The solutions where the three mass points are steadily rotating about their common center of mass are called relative equilibria. In particular, for the Hamiltonian system defined by (8) they correspond to the critical points of U, that is ∂U ∂U = 0, = 0. (9) ∂R ∂z In [6] it is studied, joint to other problems, the number of solutions of (9) on z = 0 for positive values of A, A1 , B, and B1 in terms of C. In this case, there is only one bifurcation value, C2 in the next result, that distinguish the cases with 0, 1 or 2 relative equilibria. The following result extends the study of the number of solutions of (9), also in terms of C, for arbitrary non-zero values of A, A1 , B, and B1 . We show that the effective bifurcation set is bigger. In the proof it can be seen that there are values with more than two solutions. It is not restrictive to assume C ≥ 0. Theorem 4.1. For non-vanishing values of A, A1 , B, B1 , and M the C-effective bifurca- tion values for the number of solutions of (9) belong to the set EC = {C1 , C2 , C3 } ∩ R, where s 3M 2 (4AB1 − A1 B)2 p √4 C1 = 4 − , C2 = 4 3M 2 (B + 16B1 )(A + 4A1 ), C3 = 3ABM 2 . 16A1 B1 Moreover, there exist values of the parameters such that EC is, actually, the C-effective bifurcation set. Proof. To eliminate the square roots in system (9) we write R2 + 4z 2 = r2 , for some new positive variable r. The family of planar curves Fr (R, z) = R2 + 4z 2 − r2 = 0 has genus zero and can be rationally parameterized as r(1 − t2 ) rt R= , z= , 1 + t2 1 + t2 10 ARMENGOL GASULL, J. TOMÁS LÁZARO, AND JOAN TORREGROSA for t ∈ R and r > 0. Observe that the origin (r = 0) has been removed since it is a singular point of (9). In the parameters (r, t), equations (9) are equivalent to ∂U P8 (r, t) ∂U Q7 (r, t) = = 0, = = 0, ∂r M r (t − 1)3 (t + 1)3 3 ∂t M r (t − 1)4 (t + 1)4 3 where P8 (r, t) = (Ar2 − 4A1 r2 + 3B − 48B1 )t6 − (Ar2 − 12A1 r2 − 9B − 144B1 )t4  − (Ar2 + 12A1 r2 − 9B + 144B1 )t2 + Ar2 + 4A1 r2 + 3B + 48B1 M + 2C 2 r(t4 − 1),  Q7 (r, t) =t (Ar2 + 3B)t4 − (2Ar2 − 6B)t2 + Ar2 + 3B M + 2C 2 rt(t4 − 1). To solve them we follow the procedure presented in Section 2: we compute the resultants Rr (t) = Res(P8 , Q7 , r) = 24 3M 2 t2 (t2 − 1)6 V8 (t), Rt (r) = Res(P8 , Q7 , t) = 220 36 B 6 M 10 (A1 r2 + 12B1 )4 W2 (r), where V8 (t) =3(4AB1 t4 − A1 Bt4 − 8AB1 t2 − 2A1 Bt2 + 4AB1 − A1 B)2 M 2 + 16A1 B1 C 4 (t4 − 1)2 ,  W2 (r) = (4A1 + 1)r2 + 3(B + 16B1 ) M − 2C 2 r, and take br (t) = t(t2 − 1)V8 (t), R bt (r) = (A1 r2 + 12B1 )W2 (r), R which contain all the terms in the preceding functions, but without multiplicity. Now, br (t), R Rrt = Res(R br0 (t), t) = 2118 38 A4 B 12 A25 17 16 16 4 5 4 1 B1 M C (C − C1 ) (C − C3 ) , 4 bt (r), R Rtr = Res(R b0 (r), r) = −214 33 A4 B 3 (A + 4A1 )(C 4 − C1 )2 (C 4 − C2 ), t 1 1 where C1 , C2 , and C3 are the constants given at the statement. Thus, the product of the factors of Rrt and Rtr with no multiplicity, namely C (C 4 − C1 ) (C 4 − C2 ) (C 4 − C3 ) = 0, divides the positive real line in four different intervals having the same number of real simple solutions of system (9). The case A + 4A1 = 0 becomes C2 = 0 and the set EC does not vary. Finally, taking for instance A = 2, B = 3, A1 = −1/3, B = 1/5, and M = 2 we get 0 < C2 < C3 < C1 . The number of zeros on the intervals that they define is, respectively, 0, 2, 6, and 2. In consequence, for these values of the parameters, EC is the C-effective bifurcation set.  Remark 4.2. (i) From this proof it is not difficult to see that, when all the parameters are positive and one restricts himself to the surface z = 0, the C-effective bifurcation set reduces to {C2 }, recovering the result of [6]. Indeed, C1 is not real since −3M 2 (4AB1 − A1 B)2 /(16A1 B1 ) < 0. The condition z = 0 implies t = 0 in the rational parameterization (since r > 0). This means that Q7 (r, 0) = 0 and P8 (r, 0) = (Ar2 +4A1 r2 +3B +48B1 )M − 2C 2 r. The discriminant of P8 (r, 0) is −3ABM 2 − 48AB1 M 2 − 12A1 BM 2 − 192A1 B1 M 2 + C 4 , whose unique positive real root (there is one real negative and two conjugate pure imaginary) is C2 . RATIONAL PARAMETERIZATIONS FOR SOLVING EQUATIONS 11 (ii) Taking A = 2, B = 3, A1 = −1, B = 4, and M = 2 it can be proved that the C-effective bifurcation set is {C3 }. In fact, for 0 < C < C3 the number of equilibria is 1. Even though on each of the intervals C3 < C < C1 and C > C1 this number is 3, their corresponding relative position changes between aligned and triangular-shaped. (iii) There exist values of the parameters whose C-effective bifurcation set is EC , ex- hibiting different number of equilibria. Indeed, taking A = 5, B = 3, A1 = −1, B1 = 4, and M = 2, the respective number of equilibria is 0, 4, 6, and 4. 4.2. Central Configurations. The equations governing the charged Newton’s planar n-body problem are X mi mj mi q̈i = − (qi − qj ) , i = 1, . . . , n. j6=i |qi − qj |3 Here the units have been suitably taken, |q| denotes the Euclidean norm of q, and mi ∈ R \ {0}. We will assume that the center of mass of the n charged particles is located at the origin. Then, they form a spatial charged central configuration if there exists a constant λ ∈ R such that X mj λqi = − 3 (qj − qi ) , i = 1 . . . , n, (10) j6=i |q j − qi | (see [19] and references therein). Actually, the geometric configuration (up to dilations and rotations) of the particles is such that the position vectors are always parallel to the acceleration ones. For a fixed set of masses, to find all possible central configurations is a very hard problem. Two of the main obstructions to solve system (10) are the large number of unknowns and to deal with equations involving square roots. Usual reductions of this problem come from assuming symmetries and restrictions on the dimension of the spatial positions. One of the most studied reduction is the planar non-collinear. In some cases, the method of rational parameterization transforms this problem into seeking for zeros of polynomials in several variables. This is the approach used by Leandro in [18]. He studies the finiteness of central configurations when the d masses are in a plane and two more stay in the orthogonal line to this plane passing through the origin, where is located the center of masses of the system. The restriction to four masses moving on a plane is known as a kite central configuration. The case with three equal masses has been tackled in [7, 22] and the case of two couples of equal masses in [4]. Both cases, in the more general situation of charged masses, have been considered in [3]. (0, `) (−1, 0) (1, 0) (0, −k) Figure 3. Kite configuration, where ◦ denotes the (charged) mass m. 12 ARMENGOL GASULL, J. TOMÁS LÁZARO, AND JOAN TORREGROSA The aim of this section is to describe the number of different kite central configurations in the first frame, charged or not, using rational parameterizations. Without loss of generality we can fix m1 = m2 = m3 = 1 and m4 = m ∈ R \ {0}. Moreover, we locate them as in Figure 3. That is, with the charged mass at the point with coordinates (0, `) and satisfying k + ` > 0. In this setting, following [7], equation (10) reduces to  −3/2   −3/2  −3 2 2 F1 (k, `, m) = m (k + `) (k + `) − ` + 1 + 2k k + 1 − 1/8 = 0,  −3/2   −3/2  (11) F2 (k, `) = (k + `) (k + `)−3 − k 2 + 1 + 2 ` `2 + 1 − 1/8 = 0. Solutions of (11) correspond to different kite central configurations. The non-charged problem (m > 0) is studied in [7, 22] where the authors show that the effective bifurcation set is formed by only one value. In the charged problem, as it will be seen in the next result, the bifurcation set is increased by two negative values of m. Our proof for the non-charged problem is also shorter. Theorem 4.3. The effective bifurcation set for the number of solutions of system (11) is {−1, −1/4, m∗ } where m∗ ≈ 1.00266054757261 is the closest to 1 root of a suitable irreducible polynomial with integer coefficients and degree 240. Indeed, in the intervals (−∞, −1), (−1, −1/4), (−1/4, m∗ ), and (m∗ , +∞) this number is 3, 4, 5, and 3, respec- tively. Proof. To eliminate the square roots in (11) we impose that k 2 + 1 = u2 and `2 + 1 = v 2 , for some new variables u and v. Again, both planar curves have genus zero and can be parameterized through the rational parameterizations  2   2  s − 1 s2 + 1 t − 1 t2 + 1 (k, u) = , , (`, v) = , , 2s 2s 2t 2t with s, t > 0. Substituting them into equations (11) they become P20 (s, t, m) Fe1 (s, t, m) = = 0, (s2 + 1)3 (t2 + 1)3 (t + s)2 (st − 1)2 s Q20 (s, t) Fe2 (s, t) = 2 = 0, (t + 1) (s + 1)3 (t + s)2 (st − 1)2 t 3 2 where P20 , Q20 are polynomials with integer coefficients of degree 20 in s, t. Moreover P20 has degree 1 on m. Recall that we only need to study the number of intersections of the algebraic curves Fe1 (k, `, m) = 0, Fe1 (k, `) = 0 when k + ` > 0. To do it we use the procedure explained in Section 2. First we compute the resultants Rs (t) = Res(P20 , Q20 , s) = t20 (t2 + 1)51 (3t2 − 1)U100 (t), Rt (s) = Res(P20 , Q20 , t) = s20 (s2 + 1)51 (3s2 − 1)V100 (s). Afterwards we consider the functions (with no multiple roots) bs (t) = t(3t2 − 1)U100 (t), R bt (s) = s(3s2 − 1)V100 (s), R RATIONAL PARAMETERIZATIONS FOR SOLVING EQUATIONS 13 and get the new two resultants  bs (t), R Rst (m) = Res R b0 (t), t =m16 (4m − 1)24 (4m + 1)24 (m − 1)60 (m + 1)64 (m + 2)24 s (249m2 − 162m − 23)2 P18 2 2 (m)P240 (m)P644 (m),  bt (s), R Rts (m) = Res R b0 (s), s =m1076 (4m − 1)7 (4m + 1)7 (m − 1)60 (m + 1)64 (m + 2)78 t (249m2 − 162m − 23)2 P240 (m)P408 2 (m), Pi being polynomials of degree i, irreducibles in Q(t). They are not explicitly written because of their big size. Finally, we consider a polynomial containing all the simple factors appearing in Rst and Rts , that is, RP Q (m) =m(4m − 1)(4m + 1)(m − 1)(m + 1)(m + 2) (249m2 − 162m − 23)P18 (m)P240 (m)P644 (m)P408 (m). This polynomial has 74 simple real roots. We analyze these roots taking rational values of the parameter m just before and after them and studying if an effective change in their behavior undergoes. We select those with such effective change, i.e., the effective bifurcation values. In this problem they are just three: −1, −1/4, m∗ , with m∗ being the real root of P240 which is closest to 1.  5. Zeros of Melnikov functions In this section we start with a similar question as the one posed in Section 3. Indeed, we consider a planar non-autonomous periodic differential equation with a continuum of periodic orbits and we study how many of these periodic solutions remain after small perturbations. Those which remain are isolated (limit cycles) and are closely related to the simple zeros of its corresponding Melnikov function. In [9], the following system of coupled generalized Abel equations is considered: ( ṙ = a(t)r2 + εf (t)rn sm , ṡ = b(t)s2 + εg(t)rp sq , with n + m ≥ 3, p + q ≥ 3, and a, b, f and g, T -periodic trigonometric polynomials. Using a variational approach, the authors derive the corresponding Melnikov function:  Z T  n m f (t)  ρ1 ρ2 n−2 (1 − B(t)ρ )m dt  M (ρ) =  Z0 (1 − A(t)ρ1 ) 2 , (12)  p q T g(t)  ρ1 ρ2 p q−2 dt 0 (1 − A(t)ρ1 ) (1 − B(t)ρ2 ) Rt Rt where ρ = (ρ1 , ρ2 ) ∈ R2 , A(t) = 0 a(s) ds, B(t) = 0 b(s) ds, and A(T ) = B(T ) = 0. Motivated by this problem, it seems interesting to study the number of limit cycles arising from a system of type ( ṙ = (cos t)r2 + ε (F0 + F1 sin t + · · · + Fn sinn t) rs2 , (13) ṡ = (cos t)s2 + ε (G0 + G1 sin t + · · · + Gn sinn t) r2 s, where Fi , Gi are real constants such that Fn Gn 6= 0, ε is an small parameter and n ∈ N. The case n = 2 was already considered in [9]. It is not restrictive (rescaling if necessary) to assume Fn = Gn = 1. Denote F = (F0 , F1 , . . . , Fn−1 ) and G = (G0 , G1 , . . . , Gn−1 ). We know that simple real zeros of (12) in D = {|ρi | < 1, i = 1, 2 and ρ1 ρ2 6= 0} provide limit cycles of (13) for small enough ε. 14 ARMENGOL GASULL, J. TOMÁS LÁZARO, AND JOAN TORREGROSA Several computations, that we omit for the sake of shortness, prove that system M (ρ) = 0 in D \ {ρ1 ρ2 6= 0} is equivalent to ( p PF (ρ1 , ρ2 ) + QF (ρ1 , ρ2 ) 1 − ρ22 = 0, p (14) PG (ρ1 , ρ2 ) + QG (ρ1 , ρ2 ) 1 − ρ21 = 0, for some polynomial systems PF , PG , QF , QG . In [9, Prop. 11], for n = 2 and a suitable choice of the parameters F, G, it is showed, by squaring the latter equations, that the system (14) has four simple solutions. Namely, they fix F = F∗ = (−1/2, −1/10), G = G∗ = (−3/5, 3/10) and get, using Gröbner Basis, that the squared system ( Φ1 (ρ1 , ρ2 ) = PF2∗ (ρ1 , ρ2 ) − Q2F∗ (ρ1 , ρ2 ) (1 − ρ22 ) = 0, Φ2 (ρ1 , ρ2 ) = PG2 ∗ (ρ1 , ρ2 ) − Q2G∗ (ρ1 , ρ2 ) (1 − ρ21 ) = 0. has 10 solutions in D from whose only 4 are actual solutions of (14). Let us show how our method works for n = 2 and with the same choice of the parameters F, G. Introducing the new variables si as 1 − ρ2i = s2i , i = 1, 2 we have once more two genus zero curves for which rational parameterizations exist. Plugging in system (14) the rational parameterizations ρi = 2σi /(1 + σi2 ) and si = (1 − σi2 )/(1 + σi2 ), for σi ∈ (−1, 1), i = 1, 2 we obtain an equivalent polynomial system S, with rational coefficients and degree 6 in σ1 , σ2 . This system is triangularized again by using Gröbner Basis, proving that there exist polynomials f and g, of degree 13, 14, respectively, such that the system S is equivalent to {σ1 = f (σ2 ), g(σ2 ) = 0}. The polynomial g has only 6 zeros in (−1, 1). Using interval arithmetic it is easy to prove that only 4 of them belong to D: for the values σ2 such that g(σ2 ) = 0 one checks if f (σ2 ) ∈ (−1, 1). We again find a clear advantage with respect to the previous method: no spurious solutions appear. The study of the number of limit cycles of (13) for small ε but arbitrary n presents some main difficulties. The most important one is that the complexity increases with the degree n. A second one is how to choose suitable parameters F and G giving rise to the highest number of limit cycles. Using the preceding procedure, the following result provides high lower bounds for small values of n. Theorem 5.1. For n = 1, . . . , 16, there exist values of F and G such that system (13) has n2 limit cycles, for small enough values of ε. Proof. However the result is asserted and proved for n = 1, . . . , 16, we strongly believe that the same procedure is applicable to larger values of n. The real limitation to that is its high computational cost. We describe the method for a general n without providing explicit expressions of the polynomials appearing during the process. First, consider the same change as above: ρi = 2σi /(1 + σi2 ), for i = 1, 2. So, system (14) reads    (σ 2 + 1)2  2 2 2 Λ1 (σ1 , σ2 ) = 0, (σ1 + 1)(σ2 − 1)3   (σ 2 + 1)2  2 1 2 Λ2 (σ1 , σ2 ) = 0, (σ2 + 1)(σ1 − 1)3 where Λ1 (σ1 , σ2 ) and Λ2 (σ1 , σ2 ) are polynomials of degree n + 4. Thus, our problem be- comes to find the number of real solutions of S = {Λ1 (σ1 , σ2 ) = 0, Λ2 (σ1 , σ2 ) = 0} in the open square (−1, 1)2 . Second, we look for values F∗ and G∗ such that the polynomial system S has n solutions, for instance, at the points (−1 + 2j/(n + 1), −1 + 2j/(n + 1)) with j = 1, . . . , n. Using again Gröbner basis, the system S is written as the union of systems of type {σ1 = f (σ2 ), g(σ2 ) = 0} for some polynomials f and g. A detailed study with interval RATIONAL PARAMETERIZATIONS FOR SOLVING EQUATIONS 15 arithmetic shows that S has in fact n2 −n more simple solutions, having n2 simple solutions in total.  We finish this section with a simple example of 1-parameter bifurcation in equation (14) with n = 2, whose number of limit cycles varies with the parameter. Let consider F = (−1/2, −1/3) and G = (−3/4 + 2/3λ, λ). Thus, performing the change of variables ρi = 2σi /(1 + σi2 ), for i = 1, 2, system (14) reads  2 4   3σ1 σ2 − 6σ1 σ25 − 2σ1 σ24 − 9σ12 σ22 + 12σ1 σ23 + 3σ24    −4σ 2 σ2 + 8σ1 σ 2 + 6σ1 σ2 − 9σ 2 + 2σ1 − 4σ2 = 0, 1 2 2  12λσ1 σ2 − 12σ1 σ2 + 6σ1 σ2 − 8σ12 σ22 λ − 48λσ12 σ2 + 24λσ1 σ22 + 6σ14 + 24σ13 σ2  4 5 4 2    −15σ12 σ22 − 8λσ12 + 32λσ1 σ2 − 8λσ22 + 24λσ1 − 12λσ2 − 15σ12 + 3σ22 − 8λ + 3 = 0. Using the procedure described in Section 2, it can be proved that its effective bifurcation values are λ = −3/20 and λ = 3/4. In particular, for small enough ε, if −3/20 < λ < 3/4 the system (13) has four limit cycles while only two if λ < −3/20 or λ > 3/4. Acknowledgements This work has received funding from the Ministerio de Economı́a, Industria y Compet- itividad - Agencia Estatal de Investigación (Grants MTM2015-65715-P and MTM2016- 77278-P (FEDER)), the Agència de Gestió d’Ajuts Universitaris i de Recerca (Grants 2017 SGR 1049 and 2017 SGR 1617). References [1] S. S. Abhyankar. What is the difference between a parabola and a hyperbola? Math. Intelligencer, 10(4):36–43, 1988. [2] S. S. Abhyankar and C. L. Bajaj. Automatic parameterization of rational curves and surfaces. III. Algebraic plane curves. Comput. Aided Geom. Design, 5(4):309–321, 1988. [3] F. Alfaro and E. Pérez-Chavela. Families of continua of central configurations in charged problems. Dyn. Contin. Discrete Impuls. Syst. Ser. A Math. Anal., 9(4):463–475, 2002. [4] M. Álvarez-Ramı́rez and J. Llibre. The symmetric central configurations of the 4-body problem with masses m1 = m2 6= m3 = m4 . Appl. Math. Comput., 219:5996–6001, 2013. [5] H. A. Arnold. The Crossed Ladders. Math. Mag., 29(3):153–154, 1956. [6] J. A. Arredondo, E. Pérez-Chavela, and C. Stoica. Dynamics in the Schwarzschild isosceles three body problem. J. Nonlinear Sci., 24(6):997–1032, 2014. [7] J. Bernat, J. Llibre, and E. Pérez-Chavela. On the planar central configurations of the 4-body problem with three equal masses. Dyn. Contin. Discrete Impuls. Syst. Ser. A Math. Anal., 16(1):1– 13, 2009. [8] C. Christopher and C. Li. Limit cycles of differential equations. Advanced Courses in Mathematics. CRM Barcelona. Birkhäuser Verlag, Basel, 2007. [9] B. Coll, A. Gasull, and R. Prohens. Periodic orbits for perturbed non-autonomous differential equa- tions. Bull. Sci. Math., 136(7):803–819, 2012. [10] F. Cukierman. Notas sobre integrales abelianas. Spanish, Notes of the course Geome- trı́a Algebraica, Escuela de Matemática de América Latina y el Caribe EMALCA, Perú, http://mate.dm.uba.ar/∼fcukier/pdfs/IntegralesAbelianas.pdf, 2008. [11] A. Ferragut, J. D. Garcı́a-Saldaña, and A. Gasull. Detection of special curves via the double resultant. Qual. Theory Dyn. Syst., 16(1):101–117, 2017. [12] J. D. Garcı́a-Saldaña, A. Gasull, and H. Giacomini. Bifurcation values for a family of planar vector fields of degree five. Discrete Contin. Dyn. Syst., 35(2):669–701, 2015. [13] M. Gardner. Mathematical circus. MAA Spectrum. Mathematical Association of America, Wash- ington, DC, 1992. More puzzles, games, paradoxes, and other mathematical entertainments from Scientific American, Revised reprint of the 1981 edition, With a preface by Donald Knuth. [14] A. Gasull, J. T. Lázaro, and J. Torregrosa. On the Chebyshev property for a new family of functions. J. Math. Anal. Appl., 387(2):631–644, 2012. 16 ARMENGOL GASULL, J. TOMÁS LÁZARO, AND JOAN TORREGROSA [15] A. Gasull, J. Tomás Lázaro, and J. Torregrosa. Upper bounds for the number of zeroes for some Abelian integrals. Nonlinear Anal., 75(13):5169–5179, 2012. [16] I. M. Gelfand, M. M. Kapranov, and A. V. Zelevinsky. Discriminants, resultants, and multidimen- sional determinants. Mathematics: Theory & Applications. Birkhäuser Boston, Inc., Boston, MA, 1994. [17] D. J. Jeffrey and A. D. Rich. Computer Algebra Systems: A Practical Guide, chapter Simplifying square roots of square roots by denesting, pages xvi+436. John Wiley & Sons, Chichester, United Kingdom, 1999. [18] E. S. G. Leandro. Finiteness and bifurcations of some symmetrical classes of central configurations. Arch. Ration. Mech. Anal., 167(2):147–177, 2003. [19] R. Moeckel. On central configurations. Math. Z., 205(4):499–517, 1990. [20] J. Schicho. Rational parametrization of surfaces. J. Symbolic Comput., 26(1):1–29, 1998. [21] J. R. Sendra, F. Winkler, and S. Pérez-Dı́az. Rational algebraic curves, volume 22 of Algorithms and Computation in Mathematics. Springer, Berlin, 2008. A computer algebra approach. [22] J. Shi and Z. Xie. Classification of four-body central configurations with three equal masses. J. Math. Anal. Appl., 363(2):512–524, 2010. [23] B. Sturmfels. Solving systems of polynomial equations, volume 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. Departament de Matemàtiques, Universitat Autònoma de Barcelona, 08193 Bellaterra, Barcelona, Spain Email address: gasull@mat.uab.cat Departament de Matemàtiques, Universitat Politècnica de Catalunya, Barcelona, Spain Email address: jose.tomas.lazaro@upc.edu Departament de Matemàtiques, Universitat Autònoma de Barcelona, 08193 Bellaterra, Barcelona, Spain Email address: torre@mat.uab.cat