Stability interchanges in a curved Sitnikov problem

TL;DR

This paper investigates the stability behavior of equilibrium points in a curved Sitnikov problem, revealing stability interchanges as orbital parameters change, with broader implications for the n-body problem in curved spaces.

Contribution

It introduces a general theorem on stability interchanges and applies it to a curved Sitnikov problem, linking stability phenomena to curved n-body dynamics.

Findings

One equilibrium point undergoes stability interchanges as the semi-major axis varies.

A general theorem on stability interchanges is formulated and proved.

The results connect stability behavior to the n-body problem in curved spaces.

Abstract

We consider a curved Sitnikov problem, in which an infinitesimal particle moves on a circle under the gravitational influence of two equal masses in Keplerian motion within a plane perpendicular to that circle. There are two equilibrium points, whose stability we are studying. We show that one of the equilibrium points undergoes stability interchanges as the semi-major axis of the Keplerian ellipses approaches the diameter of that circle. To derive this result, we first formulate and prove a general theorem on stability interchanges, and then we apply it to our model. The motivation for our model resides with the $n$-body problem in spaces of constant curvature.

Equations

$$\ddot{z}=-\frac{2z}{(z^2+r^2(t){)}^{3/2}},$$

$$r(t)=1-\epsilon \cos u(t),$$

$$r(t)=(1-\epsilon \cos t)+O({\epsilon }^2),$$

$$\begin{array}{rl}{r}_{\epsilon }(t;r)=& r\rho (t;\epsilon ),r>0,\\ \rho (t;\epsilon )=& \left(1-\epsilon \cos (u(t))\right)=\left(1-\epsilon \cos (t)\right)+\mathcal{O}({\epsilon }^2),\end{array}$$

$${\mathbf{x}}_1(t)$$

$${\mathbf{x}}_2(t)$$

$${\mathbf{F}}_{\epsilon }(y,t;R,r)$$

$$||\mathbf{x}-{\mathbf{x}}_1||={\left[r_{\epsilon }^2(t;r)+2R^2+2R{r}_{\epsilon }(t;r)\cos t-2y(R+r_{\epsilon }(t;r)\cos t)\right]}^{1/2},$$

$$||\mathbf{x}-{\mathbf{x}}_2||={\left[r_{\epsilon }^2(t;r)+2R^2-2R{r}_{\epsilon }(t;r)\cos t-2y(R-r_{\epsilon }(t;r)\cos t)\right]}^{1/2}.$$

$${\mathbf{F}}_{\epsilon }(q,t;R,r)$$

$${\mathbf{F}}_{\epsilon }(q,t;R,r)\cdot \mathbf{u}(q)=-\frac{(R+r_{\epsilon }(t;r)\cos (t))\sin (q)}{||\mathbf{x}-{\mathbf{x}}_1|{|}^3}-\frac{(R-r_{\epsilon }(t;r)\cos (t))\sin (q)}{||\mathbf{x}-{\mathbf{x}}_2|{|}^3}.$$

$$\dot{q}$$

$$\dot{p}$$

$$f_{\epsilon }(q,t;R,r)$$

$$||\mathbf{x}-{\mathbf{x}}_1||$$

$$||\mathbf{x}-{\mathbf{x}}_2||$$

$$H_{\epsilon }(q,p,t;R,r)=\frac{p^2}{2}+V_{\epsilon }(q,t;R,r),$$

$$V_{\epsilon }(q,t;R,r)$$

$$\begin{gathered} f_{\epsilon }(w,t;R,r)=-\frac{(R+r_{\epsilon }(t;r)\cos t)\sin \left(w/R\right)}{{\left[r_{\epsilon }^2(t;r)+2R(1-\cos \left(w/R\right))(R+r_{\epsilon }(t;r)\cos t)\right]}^{\frac{3}{2}}} \\ -\frac{(R-r_{\epsilon }(t;r)\cos t)\sin \left(w/R\right)}{{\left[r_{\epsilon }^2(t;r)+2R(1-\cos \left(w/R\right))(R-r_{\epsilon }(t;r)\cos t)\right]}^{\frac{3}{2}}} \end{gathered}$$

$$\begin{gathered} f_{\epsilon }(w,t;R,r)=-\frac{w\frac{\sin \left(w/R\right)}{\left(w/R\right)}\left(1+\frac{r_{\epsilon }(t;r)}{R}\cos t\right)}{{\left[r_{\epsilon }^2(t;r)+2w^2\frac{(1-\cos \left(w/R\right))}{{\left(w/R\right)}^2}\left(1+\frac{r_{\epsilon }(t;r)}{R}\cos t\right)\right]}^{\frac{3}{2}}} \\ -\frac{w\frac{\sin \left(w/R\right)}{\left(w/R\right)}\left(1-\frac{r_{\epsilon }(t;r)}{R}\cos t\right)}{{\left[r_{\epsilon }^2(t;r)+2w^2\frac{(1-\cos \left(w/R\right))}{{\left(w/R\right)}^2}\left(1-\frac{r_{\epsilon }(t;r)}{R}\cos t\right)\right]}^{\frac{3}{2}}}. \end{gathered}$$

$$\underset{R\to \infty }{\lim }{f}_{\epsilon }(w,t;R,r)=-\frac{2w}{{\left(r_{\epsilon }^2(t;r)+w^2\right)}^{3/2}},$$

$$\underset{r\to 0}{\lim }{f}_{\epsilon }(q,t;R,r)=-\frac{\sin (q)}{\sqrt{2}{R}^2(1-\cos (q){)}^{3/2}}.$$

$$-\frac{1}{R{q}^2}+\frac{1}{R(2\pi -q{)}^2},$$

$$V_1(q_1)=-\frac{1}{R^2\sqrt{2}(1-\cos (q_1){)}^{1/2}},$$

$$V_2(q_2)=-\frac{1}{R{q}_2}-\frac{1}{R(2\pi -q_2)},$$

$$H_i(p_i,q_i)=\frac{1}{2}{p}_i^2+V_i(q_i),$$

$$\begin{gathered} f_{\epsilon }(q,t;r)=-\frac{(1+r\rho (t;\epsilon )\cos (t))\sin (q)}{{\left[r^2{\rho }^2(t;\epsilon )+2(1-\cos q)(1+r\rho (t;\epsilon )\cos t)\right]}^{3/2}} \\ -\frac{(1-r\rho (t;\epsilon )\cos (t))\sin (q)}{{\left[r^2{\rho }^2(t;\epsilon )+2(1-\cos q)(1-r\rho (t;\epsilon )\cos t)\right]}^{3/2}}. \end{gathered}$$

$$\mathcal{X}_{\varepsilon}(q,p,s;r)=\left\{\begin{array}[]{rcl}\dot{q}&=&p\\ \dot{p}&=&{f}_{\varepsilon}(q,s;r)\\ \dot{s}&=&1\end{array}\right..$$

$${\mathbb{S}}_1(q,p,s)$$

$${\mathbb{S}}_2(q,p,s)$$

$${\mathbb{S}}_3(q,p,s)$$

$${\mathbb{S}}_4(q,p,s)$$

$$\mathbf{z}(s,t,\lambda )=\mathbf{x}(s,\lambda )-\mathbf{y}(t,\lambda ).$$

$$\frac{|\mathbf{z}(0,0,\lambda )|={\min }_{s,t}|\mathbf{z}(s,t,\lambda )|\text{buildrel}def}{=\delta (\lambda ),}$$

Full text

1. Introduction

1.1. A curved Sitnikov problem

We consider the following curved Sitnikov problem: Two bodies of equal masses (primaries) move, under mutual gravity, on Keplerian ellipses about their center of mass. A third, massless particle is confined to a circle passing through the center of mass of the primaries, denoted by $P_{0}$, and perpendicular to the plane of motion of the primaries; the second intersection point of the circle with that plane is denoted by $P_{1}$. We assume that the massless particle moves under the gravitational influence of the primaries without affecting them. The dynamics of the massless particle has two equilibrium points, at $P_{0}$ and $P_{1}$. We focus on the local dynamics near these two points, more precisely, on the dependence of the linear stability of these points on the parameters of the problem.

When the Keplerian ellipses are not too large or too small, $P_{0}$ is a local center and $P_{1}$ is a hyperbolic fixed point. When we increase the size of the Keplerian ellipses, as the distance between $P_{1}$ and the closest ellipse approaches zero, then $P_{1}$ undergoes stability interchanges. That is, there exists a sequence of open, mutually disjoint intervals of values of the semi-major axis of the Keplerian ellipses, such that, on each of these intervals the linearized stability of $P_{1}$ is strongly stable, and each complementary interval contains values where the linearized stability is not strongly stable, i.e., it is either hyperbolic or parabolic. The length of these intervals approaches zero when the semi-major axis of the Keplerian ellipses approaches the diameter of the circle on which the massless particle moves. This phenomenon is the main focus in the paper.

It is stated in [16] and suggested by numerical evidence [11, 12] that the linearized stability of the point $P_{0}$ also undergoes stability interchanges when the size of the binary is kept fixed and the eccentricity of the Keplerian ellipses approaches $1$.

Stability interchanges of the type described above are ubiquitous in systems of varying parameters; they appear, for example, in the classical Hill’s equation and in the Mathieu equation [24]. To prove the occurrence of this phenomenon in our curved Sitnikov problem, we first formulate a general result on stability interchanges for a general class of simple mechanical systems. More precisely, we consider the motion of two bodies — one massive and one massless — which are confined to a pair of curves and move under Newtonian gravity. We let the distance between the two curves be controlled by some parameter $\lambda$. We assume that the position of the infinitesimal particle that achieves the minimum distance between the curves is an equilibrium point. We show that, in the case when the minimum distance between the two curves approaches zero, corresponding to $\lambda\rightarrow 0$, there existence a sequence of mutually disjoint open intervals $(\lambda_{2n-1},\lambda_{2n})$, whose lengths approach zero as $\lambda\to 0$, such that whenever $\lambda\in(\lambda_{2n-1},\lambda_{2n})$ the linearized stability of the equilibrium point is strongly stable, and each complementary interval contains values of $\lambda$ where the linearized stability is not strongly stable. From this result we derive the above mentioned stability interchange result for the curved Sitnikov problem.

The curved Sitnikov problem considered in this paper is an extension of the classical Sitnikov problem described in Section 1.2 (also, see e.g., [18]). When the radius of the circle approaches infinity, in the limit we obtain the classical Sitnikov problem — the infinitesimal mass moves along the line perpendicular to the plane of the primaries and passing through the center of mass. The equilibrium point $P_{1}$ becomes the point at infinity and is of a degenerate hyperbolic type. Thus, stability interchanges of $P_{1}$ represent a new phenomenon that we encounter in the curved Sitnikov problem but not in the classical one. Also in the last case, it is well known that for $\varepsilon=0$ the classical Sitnikov problem is integrable. In the case of the curved one, numerical evidence suggests that it is not (see Figure 2).

The motivation for considering the curved Sitnikov problem resides in the $n$-body problem in spaces with constant curvature, and with models of planetary motions in binary star systems, as discussed in Section 1.3.

1.2. Classical Sitnikov problem

We recall here the classical Sitnikov problem. Two bodies (primaries) of equal masses $m_{1}=m_{2}=1$ move in a plane on Keplerian ellipses of eccentricity $\varepsilon$ about their center of mass, and a third, massless particle moves on a line perpendicular to the plane of the primaries and passing through their center of mass. By choosing the plane of the primaries the $xy$-plane and the line on which the massless particle moves the $z$-axis, the equations of motion of the massless particle can be written, in appropriate units, as

[TABLE]

where $r(t)$ is the distance from the primaries to their center of mass given by

[TABLE]

where $u(t)$ is the eccentric anomaly in the Kepler problem. By normalizing the time we can assume that the period of the primaries is $2\pi$, and

[TABLE]

for small $\varepsilon$.

When $\varepsilon=0$, i.e., the primaries move on a circular orbit and the dynamics of the massless particle is described by a 1-degree of freedom Hamiltonian and so is integrable. Depending on the energy level, one has the following types of solutions: an equilibrium solution, when the particle rests at the center of mass of the primaries; periodic solutions around the center of mass; escape orbits, either parabolic, that reach infinity with zero velocity, or hyperbolic, that reach infinity with positive velocity.

When $\varepsilon\in(0,1)$, the differential equation (1) is non-autonomous and the system is non-integrable. Consider the case $\varepsilon\ll 1$. The system also has bounded and unbounded orbits, as well as unbounded oscillatory orbits and capture orbits (oscillatory orbits are those for which $\limsup_{t\to\pm\infty}|z(t)|=+\infty$ and $\liminf_{t\to\pm\infty}|z(t)|<+\infty$, and capture orbits are those for which $\limsup_{t\to-\infty}|z(t)|=+\infty$ and $\limsup_{t\to+\infty}|z(t)|<+\infty$). In his famous paper about the final evolutions in the three body problem, Chazy introduced the term oscillatory motions [5], although he did not find examples of these, leaving the question of their existence open. Sitnikov’s model yielded the first example of oscillatory motions [23]. There are many relevant works on this problem, including [1, 2, 3, 17, 18, 22, 6, 9, 10].

The curved Sitnikov problem introduced in Section 1.1 is a modification of the classical problem when the massless particle moves on a circle rather than a line. Here we regard the circle as a very simple restricted model of a space with constant curvature. In Subsection 1.3 we introduce and summarize some aspects of this problem.

1.3. The $n$-body problem in spaces with

constant curvature

The $n$-body problem on spaces with constant curvature is a natural extension of the $n$-body problem in the Euclidean space; in either case the gravitational law considered is Newtonian. The extension was first proposed independently by the founders of hyperbolic geometry, Nikolai Lobachevsky and János Bolyai. It was subsequently studied in the late 19th, early 20th century, by Serret, Killing, Lipschitz, Liebmann, Schering, etc. Schrödinger developed a quantum mechanical analogue of the Kepler problem on the two-sphere in 1940. The interest in the problem was revived by Kozlov, Harin, Borisov, Mamaev, Kilin, Shchepetilov, Vozmischeva, and others, in the 1990’s. A more recent surge of interest was stimulated by the works on relative equilibria in spaces with constant curvature (both positive and negative) by Diacu, Pérez-Chavela, Santoprete, and others, starting in the 2010’s. See [7] for a history of the problem and a comprehensive list of references.

A distinctive aspect of the $n$-body problem on curved spaces is that the lack of (Galilean) translational invariance results in the lack of center-of-mass and linear-momentum integrals. Hence, the study of the motion cannot be reduced to a barycentric coordinate system.

As a consequence, the two-body problem on a sphere can no longer be reduced to the corresponding problem of motion in a central potential field, as is the case for the Kepler problem in the Euclidean space. As it turns out, the two-body problem on the sphere is not integrable [25].

Studying the three-body problem on spaces with curvature is also challenging. Perhaps the simplest model is the restricted three-body problem on a circle. This was studied in [8]. First, they consider the motion of the two primaries on the circle, which is integrable, collisions can be regularized, and all orbits can be classified into three different classes (elliptic, hyperbolic, parabolic). Then they consider the motion of the massless particle under the gravity of the primaries, when one or both primaries are at a fixed position. They obtain once again a complete classification of all orbits of the massless particle.

In this paper we take the ideas from above one step further, by considering the curved Sitnikov problem, with the massless particle moving on a circle under the gravitational influence of two primaries that move on Keplerian ellipses in a plane perpendicular to that circle. In the limit case, when the primaries are identified with one point, that is when the primaries coalesce into a single body, the Keplerian ellipses degenerate to a point, and the limit problem coincides with the two-body problem on a circle described above.

While the motivation of this work is theoretical, there are possible connections with the dynamics of planets in binary star systems. About 20 planets outside of the Solar System have been confirmed to orbit about binary stars systems; since more than half of the main sequence stars have at least one stellar companion, it is expected that a substantial fraction of planets form in binary stars systems. The orbital dynamics of such planets can vary widely, with some planets orbiting one star and some others orbiting both stars. Some chaotic-like planetary orbits have also been observed, e.g. planet Kepler-413b orbiting Kepler-413 A and Kepler-413 B in the constellation Cygnus, which displays erratic precession. This planet’s orbit is tilted away from the plane of binaries and deviating from Kepler’s laws. It is hypothesized that this tilt may be due to the gravitational influence of a third star nearby [13]. Of a related interest is the relativistic version of the Sitnikov problem [14]. Thus, mathematical models like the one considered in this paper could be helpful to understand possible types of planetary orbits in binary stars systems.

To complete this introduction, the paper is organized as follows: In Section 2 we go deeper in the description of the curved Sitnikov problem, studying the limit cases and its general properties. In Section 3 we present a general result on stability interchanges. In Section 4 we show that the equilibrium points in the curved Sitnikov problem present stability interchanges. Finally, in order to have a self contained paper, we add an Appendix with general results (without proofs) from Floquet theory.

2. The curved Sitnikov problem

2.1. Description of the model

We consider two bodies with equal masses (primaries) moving under mutual Newtonian gravity on identical elliptical orbits of eccentricity $\varepsilon$, about their center of mass. For small values of $\varepsilon$, the distance $r(t)$ from either primary to the center of mass of the binary is given by

[TABLE]

where $u(t)$ is the eccentric anomaly, which satisfies Kepler’s equation $u-\varepsilon\sin u=(2\pi/\tau)t$, where $\tau/2\pi$ is the mean motion111The mean motion is the time-average angular velocity over an orbit. of the primaries. The expansion in (4) is convergent for $\varepsilon<\varepsilon_{c}=0.6627...$; see [21].

A massless particle is confined on a circle of radius $R$ passing through the center of mass of the binary and perpendicular to the plane of its motion. We assume that the only force acting on the infinitesimal particle is the component along the circle of the resultant of the gravitational forces exerted by the primaries. The motion of the primaries take place in the $xy$-coordinate plane and the circle with radius $R$ is in the $yz$-coordinate plane. See Figure 1.

We place the center of mass at the point $(0,R,0)$ in the $xyz$-coordinate system. The position of the primaries are determined by the functions

[TABLE]

Note that $t=0$ corresponds to the passage of the primaries through the pericenter at $y=R\pm r(1-\varepsilon)$ and $t=\pi$ to the passage of the primaries through the apocenter at $y=R\pm r(1+\varepsilon)$; both peri- and apo-centers lie on the plane of the circle of radius $R$ in the $y$-axis.

The position of the infinitesimal particle is $\mathbf{x}(t)=(0,y(t),z(t))$ (taking into account the restriction of motion for the infinitesimal particle to the circle $y^{2}+z^{2}=R^{2}$). We will derive the equations of motion by computing the gravitational forces exerted by the primaries:

[TABLE]

where the distance from the particle to each primary is

[TABLE]

We note that when $r(1+\varepsilon)=2R$ the elliptical orbit of the primary with the apo-center at $y<R$ crosses the circle of radius $R$, hence collisions between the primary and the infinitesimal mass are possible. Therefore we will restrict to $r<\frac{2R}{1+\varepsilon}$; when $\varepsilon=0$, this means $r<2R$.

We write (2.1) in polar coordinates, that is $y=R\cos q$, $z=R\sin q$, and we obtain

[TABLE]

The origin $q=0$ corresponds to the point $(0,R,0)$ in the $xyz$-coordinate system. Thus, the primaries move on elliptical orbits around this point.

Next we will retain the component along the circle of the resulting force (2.1). That is, we will ignore the constraint force that confines the motion of the particle to the circle, as this force acts perpendicularly to the tangential component of the gravitational attraction force. The unit tangent vector to the circle of radius $R$ at $(0,R\cos(q),R\sin(q))$ pointing in the positive direction is given by $\mathbf{u}(q)=(0,-\sin(q),\cos(q))$. The component of the force $\mathbf{F}_{\varepsilon}(q,t;R,r)$ along the circle is computed as

[TABLE]

The motion of the particle, as a Hamiltonian system of one-and-a-half degrees of freedom, corresponds to

[TABLE]

where

[TABLE]

Hence

[TABLE]

where the potential is given by

[TABLE]

2.2. Limit cases

The curved Sitnikov problem can be viewed as a link between the classical Sitnikov problem and the Kepler problem on the circle, mentioned in Section 1.

2.2.1. The limit $R\rightarrow\infty$.

We express (7) in terms of the arc length $w=Rq$, obtaining

[TABLE]

which we can write in a suitable form as

[TABLE]

Letting $R$ tend to infinity we obtain

[TABLE]

which is the classical Sitnikov Problem.

2.2.2. The limit $r\rightarrow 0$.

When we take the limit $r\rightarrow 0$ in (7) we are fusing the primaries into a large mass at the center of mass and we obtain a two-body problem on the circle. The component force along the circle corresponds to

[TABLE]

This problem was studied in [8] with a different force given by

[TABLE]

the distance between the large mass and the particle is measured by the arc length (in that paper the authors assume that $R=1$). The potential of the force (12) is

[TABLE]

where $q_{1}$ denotes the angular coordinate, and the potential for (13) is

[TABLE]

where $q_{2}$ denotes the angular coordinate. Each problem defines an autonomous system with Hamiltonian

[TABLE]

taking $p_{1}=dq_{1}/dt$, $p_{2}=dq_{2}/dt$ and $i=1,2$. Let $\phi^{i}_{t}$ be the flow of the Hamiltonian $H_{i}$, and let $A_{i}$ denote the phase space, $i=1,\ 2$

Using that all orbits are determined by the energy relations given by (14), it is not difficult to define a homeomorphism $g:A_{1}\to A_{2}$ which maps orbits of system (12) into orbits of system (13). In the same way we can define a homeomorphism $h:A_{2}\to A_{1}$ which in fact is $g^{-1}$. This shows the $C^{0}$ equivalence of the respective flows. One can show that the two corresponding flows are $C^{0}$–equivalent.

We recall from [8] that the solutions of the two-body problem on the circle (apart from the equilibrium antipodal to the fixed body) are classified in three families (elliptic, parabolic and hyperbolic solutions) according to their energy level. The elliptic solutions come out of a collision, stop instantaneously, and reverse their path back to the collision with the fixed body. The parabolic solutions come out of a collision and approach the equilibrium as $t\rightarrow\infty$. Hyperbolic motions comes out of a collision with the fixed body, traverse the whole circle and return to a collision.

We remark that the two limit cases $R\to\infty$ and $r\to 0$ are not equivalent. Indeed, in the case $r\to 0$ the resulting system is autonomous, the point $q=0$ is a singularity for the system, and the point $q=\pi$ is a hyperbolic fixed point, while in the case $R\to\infty$ the resulting system is non-autonomous (for $\varepsilon\neq 0$), the point $q=0$ is a fixed point of elliptic type, and the point $q=\infty$ is a degenerate hyperbolic periodic orbit.

2.3. General properties

2.3.1. Extended phase space, symmetries, and equilibrium points

It is clear that, besides the limit cases $R\to\infty$ and $r\to 0$, the dynamics of the system depends only on the ratio $r/R$, so we can fix $R=1$ and study the dependence of the global dynamics on $r$ where $0<r<2$. In this case using (8) and (4) we get

[TABLE]

To study the non-autonomous system (8) we will make the system autonomous by introducing the time as an extra dependent variable

[TABLE]

This vector field is defined on $[0,2\pi]\times\mathbb{R}\times[0,2\pi]$, where we identify the boundary points of the closed intervals. The flow of $\mathcal{X}_{\varepsilon}$ possesses symmetries defined by the functions

[TABLE]

in the sense that

• •

$\mathbb{S}_{1}(\mathcal{X}_{\varepsilon}(q,p,s))=\mathcal{X}_{\varepsilon}(\mathbb{S}_{1}(q,p,s))$,

• •

$\mathcal{X}_{\varepsilon}(q,p,s)=\mathcal{X}_{\varepsilon}(\mathbb{S}_{2}(q,p,s))$,

• •

$\mathcal{X}_{\varepsilon}(q,p,s)=\mathcal{X}_{\varepsilon}(\mathbb{S}_{3}(q,p,s))$

• •

$\mathbb{S}_{4}(\mathcal{X}_{\varepsilon}(q,p,s))=-\mathcal{X}_{\varepsilon}(\mathbb{S}_{4}(q,p,s))$,

as can be verified by a direct computation.

The function $\mathbb{S}_{1}$ describes the symmetry respect to the trajectory $(0,0,s)$, $\mathbb{S}_{2}$ and $\mathbb{S}_{3}$ describe the bi-periodicity of ${f}_{\varepsilon}(q,s;r)$ and $\mathbb{S}_{4}$ describes the reversibility of the system.

System (8) has two equilibria $(0,0)$ and $(\pi,0)$, which correspond to periodic orbits for $\mathcal{X}_{\varepsilon}$.

While the classical Sitnikov equation is autonomous for $\varepsilon=0$, our equation (16) is not, and thus we expect it to be non–integrable, as is borne out by numerical simulation. Figure 2 shows a Poincaré section corresponding to $s=0$ (mod $2\pi$), with $\varepsilon=0$ and $r=1$. This simulation suggests that the invariant KAM circles coexist with chaotic regions.

In the sequel, we will analyze the dynamics around the equilibrium points $(\pi,0)$ and $(0,0)$. One important phenomenon that we will observe is that both equilibrium points undergo stability interchanges as parameters are varied. More precisely, when $\varepsilon$ sufficiently small is kept fixed and $r\rightarrow\frac{2R}{1+\varepsilon}$, the point $(\pi,0)$ undergoes infinitely many changes in stability, and when $r$ is kept fixed and $\varepsilon\to 1$, the point $(0,0)$ undergoes infinitely many changes in stability.

In the next section we will first prove a general result.

3. A general result on stability interchanges

In this section we switch to a more general mechanical model which exhibits stability interchanges. We consider the motion under mutual gravity of an infinitesimal particle and a heavy mass each constrained to its own curve and moving under gravitational attraction, and study the linear stability of the equilibrium point corresponding to the closest position between the particles along the curves they are moving on. In Section 4, we will apply this general result to the equilibrium point $P_{1}$ of the curved Sitnikov problem described in Section 1.1. The fact that in the curved Sitnikov problem there are two heavy masses, rather than a single one as considered in this section, does not change the validity of the stability interchanges result, since, as we shall see, what it ultimately matters is the time-periodic gravitational potential acting on the infinitesimal particle.

To describe the setting of this section, consider a particle constrained to a curve ${\bf x}={\bf x}(s,\lambda)$ in ${\mathbb{R}}^{3}$, where $s$ is the arc length along the curve and $\lambda$ is a parameter with values in some interval $[0,\lambda_{0}]$, $\lambda_{0}>0$. Another (much larger) gravitational mass undergoes a *prescribed * periodic motion according to ${\bf y}={\bf y}(t,\lambda)={\bf y}(t+1,\lambda)$; see Figure 3. We assume that the mass of the particle at ${\bf x}(s)$ is negligible compared to the mass at ${\bf y}(t)$, treating the particle at ${\bf x}(s)$ as massless.

To write the equation of motion for the unknown coordinate $s$ of the massless particle, let:

[TABLE]

Assume that $s=0,\ t=0$ minimize the distance between the two curves:

[TABLE]

for all $\lambda\in[0,\lambda_{0}]$, and that this minimum point is non-degenerate with respect to $t$, in the sense that

[TABLE]

Moreover, we make the following orthogonality assumption:

[TABLE]

In the sequel we will study the case when the minimum distance $\min_{s,t}|{\bf z}(s,t,\lambda)|\rightarrow 0$, that is, the shortest distance from the orbit ${\bf y}(t,\lambda)$ of the massive body to the curve ${\bf x}(s,\lambda)$ drawn by the infinitesimal mass approaches [math] as $\lambda\to\lambda_{0}$. See Figure 3. In the curved Sitnikov problem, this corresponds to the case when $r\to\frac{2R}{1+\varepsilon}$ ($r\to 2R$ when $\varepsilon=0$).

To write the equation of motion for $s$, we note that the Newtonian gravitational potential of the particle at ${\bf x}(s)$ is a function of $s$ and $t$ given by

[TABLE]

and the evolution of $s$ is governed by the Euler–Lagrange equation $\frac{d}{dt}L_{\dot{s}}-L_{s}=0$ with the Lagrangian

[TABLE]

leading to222to explain this form of the Lagrangian, we note that the equations for our massless particle are obtained by taking the limit of the particle of small mass $m$; for such a particle, in the ambient potential $U$, the Lagrangian is $\frac{1}{2}m\dot{s}^{2}-mU(s,t,\lambda)$ – the factor $m$ in front of $U$ is due to the fact that $U$ is the potential energy of the unit mass. Dividing the Euler–Lagrange equation by $m$ gives (22).

[TABLE]

Note that $s=0$ is an equilibrium for any $\lambda$, since $U^{\prime}(0,t,\lambda)=0$ for all $t$ and for all $\lambda$, according to (18). Linearizing (22) around the equilibrium $s=0$ we obtain

[TABLE]

where $a(t,\lambda)=U^{\prime\prime}(0,t,\lambda)$.

We have the following general result:

Theorem 1**.**

Assume that (18), (19), (20) hold, that ${\bf x}$ and ${\bf y}$ are both bounded in the $C^{2}$–norm uniformly in $\lambda$, and

[TABLE]

Then there exists an infinite sequence

[TABLE]

such that the equilibrium solution $s=0$ of (22) is linearly strongly stable333The equilibrium solution is linearly strongly stable if the Floquet multipliers of (23) lie on the unit circle and are not real, or equivalently, if the linearized system lies in the interior of the set of stable systems. for all $\lambda\in(\lambda_{2n},\lambda_{2n-1})$. Furthermore, each complementary $\lambda$–interval contains points where the linearized equilibrium is not strongly stable, i.e. is either hyperbolic or parabolic.

The proof of Theorem 1 relies on Lemmas 2, 3 and 4 stated below.

Lemma 2**.**

Consider the linear system

[TABLE]

where $a(t,\lambda)$ is a continuous function of $t\in[0,1]$, $\lambda\in(0,\bar{\lambda}]$, where $\bar{\lambda}>0$. Let $z(t,\lambda)=x+i\dot{x}$, where $x=x(t,\lambda)$ is a nontrivial solution of (26). Assume that there exists an interval $[t_{0}(\lambda),t_{1}(\lambda)]\subset[0,1]$, possibly depending on $\lambda$, such that

[TABLE]

Here $\arg z(t,\lambda)$ is defined as a continuous function of $t$. Although this choice of $\arg$ is unique only modulo $2\pi$, its increment as stated in equation (27) is uniquely defined.

Then there exists a sequence $\{\lambda_{k}\}_{k=0}^{\infty}$ satisfying (25) such that the Floquet matrix of (23) is strongly stable for all $\lambda\in(\lambda_{2n},\lambda_{2n-1})$, for any $n>0$. Furthermore, every complementary $\lambda$–interval contains values of $\lambda$ for which the Floquet matrix is not strongly stable.

Proof of this lemma can be found in [15].

Lemma 3**.**

Consider the linear system

[TABLE]

and assume that $a(t)>0$ on some interval $[t_{0},t_{1}]$. For any (nontrivial) solution $x(t)$ the corresponding phase vector $z(t)=x+i\dot{x}$ rotates by

[TABLE]

Proof.

Writing the differential equation $\ddot{x}+a(t)x=0$ as a system $\dot{x}=y,\ \dot{y}=-a(t)x$, we obtain (using complex notation):

[TABLE]

We conclude that for any solution of (28), the angle $\theta=\arg(x+iy)$ satisfies

[TABLE]

To invoke comparison estimates, consider $\bar{\theta}(t)$ which satisfies

[TABLE]

By the comparison estimate, we conclude:

[TABLE]

and the proof of the lemma will be complete once we show that $\bar{\theta}$ satisfies the estimate (29). To that end we consider a solution of

[TABLE]

with the initial condition satisfying

[TABLE]

This solution is of the form

[TABLE]

where $\varphi$ is chosen so as to satisfy (32).

Since $\arg(\bar{x}+i\dot{\bar{x}})$ satisfies the same differential equation as $\bar{\theta}$, and since the initial conditions match, we conclude that $\bar{\theta}(t)=\arg(\bar{x}+i\dot{\bar{x}})$, so that

[TABLE]

where $\widehat{X}$ denotes a quantity whose absolute value does not exceed $X$. In other words, $\bar{\theta}(t)$ is given by a linear function with coefficient $-\sqrt{a_{m}}$, up to an error $<\pi/2$. The last inequality is due to the fact that the complex numbers $\cos(\sqrt{a_{m}}t-\phi)-i\sin(\sqrt{a_{m}}t-\phi)$ and $\cos(\sqrt{a_{m}}t-\phi)-i\sqrt{a_{m}}t\sin(\sqrt{a_{m}}t-\phi)$ lie in the same quadrant, so the difference between their arguments is no more than $\pi/2$.

Therefore, over the interval $[t_{0},t_{1}]$ the function $\bar{\theta}$ changes by the amount $\sqrt{a_{m}}(t_{1}-t_{0})$ with the error of at most $\frac{\pi}{2}+\frac{\pi}{2}=\pi$:

[TABLE]

restating this more formally, (33) implies

[TABLE]

Subtracting the second inequality from the first gives (33).

Substituting (34) into (31) yields (29) and completes the proof of Lemma 3. $\diamondsuit$ ∎

Lemma 4**.**

Consider the potential $U$ defined by equation (21). Assume that the functions ${\bf x}={\bf x}(s,\lambda)$, ${\bf y}={\bf y}(t,\lambda)$ satisfy

• (i)

the minimum $\min_{s,t}|{\bf z}(\cdot,\cdot,\lambda)|=\delta(\lambda)$ is non-degenerate and is achieved at $s=t=0$, where ${\bf z}(s,t,\lambda)={\bf x}(s,\lambda)-{\bf y}(t,\lambda)$,

• (ii)

$||{\bf z}(\cdot,\cdot,\lambda)||_{C^{2}}\leq M$* uniformly in $0<\lambda<\lambda_{0}$, and*

• (iii)

$\delta(\lambda)\rightarrow 0$* as $\lambda\rightarrow 0$.*

Then there exists time $\tau=\tau(\lambda)$ (approaching zero as $\lambda\rightarrow 0$) such that

[TABLE]

Proof.

Differentiating (21) with respect to $s$ twice, we get

[TABLE]

We now estimate all the dot products in the above expression to obtain a lower bound.

First,

[TABLE]

since $s$ is the arc length, and from here ${\bf z}^{\prime}\cdot{\bf z}^{\prime\prime}=0$. Now, to estimate ${\bf z}\cdot{\bf z}$ and ${\bf z}\cdot{\bf z}^{\prime}$ we observe that ${\bf z}\cdot{\bf z}^{\prime}=\frac{1}{2}({\bf z}\cdot{\bf z})^{\prime}$ and we note that the first expression, as a function of $t$ with $s=0$ fixed has a minimum at $t=0$ that we call $\delta^{2}$, and that the second function vanishes at $t=0$. Applying Taylor’s formula with respect to $t$ we then have

[TABLE]

where the constant $k$ is determined by the $C^{2}$–norm $M$ of ${\bf z}$. For the remaining dot product we have (still keeping $s=0$ and $t$ arbitrary):

[TABLE]

Using the above estimates in (36), we obtain

[TABLE]

Now we restrict $t$ to have $\delta^{2}+kt^{2}\leq 2\delta^{2}$; this guarantees that the denominator in (40) does not exceed $(2\delta)^{5/2}$; to bound the numerator, we further restrict $t$ so that the dominant part $\delta^{2}-3k^{2}t^{2}\geq\frac{1}{2}\delta^{2}$, thus bounding the numerator from below by

[TABLE]

if $\delta$ is sufficiently small. Summarizing, we restricted $t$ to

[TABLE]

and showed that for all such $t$ and for $\delta$ small enough

[TABLE]

where $c_{1}=2^{-9/2}$. With $\tau$ defined in (41) we obtain $\lim_{\tau\to 0}\left(\tau\cdot\min_{|t|\leq\tau}\sqrt{a(t,\lambda)}\right)=\infty$, thus completing the proof of the lemma. $\diamondsuit$ ∎

Proof of Theorem 1 Consider the linearized equation

[TABLE]

and consider the phase point $z=x+i\dot{x}$ of a nontrivial solution. Lemma 3 gives us the rotation estimate (29) for any time interval $[t_{0},t_{1}]$; let this interval be $[-\tau(\lambda),\tau(\lambda)]$ where $\tau(\lambda)$ is taken from the statement of Lemma 4. We then have from (29): (29)

[TABLE]

According to the conclusion (35) of Lemma 4, $\theta\rightarrow-\infty$ as $\lambda\rightarrow 0$. This satisfies the condition (27) Lemma 2, which now applies and its conclusion comples the proof of Theorem 1. $\diamondsuit$

4. Stability of the equilibrium points in the curved Sitnikov problem

The system (8) has two equilibria $(0,0)$ and $(\pi,0)$, which correspond to periodic orbits for $\mathcal{X}_{\varepsilon}$. The associated linear system around the fixed point $(q_{*},p_{*})$ can be written as

[TABLE]

with

[TABLE]

Let $X(t)$ be a fundamental matrix solution of system (45) given by

[TABLE]

with the initial condition $X(0)=I$, the identity matrix; $x_{1}$ is an even function and $x_{2}$ and odd one since $\partial{f}_{\varepsilon}/\partial q$ is an even function with respect to $t$. The monodromy matrix is given by $X(2\pi)$ and we denote $\lambda_{1},\lambda_{2}$ its eigenvalues, the Floquet multipliers associated to (45). These are given by

[TABLE]

the trace $\textrm{Tr}(X(2\pi))=x_{1}(2\pi)+y_{2}(2\pi)$ determines the linearized dynamics around the fixed point. Moreover, since the function $(\partial{f}_{\varepsilon}/\partial q)_{\mid q=q^{*}}$ is an even function, we know $x_{1}(2\pi)=y_{2}(2\pi)$ (see [24]) and then

[TABLE]

To emphasize the dependence on the parameters $\varepsilon,r$ we write

[TABLE]

Thus the linear stability of the equilibrium is

1. (1)

Elliptic type: $|y_{2}(2\pi;\varepsilon,r)|<1$. 2. (2)

Parabolic type: $|y_{2}(2\pi;\varepsilon,r)|=1$. 3. (3)

Hyperbolic type: $|y_{2}(2\pi;\varepsilon,r)|>1$.

Since the Wronskian is equal to $1$ for all $t$, then

[TABLE]

and we have:

[TABLE]

From this last expression it follows that, in the parabolic case, the periodic orbit corresponding to the equilibrium point $(\pi,0)$ is associated to $x_{1}(t)$ or to $x_{2}(t)$ (or to both).

4.1. Stability of the equilibrium point $(\pi,0)$. Case $\varepsilon=0$.

In this subsection we consider the system (8)-(15) for the case $\varepsilon=0$, namely when the primaries are following circular trajectories. We also fix $R=1$. The extended vector field under study is

[TABLE]

where

[TABLE]

[TABLE]

with $0<r<2$. We observe that in this case, function (51) is $\pi$–periodic (remember that ${f}_{\varepsilon}(q,t;r)$ is $2\pi$–periodic if $\varepsilon>0$).

The linear system associated to (50) around the fixed point $(\pi,0)$ is defined by the function

[TABLE]

which is $C^{1}$ respect to $t$ and $r$.

Lemma 5**.**

The function $\frac{\partial{f}_{0}}{\partial q}(\pi,t;r)$ is monotone decreasing with respect to $r$, that is, for all $t\in[0,\pi)$

[TABLE]

if $r_{1}<r_{2}$ (see Figure 4).

Proof.

Let be $F(t,r)=\frac{\partial{f}_{0}}{\partial q}(\pi,t;r)$, then by straightforward computation we get

[TABLE]

$\diamondsuit$ ∎

Proposition 6**.**

The equilibrium point $(\pi,0)$ is of hyperbolic type if $r\leq\left(\sqrt{17}-3\right)^{1/2}=1.059\dots$.

Proof.

We first show that $F(t,r)\geq 0$ for all $t$ if and only if $r\leq\left(\sqrt{17}-3\right)^{1/2}$. We compute

[TABLE]

We have $F(t,r)=F\left(\pi-t,r\right)$ so it is enough to restrict $t\in[0,\frac{\pi}{2}]$. Note that $\frac{\partial F}{\partial t}(t,r)=0$ for $t=0,\pi/2$.

For $0<r<2$, we have

[TABLE]

and

[TABLE]

It follows immediately that $F(0,r)<0$ if $r>(\sqrt{17}-3)^{1/2}$. If $r\cos t\leq 1$ then (52) implies $F(t,r)\geq 0$. Hence $F(t,r)\geq 0$ for all $t$ provided $r\leq 1$. Let $1<r\leq(\sqrt{17}-3)^{1/2}$. If $t\in[0,\cos^{-1}(1/r)]$, which is equivalent to $r\cos t>1$, then

[TABLE]

Therefore (53) implies that $F(t,r)$ is increasing in $t$ for $t\in[0,\cos^{-1}(1/r)]$, and, since $F(0,r)\geq 0$, it follows that $F(t,r)\geq 0$ for all $t$.

Now let $x(t)=x_{1}(t)+x_{2}(t)$ be a particular solution of system (45); by hypothesis it satisfies $x(0)=1,\dot{x}(0)=1$. From (45) with $\varepsilon=0$ we obtain

[TABLE]

Then, since $F(t,r)\leq 0$ for all $t$, we can apply directly Lyapunov’s instability criterion (see for instance page 60 in [4]) to show that $(\pi,0)$ is of hyperbolic type. $\diamondsuit$∎

We can in fact estimate the first value $r_{1}$ of $r$ at which the equilibrium point $(\pi,0)$ becomes of parabolic type for the first time. The idea is to extend slightly the result beyond $r=\left(\sqrt{17}-3\right)^{1/2}$, and then find the maximum of the $r$ for which Proposition 6 holds. In the way we have to do straightforward analytic but tedious computations that we decide to avoid in this paper. Finally we get that the value of $r$ to have the first parabolic solution is $r_{1}\approx 1.2472\cdots$.

In Figure 5 we show a couple of numerical simulations which illustrate the stability interchanges of the equilibrium point $(\pi,0)$ for $r\in(0,2)$ and $\varepsilon=0$. The figure on the right hand side is a plot for $r$ close to 2.

4.2. Stability of the equilibrium point $(\pi,0)$. Case $\varepsilon\neq 0$.

In this case, for every $0<r<(\sqrt{17}-3)^{1/2}$, $\partial f_{0}/\partial q(\pi,t;r)>0$, hence $\partial f_{\varepsilon}/\partial q(\pi,t;r)>0$ for all $\varepsilon>0$ sufficiently small (depending on $r$). Therefore $(\pi,0)$ remains an equilibrium point of hyperbolic type in the case when the primaries move on Keplerian ellipses of sufficiently small eccentricity $\varepsilon>0$.

We remark that Theorem 1 does not depend of the shape of the curves where ${\bf x}(x)$ or ${\bf y}(t)$ are moving, in other words we can apply Theorem 1 independently if the primaries are moving on a circle or on ellipses of eccentricity $\varepsilon>0$. Thus we obtain the following result on stability interchanges:

Theorem 7**.**

In the curved Sitnikov problem, let us fix any $\varepsilon\in[0,1)$, let $R=1$, and consider $r$ (the semi-major axes of the Keplerian ellipses traced out by the primaries) as the parameter. As $r$ approaches $\frac{2}{1+\varepsilon}$, the distance between a Keplerian ellipse and the circle of the massless particle approaches zero. There exists a sequence $r_{n}\uparrow\frac{2}{1+\varepsilon}$ satisfying

[TABLE]

such that the equilibrium point $(\pi,0)$ of equation (45) is strongly stable for $r\in(r_{2n-1},r_{2n})$ and not strongly stable for some $r$ in the complementary intervals $(r_{2n},r_{2n+1})$. In other words, the equilibrium point $(\pi,0)$ loses and then regains its strong stability infinitely many times as $r$ increases towards $\frac{2}{1+\varepsilon}$.

Proof of Theorem 7.

We verify that Theorem 1 applies. The curve ${\bf x}(s,\lambda)$ of Theorem 1 is represented by the circle of radius $R=1$, the curve ${\bf y}(t,\lambda)$ is represented by the orbit of the primary that gets closer to $(\pi,0)$ (when $\varepsilon=0$ the primaries are co-orbital), and the parameter $\lambda$ corresponds to $r$. The planes of the two curves are perpendicular, as in (20), and the minimum distance between the curves, given by (18), is non-degenerate as in (19), and it corresponds to the infinitesimal mass being at $y=-1$ and the closest primary to this point being at $y=1-r(1+\varepsilon)$. Thus, the minimum distance is $\delta(r)=2-r(1+\varepsilon)$, and it approaches [math] when $r\to\frac{2}{1+\varepsilon}$. To apply Theorem 1 we only need to verify that ${\bf x}$ and ${\bf y}$ are bounded in the $C^{2}$-norm uniformly in $r$. This is obviously true for ${\bf y}(t,r)$ since the motion of the primary is not affected by the motion of the infinitesimal mass; it is also true for ${\bf x}(s,r)$ since $s=$arclength and the motion lies on the circle of radius $R=1$, hence $\|{\bf x}(s,r)\|=1$ (the radius of the circle), $\|{\bf x}^{\prime}(s,r)\|=1$ (the unit speed of a curve parametrized by arc-length), and $\|{\bf x}^{\prime\prime}(s,r)\|=1$ (the curvature of the circle). Hence the conclusion of Theorem 1 follows immediately. $\diamondsuit$∎

4.3. Stability of the equilibrium point $q=(0,0)$.

The linear stability of $(0,0)$ is the same as of the barycenter in the classical Sitnikov problem, since

[TABLE]

The stability of this point can be treated very similarly to that of the origin for Hill’s equation, so in this analysis we use results from that theory.

As in the study of the other equilibrium point we start with the case $\varepsilon=0.$ Here the function $f_{0}(q,t;r)$ defined in equation (8) around $q=(0,0)$ is given by

[TABLE]

The local dynamics is determined by the linear part. The eigenvalues are on the unit circle, given by $\sigma_{1,2}=\pm i\sqrt{2/r^{3}}$, and the Floquet multipliers, which come from the monodromy matrix

[TABLE]

are $\lambda_{1,2}=e^{\pm i\sqrt{2/r^{3}}\pi}$.

Hence $q=(0,0)$ is of elliptic type if $\sqrt{2/r^{3}}\neq k$ for any $k\in\mathbb{Z}$, and it is of parabolic type if $\sqrt{2/r^{3}}=k$ for some $k\in\mathbb{Z}$. In fact, in the parabolic case, if $\sqrt{2/r^{3}}=2m$, $m\in\mathbb{Z}$, then there exists a $\pi$-periodic solution, and if $\sqrt{2/r^{3}}=2m+1$, $m\in\mathbb{Z}$, then there exists a $2\pi$-periodic solution.

Proposition 8**.**

The equilibium point $q=(0,0)$ of the system defined by (55) is stable for all $r\in(0,2)$.

Proof.

The linear part possess Floquet multipliers which place the system in either the elliptic or the parabolic case. In the last case, when $\sqrt{2/r^{3}}=k$ for some $k\in\mathbb{Z}$, there exists two independent eigenvectors associated to the eigenvalues $\lambda_{1,2}$. Then the conjugacy of $A$ and $\pm I$ implies that the monodromy matrix is $A=\pm I$ and $q=(0,0)$ is stable (see [20] for more details).

When the origin is of elliptic type for the linear system, we observe that the coefficient of the term of order $3$ in equation (55) is greater than zero for all $r\in(0,2)$, then we can apply directly Ortega’s theorem (see Appendix and [19]), and therefore we obtain that the equilibrium point $(0,0)$ is stable for the whole system, that is, including the nonlinear part. $\diamondsuit$ ∎

We note that, in the case when $\varepsilon=0$, stability interchanges in a weak sense appear as $r\to 0$, since $(0,0)$ switches between elliptic type when $r\neq(2/k^{2})^{1/3}$ and parabolic type when $r=(2/k^{2})^{1/3}$, $k\in\mathbb{Z}^{+}$, as noted before.

In the case when $\varepsilon\neq 0$, as we mentioned earlier, the linear stability of $(0,0)$ is the same as in the classical Sitnikov problem, so it only depends on the eccentricity parameter $\varepsilon$. The papers [16, 11, 12] state that there are stability interchanges when the size $r$ of the binary is kept fixed and the eccentricity $\varepsilon$ of the Keplerian ellipses approaches $1$. We should point out that Theorem 1 does not apply to this case since the function $r_{\varepsilon}(t,r)$ describing the motion of the primaries — corresponding to $y(t,\lambda)$ in Theorem 1 — does not remain bounded in the $C^{2}$ norm uniformly in $\varepsilon$.

Appendix A Floquet Theory

In order to have a self contained paper we add this appendix with the main results on Floquet theory, must of them are very well known for people in the field.

Consider the linear system

[TABLE]

where $A(t)$ is a $T$-periodic matrix-valued function. Let $X(t)$ be the fundamental matrix solution

[TABLE]

with the initial condition $X(0)=I$, the identity matrix. Let $\lambda_{1}$ and $\lambda_{2}$ be the eigenvalues (Floquet multipliers) of the monodromy matrix matrix $X(T)$ and let $\mu_{1},\mu_{2}$ (Floquet exponents) be such that $\lambda_{1}=e^{\mu_{1}T},\lambda_{2}=e^{\mu_{2}T}$.

Theorem 9** (Floquet’s theorem).**

Suppose $X(t)$ is a fundamental matrix solution for (56), then

[TABLE]

for all $t\in\mathbb{R}$. Also there exists a constant matrix $B$ such that $e^{TB}=X(T)$ and a $T$-periodic matrix $P(t)$, so that, for all $t$,

[TABLE]

Lemma 10**.**

Let be $\lambda$ a Floquet multiplier for (56) and $\lambda=e^{\mu T}$, then there exists a nontrivial solution $x(t)=e^{\mu t}p(t)$, with $p(t)$ a $T$-periodic function. Moreover, $x(t+T)=\lambda x(t)$.

Thus, the Floquet multipliers lead to the following characterization:

• •

If $|\lambda|<1\Leftrightarrow\textrm{Re}(\mu)<0$ then $x(t)\rightarrow 0$ as $t\rightarrow\infty$.

• •

If $|\lambda|=1\Leftrightarrow\textrm{Re}(\mu)=0$ then $x(t)$ is a pseudo-periodic, bounded solution. In particular when $\lambda=1$ then $x(t)$ is $T$-periodic and when $\lambda=-1$ then $x(t)$ is $2T$-periodic.

• •

If $|\lambda|>1$ then $\textrm{Re}(\mu)>0$ and therefore $x(t)\rightarrow\infty$ as $t\rightarrow\infty$, an unbounded solution.

Lemma 11**.**

If the Floquet multipliers satisfy $\lambda_{1}\neq\lambda_{2}$, then the equation (56) has two linearly independent solutions

[TABLE]

where $p_{1}(t)$ and $p_{2}(t)$ are $T$-periodic functions and $\mu_{1}$ and $\mu_{2}$ are the respective Floquet exponents.

In this way, the stability of the solution to (56) is

• •

Asymptotically stable if $|\lambda_{i}|<1$ for $i=1,2$.

• •

Lyapunov stable if $|\lambda_{i}|\leq 1$ for $i=1,2$, or if $|\lambda_{i}|=1$ and the algebraic multiplicity equals the geometric multiplicity.

• •

Unstable if $|\lambda_{i}|>1$ for at least one $i$, or if $|\lambda_{i}|=1$ and the algebraic multiplicity is greater than the geometric multiplicity.

A particular case for the equation (56) is the so-called Hill’s equation, namely the periodic linear second order differential equation:

[TABLE]

where $f(t)$ is a $\pi$-periodic function. Equation (58), as a first order system, is

[TABLE]

Let (57) be a fundamental matrix solution of (59). The monodromy matrix corresponds to $X(\pi)$ and, as before, let $\lambda_{1}$ and $\lambda_{2}$ be the Floquet multipliers associated to (59) and $\mu_{1}$ and $\mu_{2}$ be the corresponding Floquet exponents.

In this paper we assume that $f(t)$ is an even function; from Floquet theory we know that $x_{1}$ is an even and $x_{2}$ is an odd function. Also, the trace of $A(t)$ vanishes and

[TABLE]

Assuming $\mu_{1}=a+ib$ is the Floquet exponent corresponding to the Floquet multiplier $\lambda_{1}$, the general solution to (59) is characterized as follows:

• •

Elliptic type: $\lambda_{1}\in\mathbb{C}\setminus\mathbb{R}$, with $|\lambda_{1}|=1$ (and $\lambda_{2}=\bar{\lambda}_{1}$). The general solution is pseudo-periodic and can be written as

[TABLE]

The origin is Lyapunov stable.

• •

Parabolic type: $\lambda_{1}=\lambda_{2}=\pm 1$.

• –

If $\lambda_{1}=1$ and there are two linearly independent eigenvectors of the monodromy matrix, the general solution is

[TABLE]

and is $\pi$-periodic and Lyapunov stable.

• –

If $\lambda_{1}=1$ and there is just one eigenvector associated to this eigenvalue, thus the general solution is

[TABLE]

and is unstable.

• –

If $\lambda_{2}=-1$ and there are two linearly independent eigenvectors of the monodromy matrix, the general solution is

[TABLE]

and is $2\pi$-periodic and Lyapunov stable.

• –

If $\lambda_{1}=-1$ and there is just one eigenvector associated to this eigenvalue, thus the general solution is

[TABLE]

and is unstable.

• •

Hyperbolic type: $\lambda_{1}\in\mathbb{R}$, but $|\lambda_{1}|\neq 1$ (and $\lambda_{2}=1/\lambda_{1}$).

• –

If $\lambda_{1}>1$, then the solution is

[TABLE]

• –

If $\lambda_{1}<-1$, then the solution is

[TABLE]

Thus, the origin is unstable.

A useful result from Rafael Ortega [19] considers the nonlinear Hill equation

[TABLE]

with $n\geq 2$, where the functions $a,c:\mathbb{R}\rightarrow\mathbb{R}$ are continuous, $T$-periodic, and $\int_{0}^{T}|c(t)|\neq 0$, and the function $d:\mathbb{R}\times(-\epsilon,\epsilon)\rightarrow\mathbb{R}$, for $\epsilon>0$, is continuous, has continuous derivatives of all orders respect to $y$, is $T$-periodic function respect to $t$, and $d(t,y)=\mathcal{O}\left(|y|^{2n}\right)$ as $y\rightarrow 0$ uniformly with respect to $t\in\mathbb{R}$.

The linear part around the solution $y=0$ of (60) is

[TABLE]

Theorem 12** (Ortega’s Theorem).**

Assume the following:

1. (1)

The equation (61) is stable. 2. (2)

$c\geq 0$* or $c\leq 0$.*

Then $y=0$ is a stable solution of (60).

Acknowledgements

We thank to the anonymous referees, their remarks and suggestions help us to improve this paper. Research of M.G. was partially supported by NSF grant DMS-1515851. M. L. gratefully acknowledges support by the NSF grant DMS-1412542. The fourth author (EPC) has received partial support by the Asociación Mexicana de Cultura A.C. Parts of this work have been done while the authors visited CIMAT, Guanajuato, LFP and EPC visit Yeshiva University and M.G. visited UAM-I in Mexico City. All authors are grateful for the hospitality of these institutions.