Comparison on Sufficient Conditions for the Stability of Hill Equation : An Arnold ’ s Tongues Approach

It is known that the solutions of a second order linear differential equation with periodic coefficients are almost always analytically impossible to obtain and in order to study its properties we often require a computational approach. In this paper we compare graphically, using the Arnold Tongues, some sufficient criteria for the stability of periodic differential equations. We also present a brief explanation on how the authors, of each criterion, obtained them. And a comparison between four sufficient stability criteria and the stability zones found by perturbation methods is presented.


Introduction
The second order differential equations are often encountered in engineering and physical problems, and they have been studied for more than a hundred years.The Hill equation is a particular and representative equation among the linear periodic equations, and it receives its name after the work of W. Hill on the lunar perigee.
The Hill equation can be used to describe from the simples dynamical systems to the more complex systems; from a child playing on a swing or a spring mass system [1] to the suspension bridges [2] or the flow of the light in a 1D photonic crystal [3].The complexities of studying the properties of periodic systems often demand a computational approach.However we are not so much interested in the exact shape of the solutions but in the question whether the solution is stable or unstable.By stable solution we mean one that is bounded in the whole real line.
For any linear T-periodic ODE, the knowledge of the state transition matrix in one period; [ ] 0, t T ∈ , gives us sufficient information to know the solution in any t ∈  .The main problem with the determination of the stability of a periodic sys- tem solutions lies in the calculation of the state transition matrix, which is almost always impossible to obtain analytically.Hence it is desirable to find some criteria for determining the stability of solutions without the need of obtain them.Zhukovskii [4] gave one of the most known sufficient condition for stability, namely; let the differential equation then the solutions of the periodic system are stable.
Lyapunov in his celebrated work "The general problem of the stability of motion" [5] developed an approximation of the discriminant 1 of the periodic system (see section 2) and obtained a variety of sufficient stability conditions, being the most known the one here presented.Authors such as Borg [4] Yakubovich [6] [7] and Tang [8] had taken advantage of the properties of the Hamiltonian structure in order to obtain sufficient conditions for stability.Hochstadt [9] [10] and Xu [11] had made use of the properties of the Sturm-Liouville equation and by means of successive approximations, respectively, in order to develop a new discriminant approximation and then finding sufficient criteria for stability.
The aim of the present work is to collect and graphically display some of the most known and efficient sufficient criteria for the stability of the periodic differential equation ( ) ) ( ) p t T p t + = , known as Hill's equa- tion.There is a vast number of stability criteria, in the literature, they are based on different approaches, here we present an easy explanation of four of those approaches starting with: a) the approximation of the discriminant due to Lyapunov followed by; b) canonical forms (Hamiltonian structure); c) Sturm-Liouville equation properties and ending with d) another discriminant approximation due to Shi [12].Each stability criterion will be used in order to find zones in the x and 2 x be two linearly independent solutions of a periodic differential equation subject to the initial conditions ( ) ( ) , then the discriminant is equal to ( ) ( ) x T x T +  , T is the minimum period of the differential equation.

Preliminaries
Consider the linear differential equation with periodic coefficients where ( ) ( ) p t T p t + = is a periodic piece-wise continuous real function with zero average, i.e.
where ( ) q t is a T periodic real function and λ is a real constant.Throughout the document we will use system (1) or (2) indistinctly.
We say that a system is stable if and only if all its solutions are bounded in the whole real line.It is known that for some values of the parameters α and β (or λ ) the Equation (1) has bounded solutions and for some others the solu- tions grow without bounds.The plane of parameters α β − can be splited into stable zones (where all solutions are stable) and unstable zones (where at least one solution is unstable).Stable and unstable regions are separated by some curves known as transition curves, these transition curves are defined by having at least one periodic or anti-periodic solution [14]. where B t B t T × = + ∈ .Notice that (3) may also be written as the Hamiltonian system ( ) where ( ) H t is the symmetric and T periodic matrix . The Equation ( 4) may be taken as a definition of Hamiltonian system and it will be used in section 4.

Let
( ) 0 , t t Φ be the state transition matrix of the system (3), from Floquet Theorem it is known that the matrix ( ) 0 , t t Φ may be written as a multiplication of three matrices, two of them are time dependent matrices and one constant matrix; one of the time independent matrices is bounded and periodic, and the other is an exponential one, the former gives us information about the phase of the solutions and the latter contains information about the growth of the solutions, see [16] or [14].Such factorization is where ( ) ( ) so, the growth of the system solutions are related with monodromy matrix M .
For better understanding notice that for all 0 t ≥ , t kT τ = + , 0, ) T τ ∈ and for some positive integer k any solution of (3) can be written as Let ( ) f t and ( ) g t be real solutions of (2) subject to the initial conditions ( ) ( ) ( ) ( ) Then the state transition matrix associated to ( 2) is and the characteristic multipliers 3 are solutions of the characteristic equation where the independent term is equal to one because of the Liouville theorem 4   [19] and the function is known as the characteristic constant or the discriminant associated to (2) and it plays a fundamental roll on the determination of the stability of the system, the Trace operator in (11) is the sum of the main diagonal entries of a square matrix.By Floquet theory, the condition for the solutions, ( ) f t and ( ) g t , to be stable may be restated as: If ( ) A λ > one solution is stable and one unstable; if ( ) A λ = then one solu- tion is periodic when ( ) 1  A λ = , (or anti-periodic when ( ) A λ = − ); and the second solution may or may not be periodic (anti-periodic).
The Haupt oscillation Theorem asserts that the λ real line can be split into alternating intervals known as stability and instability intervals, the former are characterized by ( ) A λ < and the latter by ( ) A λ > , the endpoints of the 2 We say that a matrix B is symplectic if the condition B JB J ′ = is fulfilled.Symplectic matrices are of even order and they form a group [18] 3 The monodromy matrix eigenvalues i µ are usually called characteristic multipliers or just multip- liers.4 The fact that the independent term be equal to one follows from the symplectic matrix property that states that the determinant of a symplectic matrix is equal one.intervals ( ) 1 A λ = are characterized by values of λ for which the system (2)   has at least one periodic solution.The Haupt oscillation Theorem [20] is as follows Theorem 2.2.For the differential equation ( )
Notice that the Theorem 2.2 establishes conditions, in terms of λ , for the solutions of (2) to be stable; Theorem 2.2 is easily restated so that the stability conditions depend on the parameters α and β of the system (1).It is well known that the stability of solutions of any Hill equation, of the form (1), can be represented as a stability chart in the plane of parameters α β − , unstable zones are called Arnold tongues.Figure 1 shows the stable and unstable zones, in the α β − plane, of the Mathieu, Meissner and Lyapunov equations.

Stability Criteria Based on Lyapunov Approximation
In [5] Lyapunov proposed a method to approximate the discriminant ( ) p t ≥ , and in [23] he delved into the study of the approximation properties.In order to be consistent with Theorem 2.2 we will consider ( ) ( ) where each coefficient is defined as a definite n-multiple integral, that is ( )  ( )( ) ( ) ( ) ( ) notice that the sub-index of each coefficient n A is equal to the order of the n-multiple integral, that is, the coefficient 3 A requires a triple definite integral to be calculated, 4 A requires a forth order definite integral to be calculated, and so on.
In [23] Lyapunov studied in depth the properties of the approximation for ( ) A λ such as the convergence of the series and the monotonic decrease of the series coefficients.One can prove that for For details of the inequality (20) then the solutions of ( ) 0 As far as knowledge of the authors, the only reference in which this criterion is graphically shown is [ [24], pp.205] Before proceeding, we shall call the interval ( ) 0 ,λ −∞ the zeroth instability zone, the intervals ( ) , λ λ the first and the second instability zone and so on.Similarly the intervals ( ) ( ) ( ) , , , , , zeroth, first and second stability zones respectively.
Next criterion can be proved by using the above described Lyapunov method [14]. and where . Then, ( ) p t belongs to the zero stability domain.The blue and red areas are the stability zones found by criterion 0.3 and criterion 0.5 respectively.From the Lyapunov 2 criterion, one can see that the parameter a must fulfill the inequality 1 0 2 a < < .The stability area found by the Criterion 0.5 is obtained by calculating the stability area defined by the equation (23) and 500 different values of a, equispaced in the interval [ [ 0,0.5 , and then, merging the 500 resulting areas.
From Figure 2 one can notice that the blue and red regions are the same for the three equations.This is because the integrals of ( ) ( ) ( ) are all equal to zero, so both criteria depends only on the real constants α and a .
For more criterion based on Lyapunov method see [14] where a whole chapter is dedicated to study stability of the characteristic constant ( )

Stability Criteria Based on Properties of Canonical Forms
Consider a second order system in canonical (Hamiltonian) form ( ) where and J is the skew symmetric matrix defined in sec- tion 2.
The state transition matrix of (24) may be expressed as where ( ) ( ) F t is a real matrix function, and K is a real matrix with ( ) 0 Tr K = , see appendix B. It can be proved that matrix K could be defined as 5( ) where the sign ± is chosen so that K be real.From the factorization (25) we can notice that the stability of the canonical system (24) depends on the exponential matrix KT e and therefore on the matrix K; it can be proved that the matrix K is similar to one of the three matrices shown in Table 1.The set of all possible matrices K could be divided into subsets depending on the determinant sign of each matrix K, following the nomenclature of [14], we will say that: a) (24) has one stable solution and one unstable solution (the characteristic multipliers are real and distinct from 1 ± ); if K ∈ Π Table 1.Relation between the subsets  , Π and  and the stability of canonical system solutions (24),  , λ and φ are real positive constants.
K is similar to K ∈ µ ( ) x t and ( ) 1, then (24) has one periodic solution ( ) x t and one un- stable solution ( ) tx t , with 0 =  both solutions are stable and periodic (the characteristic multipliers are 1 + or 1 − ), these solutions corresponds to coexistence 6 ; and, if K ∈ then both solutions are bounded (the characteristic multipliers are complex, lie on the unit circle and are distinct from 1 ± ).All of these properties are summarized in the Table 1.
Remark 4.1.The subset Π defines the transition boundaries, i.e., Π defines lines on the α β − plane separating the stable zones from the unstable ones.Moreover if K ∈ Π and 0 ≠  then (24) has one periodic stable solution ( ) x t and one unstable solution of the form ( ) ( ) ( ) Let Ω be the set of real continuous function matrices ( ) Hamiltonian system (24), where a is a non-zero arbitrary vector and let x ϕ de- note the rotation of the solution x in time T, since ( ) Ω denote the set of matrices ( ) F t such that the rotation over a period T is π x n ϕ = , and


, each of the subsets n Ω are disjoint sets [14].
By the above mentioned properties, of the subsets Π ,  ,  and n Ω , we can say that the symmetric matrix ( ) H t in (24) belongs to one of the subsets , be two linear independent solutions of (24), and [ ] so the area of the parallelogram defined by x and y , do not change do to the influence of ( ) Remember that every linear Hamiltonian system (24) may be defined as where ( ) the matrix H is the symmetric matrix associated to (24) and ( ) ( ) ( ) is a solution of the same equation.Defining

( )
x t ω as the argument of ( ) it follows from the multiplication [ ] [ ] . Integrating (30) we get then the Hamiltonian system (24) belongs to the n-th stability region n  ( 0,1,2,3, n =  ).One must notice that Hill's Equation ( 1) could be written as in (24), setting thus Hill's equation could be seen as a Hamiltonian system.
It is easy to verify that (34) may be written as ( ) where c is a positive real constant.Then the smallest and largest eigenvalues of ( ) Notice that the Yakubovich 1 criterion inequalities could be reformulated as follow.For simplicity consider the Meissner equation, the integral over one period of the functions ( ) where the function ( ) p t is the excitation function of each periodic differential equation, see section 1.
Figure 3 shows the stability zones obtained by applying Yakubovich 1 q t c ≡ is ( ) Π .Suppose that for some function ( ) in words, this condition establishes that the distance from 2 n c ∈ to the func- tion ( ) q t is less than the distance from 2 c to the boundary of n  , so the function ( ) n q t ∈ .Thus inequality (44) is a test for the stability of ( 24).This test requires the explicit expressions for sufficient conditions for the stability of the Hamiltonian system (24).Now consider the differential equation where ( ) q t is a non-negative, non-identically equal to zero, piecewise conti- nuous periodic function with period T.This and the following results were obtained by V. A. Yakubovich [7].
Remark 4.5.Notice that in these criterion, the regions which are guaranteed stable are not convex.
For Meissner equation, the criteria 0.11 and 0.12 may be rewritten just as we did in the criterion 0.9.Table 2 shows the sufficient conditions for stability of

Meissner equation solutions
From Table 2 one can notice that the n-th stability zone ( on the parameter c.The Yakubovich 2 criterion assumptions, 1 2 2 , define parallelograms whose vertexes depends on the current value of c, for example, for the stability zone 1  the vertexes are: ( ) ( ) . And the inequality From the Figure 4 we notice that the stability zone obtained by the Yakubovich 2 criterion belongs to the cone ( )    ( ) By doing a similar procedure for Yakubovich 3 criterion we obtain the Figure 5.
Notice that the cone ( ) the zones defined by the Yakubovich 3 criterion.This is not the same for the cases of Mathieu and Lyapunov equations, see Figure 6.
Next criterion was developed by Borg [4].An alternative proof was made in [14].
Criterion 4.8 (Borg).Consider the system (45) and let ( ) Suppose that for some integer n ( ) Then all solutions of Equation ( 45) are bounded on ( ) , −∞ +∞ and the corresponding Hamiltonian in (24) is in the stability domain n  .
Figure 6 shows the stable zones obtained by the Yakubovich 2 (red), Yakubovich It is worth to notice that Borg criterion uses almost the same statements of the criteria 0.11 and 0.12 but av q is written instead of 2 c and the distance be- tween the function ( ) q t and the constant c is divided by 2. Moreover for Mathieu, Meissner and Lyapunov equations the constant av q is equal to α , so Borg criterion may be rewritten as: The solutions of Mathieu, Meissner and Lyapunov equations belongs to the stability domain n  if for some integer n the inequalities ( ) for Mathieu, Meissner and Lyapunov equations are equal to 4 2 , 2π and 5.40537 respectively.

Stability Criterion Based on Properties of the Sturm Liouville Equation
In [13] Hochstadt proved that given the differential equation notice that ϕ could be seen as a simile of the rotation of the solutions, see sec- tion 4.
These class of systems are called reversible by V. I. Arnold, see [26].
C. A. Franco, J. Collado DOI: 10.4236/am.2017.8101091499 Applied Mathematics We know that for the solutions of (55) to be stable the inequality ( ) ( ) ( ) x t and ( ) x t are linearly independent solutions of (55) then ( ) ( ) x t T + are also solutions and can be written as x t T x t x T x t x T x t T x t x T x t x T x

T t = −
and noticing that the solutions ( ) x t and ( ) x t are even and odd respectively, that is, ( ) ( ) x t x t = − and ( ) ( ) x t x t = − − .Solving the equations for ( ) x T , ( ) x T and ( ) x T  we arrive to (57), see [20].
From (57) one can easily rewrite the stability inequality ( ) ( ) The stability of the solution is then determined by the examination of the , and the number of zeros of ( ) are positive and ( ) x t , ( ) x t have no zeros in 0, 2 and 2 x belongs to the first stability zone [27].So, if ( ) then the solutions of (55) are bounded.On the other hand, by Equation (56), if ( ) . In [9] it is proved that the solutions "rotation" is bounded by and we get the following Criterion 5.1 (Hochstadt 1).A sufficient condition for the boundedness of all solutions of the periodic differential equation It is well known that, for ( ) ( ) , the transition curves are defined by points in the α β − plane for which there is at least one periodic (anti-periodic) solution of (55).In [13] Hochstadt shows that for periodic solutions of (55) the condition for some positive integer n, must be satisfied.And for anti-periodic solutions the condition ( ) ( ) must be satisfied.
If the solutions of (55) are stable then, the condition ( ) ( ) must be satisfied.By noticing the above Hochstadt generalized criterion 13 as follows [10] Criterion 5.2 (Hochstadt 2).A sufficient condition for the boundedness of all solutions of the periodic differential equation where w λ = , ( ) ( ) ( ) . Then the solution of (68), with suppose that are solutions of (68) and they are subject to the initial conditions ( ) ( ) and one can obtain ( ) 2 A λ by means of successive ap- proximations.For details see [12].
From (78) it is not so hard to prove the next Criterion 6.1 (Xu).Suppose where ( ) Xu criterion, as the Hochstadt 1 and Hochstadt 2 criteria, needs the derivative of the function ( ) q t .In order to avoid troubles we will do the same as in the previous section, i.e. take the first 20 elements of the expansion of ( ) ( ) sin sign t and then substitute them instead of ( ) ( ) Applying Xu criterion to Mathieu, Meissner and Lyapunov equations we obtain the Figure 8.
There is a vast number of sufficient conditions for the stability of periodic differential equations, we have just mentioned and explained some of them.For more stability criteria see for example: [4] where Starzhinskii not only collect sufficient stability criteria for Equation (1) but he consider second order periodic differential equations with dissipation, n-th order systems and some particular cases of the vector equation (Hamiltonian) systems.

Conclusions
We have reviewed some of the most known stability criteria for second order differential equations; these criteria were obtained by four different approaches; Perturbation methods, such as strained parameters (Lindstedt-Poincare method) and multiple scales methods, are frequently used for analysing the stability of the periodic differential equations.They are based on the assumption that the variable-coefficient terms are small in some sense.The stability boundaries associated to a Hill equation may be determined by the strained parameters method, that is, assuming that 1 β  , and then seeking the value of α such that the solutions be T or 2T periodic.The general solution and the coefficient α are written in terms of powers of β (perturbation expansion)   .The accuracy of the method depends on how many elements of the series (80) and (81) are obtained.Notice that the same procedure must be done for each transition curve, see [29].
In [30]  Even though the stability areas found with strained parameters method are larger than the ones obtained with the stability criteria, for the Mathieu equation case (82), the complexities of the method and the set of differential equations that has to be solved, for obtaining the stable zones, make the strained parameters method less suitable than the simplicity of the criteria statements for the stability analysis of periodic differential equations.Besides, we must remember that the strained parameters method has the limitation that it just the solutions are stable, for the best known forms of Hill's equation, Mathieu, Meissner and Lyapunov equations i.e., for equations of the form

Figure 1 .
Figure 1.Stability chart for (a) Mathieu, (b) Meissner and (c) Lyapunov equations, stable areas in white and unstable areas in gray.

Figure 2 Figure 2 .
Figure2shows the stability zones obtained by applying the criteria Lyapunov

Figure 4
Figure 4 Shows the parallelograms defined by different values of c and the triangular sector defined by the inequality

Figure 4 .
Figure 4. Red continuous line represents the boundary of the stability zone obtained by Yakubovich 2 criterion.Parallelograms in blue discontinuous lines are defined by the inequalities 1 1 2 c ≤ < and π 4 cot 2 c c

Figure 5 .
Figure 5. Red continuous line represents the boundary of the stability zone obtained by Yakubovich 3 criterion.Parallelograms in blue discontinuous lines are defined by the inequalities 1 1 2 c ≤ < and
the second term of the right hand side of (72) is

Figure 8 .
Figure 8. Stability zones obtained by Xu criterion for a) Mathieu b) Meissner and c) Lyapunov equations.
series (80) and (81) are substituted into the Hill equation and by grouping terms of like powers of β , one obtains a set of recursive differential equations.The initial condition of α , i.e. 0 α , is the value of α where the Arnold tongues rise, and depends on the excitation function period,

Figure 9
Figure 9(a) shows the transition curves, associated to the first tongue, obtained by the strained parameters method (Lindstedt-Poincare method) and the actual one.Figures 9(b)-(e) show the stability zones obtained by Yakubovich 1, Hochstadt 2, Borg and Xu criteria respectively, in the ranges
∈Ω , and K belongs to  , Π or  respectively, that is, ∈ , i.e. the system (24) is stable and ( ) n F t ∈Ω and K ∈ , then the rotation ϕ of all the solutions of (24) will be t and vectors x , y cannot overlap each other neither π θ ≥ .Then the rotation number of x and y must coincide.Remark 4.3.The definition of the subsets n Π , n  and n  allows us to discriminate between stable (unstable) zones where the solutions rotation have similar properties, see appendix B.

Table 2 .
moreover, by following the same procedure we can prove that each subset of n Yakubovich 2 and Yakubovich 3 criteria for the Meissner equation, where Hochstadt 1 and Hochstadt 2, require the derivative of the func-

Table 3 .
Table4shows us the percentage of the stability zone obtained by using each

Table 3 .
Percentage of stable area obtained by numerical calculation in the plane α β − ,

Table 4 .
Percentage of stability zone obtained by applying each criterion on the Mathieu, Meissner and Lyapunov equation.