Global and local stability for ghosts coupled to positive energy degrees of freedom

Negative kinetic energies correspond to ghost degrees of freedom, which are potentially of relevance for cosmology, quantum gravity, and high energy physics. We present a novel wide class of stable mechanical systems where a positive energy degree of freedom interacts with a ghost. These theories have Hamiltonians unbounded from above and from below, are integrable, and contain free functions. We show analytically that their classical motion is bounded for all initial data. Moreover, we derive conditions allowing for Lyapunov stable equilibrium points. A subclass of these stable systems has simple polynomial potentials with stable equilibrium points entirely due to interactions with the ghost. All these findings are fully supported by numerical computations which we also use to gather evidence for stability in various nonintegrable systems.


I. INTRODUCTION
Ghosts are dynamical degrees of freedom (DoF), described by a Hamiltonian with negative kinetic energies, as opposed to the kinetic energies of standard degrees of freedom.Such ghosts generically appear in the context of theories with higher derivatives [1] (see also [2][3][4] for reviews) and have been considered in a variety of physical scenarios: they have been advocated as a way to regulate the ultraviolet behaviour of quantum field theories [5][6][7][8]; they have been suggested as a way to address the cosmological constant problem [9][10][11]; they have been used to obtain bouncing cosmologies [12]; and they have been proposed in the context of dark energy [13] or, more recently, as a way to ameliorate the Hubble tension [14].
In spite of these applications, ghosts are typically disregarded since they are expected to generically lead to catastrophic instabilities, both at the classical and at the quantum level (see e.g [15]) -for the purpose of this work, we will only be concerned with the classical case.
Whenever positive-and negative-kinetic-energy degrees of freedom do not interact, the presence of a ghost is trivially harmless.Without coupling, no energy exchange can take place and the respective ghostly and non-ghostly Hamiltonians can each be stable due to suit-able self-interactions.
In contrast, it is widely expected that a generic coupling between the positive-and negative-kineticenergy degrees of freedom, leads to an energy exchange between them that instigates an unbounded growth of absolute values of both energies.Hence, one expects a runaway evolution to larger and larger phase-space variables.However, in previous work [16], coauthored by some of the authors of the current paper, a specific counterexample to such a catastrophe was presented.Besides completing a stability proof for a class of models introduced in [16], the key physical motivation for this work can be understood in terms of the following simple question: Is the stable model from [16] too specific or can stability be extended to a wider class of integrable and nonintegrable systems with ghosts?
The potential for stable motion in the presence of interacting ghosts has, of course, been previously considered.(While we will be concerned only with classical theories, we refer to [6,8,[17][18][19][20][21][22][23][24][25][26][27] for proposals at the quantum level.)At the classical level, it has been argued that theories with ghosts can be considered "benign", provided that the runaway is sufficiently mild and the motion can be extended to the infinite future [25,[28][29][30][31].At the same time, indications have been found that classical motion in the presence of ghosts can exhibit so-called "islands of stability" (i.e., bounded motion for a restricted set of initial conditions) [21,[28][29][30][31][32][33][34][35][36].All of these indications are numerical (and, as such, cannot be fully conclusive, as they cannot cover the whole evolution of the system to the infinite future along even a single trajectory, not to mention all the Hamiltonian trajectories) and/or hold only for restricted sets of initial conditions, i.e., do not conclusively advocate for global stability.Moreover, we still lack a physical rationale or criterion to distinguish potentially stable ghostly interactions from unstable ones.
As mentioned before, the recent work [16] has established the first conclusive proof of globally stable motion for a positive-energy harmonic oscillator interacting through a specific non-polynomial potential with a ghost in the form of a negative-energy oscillator.It was shown that, for arbitrary initial conditions, the motion remains bounded.The proof relies on the integrable nature of the model.We stress that this kind of global stability (also referred to as the Lagrange stability) is different from the above notions of "benign" ghosts and/or "islands of stability".We also stress that there is no direct connection of global stability to local (Lyapunov) stability [37].In this context, we would like to mention that [16] has also proven Lyapunov stability of the equilibrium configuration of the model.It is global (Lagrange) stability that we are mostly concerned with and we will, from here on, refer to it simply as stability, everywhere where this would not cause confusion.One of the reason for this concern is that quantum mechanics effectively probes the whole configuration space.Thus, regions beyond the, usually small, "islands of stability" can force the wave function and probability distributions to runaway to large values of phase-space coordinates in the quantized theory.
The purpose of this paper is the following: First, we extend the proof in [16] to a much larger class of integrable models with kinetic terms quadratic in canonical momenta and with two free functions of (particular combinations) of the coordinates canonically conjugated to them.This class of systems encompasses, in particular, a tower of stable polynomial interactions.Our proof also covers the subclass that has been announced in [16].Second, we prove that a wide subclass of these models possesses locally stable vacua.Third, we propose physical criteria for stable motion of more generic systems with ghosts interacting with positive energy degrees of freedom.Fourth, we investigate the above numerically, both for integrable and non integrable models.
With these goals in mind, we start in Sec.(II) by introducing a generic procedure to relate classical mechanical systems of N degrees of freedom with and without ghosts by means of suitable complex canonical transformations.We then consider a particular class of integrable models in Sec.(III), for which in Subsec.(III A) we provide a fully analytic proof of the boundedness of motion, identify a polynomial subclass in Subsec.(III C), and discuss Lyapunov stability of equilibrium points in Subsec.(III E).In particular, we find that new Lyapunov stable vacua can be instigated by polynomial interactions with a ghost.Incidentally, we stress that the polyno-mial models introduced in Subsec.(III C) stay polynomial and integrable when the degrees of freedom all have positive kinetic energies and appear as such to be novel to this work.Finally, we extend the discussion to nonintegrable models in Sec.(IV), for which we gather numerical evidence that stable (or at least longlived) motion can persist in absence of integrability.We close with our conclusions in Sec.(V).We provide additional technical material in Appendices.In particular, in App.(A) we list various classes of integrable models; in App.(B) we derive equations defining the models with the fist integral quadratic in canonical momenta and demonstrate that Lyapunov stable equilibrium points are located in the saddle points of the potential; in App.(C) we discuss stability in alternatively complexified Liouville models; in App.(D) and App.(E) we provide more details on integrable models with polynomial potentials.

II. INTEGRABLE INTERACTING GHOSTS VIA COMPLEX CANONICAL TRANSFORMATIONS
We start by considering N classical mechanical degrees of freedom (DoFs) with a Hamiltonian where p i denotes the canonical momenta conjugate to the coordinates x i .Together, these denote the phase-space variables ξ i = (x i , p i ).For now, all N DoFs have samesign kinetic terms and thus correspond to modes with positive kinetic energy, as the masses m i are assumed to be positive.We notice that a rotation in the complex plane to the imaginary axis of any of the momentum variables, i.e., p n → ±ip n , flips the sign of the respective kinetic term, i.e., p 2 n /m n → −p 2 n /m n , and thus turns the n-th degree of freedom (DoF) with positive kinetic energy into a DoF with negative kinetic energy.Hence, this transformation effectively flips the sign of the mass m n .This rotation in the complex plane can be completed into a complex canonical transformation by simultaneously transforming the respective x n → ∓ix n .The resulting (univalent) transformation x j → ±i x j for j = n p j → ∓i p j for j = n is canonical in the sense that it preserves the Poisson bracket {} PB , i.e., Given that C ± x n is a canonical transformation, it moreover preserves the Hamilton equations, the time-independence of the Hamiltonian itself, as well as all the timeindependence of any other existing constant of motion.We note also that Moreover, any function F ( ⃗ ξ) of the phase-space variables ⃗ ξ which is even under the simultaneous change of the sign of the variables x n and p n , i.e., is invariant under the successive application of two identical C ± x n .We will also refer to such functions as '(x n ,p n )-parity even'.On the other hand, we refer to a potential satisfying V (x i̸ =n , x n ) = V (x i̸ =n , −x n ) as 'x n -parity even'.Obviously, the Hamiltonian in Eq. ( 1) is (x n ,p n )-parity even if and only if the potential is x n -parity even.
In general, C ± x n does not preserve real-valuedness of the Hamiltonian.(The same, and all of what follows, also holds for other constants of motion.)However, realvaluedness is obviously preserved if the potential is polynomial (or Taylor expandable with respect to x n ) and x n -parity even.This is the case, e.g., for the canonical Hamiltonian of N non interacting positive kinetic-energy modes (where here and henceforth we normalize the kinetic terms just to be ±p 2 i /2, i.e., we set the masses to one, but keep the possibility to have a different frequency ω i for each mode) In this case, the canonical transformation C ± x n transforms a real-valued Hamiltonian with N positive-energy modes to a dual real-valued Hamiltonian with N − 1 positive kinetic-energy DoFs and one negative kinetic-energy DoF, i.e., a ghost.This can easily be extended by performing several distinct such canonical transformations to introduce multiple ghost degrees of freedom.For the general Hamiltonian in Eq. ( 1), if the potential is parity-even in all of the respective variables, then all of these Hamiltonians are real-valued as far as the original one is.All of the above also holds if further parity-even interactions are included.
In the rest of this work, we will specify to systems with two degrees of freedom x and y and with a Hamiltonian where σ = ±1 decides between the case of two positive kinetic terms (σ = +1) and the case of two opposite-sign kinetic terms (σ = −1).We refer to σ = +1 as the "PP" case, denoting two positive kinetic-energy DoFs, and to σ = −1 as the "PG" case, where a positive kinetic-energy DoF is coupled to a negative kinetic-energy (ghost) DoF.
Starting with a PP model, and applying a complex canonical transformation (cf.Eq. ( 2)) on the variable y, we obtain a new model of PG type with Hamiltonian 1 2 p 2 x −1 2 p 2 y + V (x, ±iy) which can be real depending on the V (x, y) one starts out with.
In the following, it will be sometimes convenient to split the potential V (x, y) into three parts as follows for a PG system with σ = −1, such that V = V P + V G + V int .We refer to V P (x) and V G (y) as the decoupled potentials, to V int (x, y) as the interaction potential, and to V (x, y) as the full potential.This makes sense, in particular, if the structure of the potential V (x, y) is such that y = 0 is a solution of y-equation of motion for arbitrary x(t) and analogously x = 0 is a solution of the x-equation of motion for arbitrary y(t).If this property holds 1 one can separate self-interactions from cross-coupling interactions.
More specifically, for what concerns the analytic part of this work, our starting point will be integrable PP systems with a time independent Hamiltonian H of the form in Eq. ( 6) and with an extra constant of motion I(ξ) such that The study of such systems has been pioneered by Liouville [38,39] who showed that if an N dimensional Hamiltonian system possesses N functionally independent and Poisson-commuting constants of motion, then the motion is integrable (by quadrature).In principle, there can be more than N (i.e., up to (2N − 1)) functionally independent constants of motion.If so, the system is referred to as superintegrable.Some Hamiltonian systems are trivially integrable.For instance, the canonical Hamiltonian of N uncoupled harmonic oscillators, see Eq. ( 5), possesses N conserved quantities: the individual Hamiltonians of each oscillator.For other Hamiltonian systems, integrability is, of course, nontrivial and takes more effort to be revealed.
One way to systematically search for constants of motion is by direct methods.Direct methods (i) make an ansatz for I(ξ) -typically a polynomial in momenta, dressed with undetermined functions as prefactorsand then (ii) classify solutions to the resulting (highly overdetermined) set of partial differential equations (PDEs) implied by {I, H} PB = 0.When restricting to two particle systems, such direct searches for constants of motion (and the corresponding integrable potentials) have been performed at fixed polynomial order in the momenta.Up to quadratic order, the classification has been pioneered in [40], formalized in [41], and subsequently completed in [42][43][44][45][46][47][48][49].We summarize the complete set of resulting pairs of integrable potentials V (x, y) and constants of motion in App.(A).At cubic order in the momenta and beyond, no full classification is known, cf.[49] for a comprehensive review of partial results and for references to related methods.
As stressed above, the on-shell constancy of some quantity I is preserved by the action of C ± x n as a mere consequence of the algebraic nature of equation Eq. ( 8) and the canonical nature of the transformation C ± x n .Hence, acting on a PP integrable system with a transformation C ± x n , we obtain a PG integrable one.Upon this action, the Hamiltonian can, in general, become complex.However, an inspection of the potentials given in App.(A) shows that there exist large classes of PP integrable systems yielding real-valued PG integrable systems upon application of the complex canonical transformation C ± y (or C ± x ).This is the case, e.g., for class 1 of App.(A).In the next section, we use this real-valued PP ⇔ PG pair to analytically discuss the conditions for bounded motion in the presence of a ghost.
In passing, we also note that there can be realvalued integrable PG systems, for which the PP side is complex-valued, and which may thus have not yet received much attention.An integrable example of this kind is given by class 5 in App.(A).In this case the application of the complex canonical transformation C ± y (or C ± x ) transforms a complex V to a real one.
To avoid confusion, we stress that complex canonical transformations C ± x n preserve the physics only when accompanied by the corresponding rotations in the complex planes of the initial data.Otherwise, the transformed theory is physically different and equivalent to probing the original theory with purely imaginary initial data.Therefore, the dynamical properties of such "ghostified" theories can be completely different from those of the original theory when all quantities considered are real, which we will assume in our PG models.
Let us illustrate this by a simple example of a theory with one real degree of freedom x (and real p) Clearly, for arbitrary, we stress real, initial data, the solution perpetually evolves in a finite region of phase space bounded by the potential well cosh x and the positive kinetic energy.This solution is a closed trajectory in the phase space (x, p), see the upper panel of Fig. (1).Now let us perform the canonical transformation C − x , so that p = ip and x = −ix, transforming Eq. ( 9) to the new We show the phase space (x, p) of the original Hamiltonian in Eq. ( 9) (upper panel) as well as the phase space (x, p) of the transformed Hamiltonian in Eq. ( 10) (lower panel).In the ghostified case (lower panel), all trajectories with initial momentum |p0| > 1 (cf.red horizontal lines) evolve over the potential barrier and run away to plus or minus infinity in x, irrespective of the initial value of x0.Two such runaway trajectories are depicted in blue.In contrast, all phase space trajectories of Eq. ( 9) (upper panel) are closed curves.
theory of real variables (x, p) and Hamiltonian, i.e., We could consider the solution (x (t) , p (t)) for the system in Eq. (10), with purely imaginary initial data (x 0 , p0 ) = (ix 0 , −ip 0 ) to obtain just the real solution (x (t) , p (t)) = (−ix (t) , ip (t)) of Eq. ( 9) with the real initial data (x 0 , p 0 ).In particular, H = H for all the solutions which are related in this way.However, a solution of Eq. (10) with real initial data (x 0 , p0 ) corresponds to a purely imaginary solution of Eq. ( 9) with initial data (−ix 0 , ip 0 ).Thus, purely real solutions of Eq. ( 9) are very different from purely real solutions of Eq. ( 10) even though both systems are related by a simple canonical transformation.In particular, for arbitrary real initial coordinates x0 and initial momenta p0 > 1 (or p0 < −1), the trajectory of Eq. ( 10) evolves over the potential barrier and is unbounded, with a runaway of x(t), see the right panel of Fig. (1).As we mentioned earlier, there are no such real-valued trajectories in Eq. (9).Clearly, in the phase space (x, p) of Eq. ( 10) such real-valued trajectories are not closed, whereas all real-valued trajectories in (x, p) of Eq. ( 9) are closed, thus even the topology of the phase space trajectories is changed.In this regard it is worth mentioning that with real (x, p) there are not only an infinite number of new stable equilibrium points, but also new unstable equilibrium points in-between the former, see Fig. (1).While the above refers to local (Lyapunov) stability, the "ghostified" real theory in Eq. ( 9) is globally (Lagrange) unstable, in contrast to the original real theory Eq. (9).Two of the respective runaway trajectories are marked as thick blue lines in Fig. (1).

III. THE INTEGRABLE LIOUVILLE MODEL & A PROOF OF BOUNDED MOTION
In the following, we consider an integrable class of theories, cf.class 1 in App.(A), for which we can rigorously establish boundedness of motion.Specifically, we focus on with where f and g are arbitrary functions of the coordinates (u, v).The coordinates (u, v) are related to (x, y) via where c is an arbitrary constant.Besides the Hamiltonian, this theory exhibits a second constant of motion, i.e., where we have (for later purpose) defined the momentum-independent part of the constant of motion V = I LV | px,y→0 .Although this model was initially studied in the PP case (i.e. with σ = +1) the discussion of the previous section shows that it stays integrable also for σ = −1, i.e., in the PG case.The corresponding extra integral of motion (i.e.Eq. ( 16) with σ = −1) is obtained from the PP one (i.e.Eq. ( 16) with σ = +1), as explained in the previous section.For σ = 1 and c ⩾ 0, u and v define so-called elliptic coordinates, i.e., curves of constant u ∈ ( √ c, ∞) trace out ellipses and curves of constant v ∈ (0, √ c) trace out hyperbolas in the (x, y) plane.For σ = −1 and c ⩽ 0, curves of constant u ∈ (0, ∞) and curves of constant −i v ∈ ( √ c, ∞) both trace out hyperbolas in the (x, y) plane, cf.Fig. (2).
We note that, even though this class of integrable systems was derived, and is usually presented, with an arbitrary constant c, it is only the sign of this constant which influences the dynamics.Indeed, for c ̸ = 0 one can first perform a canonical transformation Then one can complete it with rescaling of time These transformations completely absorb |c| into the arbitrariness of f (u) and g(v), except for the limit cases c → ∞ or c → 0. For instance, the latter limit transforms the Liouville model to a subclass of class 2, see Eq. (A9), with f = 0.The above integrable Hamiltonian in Eq. ( 11) follows as a subclass of a system first investigated by Liouville [50] (in the PP case).Liouville starts out from a Hamiltonian with non-standard kinetic terms, i.e., where (α, β) are yet another set of coordinates, (p α , p β ) are the respective canonical momenta, and f (α), F (α), g(β), G(β) are arbitrary functions.Following Liouville, we found that Poisson-commutes with the Hamiltonian H LV and thus defines a second constant of motion2 .Hence, the model is integrable, irrespective of the choice of f (α), F (α), g(β), and G(β).
The coordinates (α, β) are related to (u, v) via In coordinates (u, v), the Liouville Hamiltonian H LV and constant of motion I LV read For the specific choice of F (u) = u 2 and G(v) = v 2 , and with the transformation from (u, v) to (x, y), these reduce to the Hamiltonian and constant of motion previously given in Eqs. ( 11) and ( 16) with σ = +1.For more general F (u) and G(v), the transformations to (x, y) coordinates does not bring the Hamiltonian to a form with standard kinetic terms.In the following, we will thus restrict to the case of

A. Proof of boundedness of motion
We show here that the PG Liouville model with σ = −1 (in the form of Eq. (11) and with the additional constant of motion in Eq. ( 16)) is such that -under certain conditions on f (u) and g(v) detailed below -the phasespace motion is always bounded, irrespective of the choice of initial conditions.This applies to the model considered in [16], which is a special class of PG Liouville models, but also extends the result of [16] to a much larger class of ghostly oscillator interacting with a positive energy oscillator where the motion can hence be proven analytically not to run away.
We start by defining three further helpful coordinate combinations with which we may re-express W = √ W 2 + w and write We choose here and henceforth a negative c < 0 such that the momentum-dependent part of the constant of motion (cf. the first two terms in Eq. ( 24)) are positive.Furthermore, a negative c < 0 implies and while the above v 2 (given in Eq. ( 28)) is a real expression, the so defined v is purely imaginary.In order to manipulate real and positive variables, we define the real and positive quantities ṽ and c as such that v 2 = −ṽ 2 , and thus we have all these quantities being real, as well as, by choice, the functions f and g.
We now prove that phase-space motion is bounded if with constants f 0 , g 0 ∈ R; and (ii) at large |u| and |v|, these lower bounds sharpen to with positive constants F 0 , G 0 ∈ R + as well as ζ > 2 and η > 2 (and the factor 4 is just here for later convenience).
The first step of our proof is the introduction of a new first integral J LV defined as (using Eqs. ( 31), ( 34) and ( 36)) where and γ ± are obtained to be (using Eqs. ( 26), ( 28) and ( 29)) Using Eq. ( 26), we get This, using the condition (i) above implies that where the equality may hold if g 0 = 0 or/and f 0 = 0, and hence, using Eq. ( 42), we get that showing that U is bounded below.As J LV defined in Eq. ( 41) is conserved, this implies that |xp y + yp x |, |p x | and |p y | are bounded on shell.This is the first step of our proof.As it should be clear, this first step puts forward an analogy between a conserved and bounded from below Hamiltonian and the conserved J LV , where all the "kinetic energies" are positive definite and the "potential energy" U (u, v) has just been shown to be bounded from below.
To proceed further, we use again this analogy, showing that conditions (ii) imply that U (x, y) is growing without bound at large radius R = x 2 + y 2 , which excludes runaways of the dynamical variables x and y.To do so we introduce polar coordinates (R, ϕ) such that x = R cos ϕ and y = R sin ϕ, yielding r 2 = x 2 − y 2 = R 2 κ where we have defined κ = cos 2ϕ.As a consequence, one has As we are interested in large R, we can assume that implying, in particular, that the last term in the squareroot in Eq. ( 48) is much smaller than the second one.We also have from Eq. ( 26) Using this we conclude that one has either Distinguishing the above, and moreover between or we will now show that, in all the respective four cases, U grows at least as fast as a suitable power of R.
Let us first assume that Eq. ( 53) holds.In this case we have This implies a lower bound on γ ± given by3 which holds for either sign of κ (and including κ = 0).Moreover, we have from Eqs. ( 51) and ( 52) that either In the first subcase, we can use the bounds Eqs. ( 39), ( 45) and ( 56) to get that at large R In the second subcase, we get similarly Hence, at least one of these two bounds holds whenever Eq. ( 53) does.
We then turn to the second case in Eq. ( 54), where obviously κ ̸ = 0, and we can hence rewrite Eq. ( 48) as with ξ given by The inequalities in Eqs. ( 49) and ( 54) (and the assumption c > 0) imply in turn that4 Let us then first assume that κ > 0. In this case we can obtain a useful lower bound on γ − noticing that, for ξ in the range in Eq. ( 61), one has Using this, we conclude that, for positive κ and whenever Eq. ( 54) holds, we have We also note, following the same logic, that γ + obeys exactly the same bound for negative κ whenever Eq. ( 54) holds.We can then, as previously, distinguish between the two subcases in Eqs. ( 51) and ( 52) which implies that either u 2 ≥ |κ|R 2 /2 or ṽ2 ≥ |κ|R 2 /2.In the first subcase, and for positive κ, we can use the bounds in Eqs. ( 39) and (63) to get that at large R where to go from the first line to the second line we used once more Eq.( 54) and to go from the second line to the third line we used the fact that 1 ≥ κ > 0. Hence, using this last inequality as well as Eq. ( 45), we get that, for positive κ and whenever Eq. ( 54) and u 2 ≥ κR 2 /2 hold, we have Let us then turn to discuss the subcase where ṽ2 > |κ|R 2 /2.In this subcase, still assuming that κ > 0 we notice that γ + ≥ 1 and hence we get simply using Eqs.( 40) and ( 46) that where to go from the first line to the second line we used again Eq. ( 54).Similar bounds can be obtained for negative κ, mutatis mutandis.Indeed, we find in this case, following the same logic, one has either or Hence, we conclude from the considerations of Eqs. ( 57), ( 58) and ( 64) to (67) that the hypotheses (i) and (ii) above are enough to conclude that This concludes our proof: in all cases, the "potential energy" U grows without bound with R implying that |x| and |y| stay bounded on shell.
To summarize, we have shown that the phase-space motion stays bounded at all times, irrespective of the initial conditions, for functions f and g which (i) are bounded below, and (ii) at large values of their arguments, obey respectively Eq. ( 39) (where F 0 > 0 and ζ > 2) and Eq. ( 40) (where G 0 > 0 and η > 2).We stress that the behaviour of the functions f and g, entering in condition (ii), only needs to be specified at large |u| and |v|, implying, in particular, that the lower bounds f 0 and g 0 can be negative.
The above proof can also be extended to the case where, starting from the Liouville PP model, one ghostifies x by applying a canonical transformation C ± x .The obtained models are in fact identical (up to a global sign in the Hamiltonian) to the previously considered PG models obtained by ghostifying y from the PP model via a canonical transformation C ± y .As a result, the two kinds of models with a ghost are stable under the same conditions (i) and (ii).This is shown explicitly in App.(C).

B. A first application of the proof
As a first application of the above proof, we consider potentials V LV (x, y) generated by functions f and g which are polynomials in u 2 and v 2 respectively, i.e. of the form with and the potential V LV (x, y) is given by where W , W and w are defined in the previous subsection.Note that, in the above, the first term on the righthand side is a mere constant and, as such, can be omitted.Moreover, the second term proportional to W is a quadratic term for x and y (plus an irrelevant constant proportional to c).Choosing the functions f and g defined by Eqs. ( 70) and (71) fulfill conditions (i) in Eqs. ( 37) and ( 38) and (ii) in Eqs. ( 39) and ( 40 and c = 1 yields exactly the main model of reference [16] (whose Hamiltonian is given in equation ( 1) of this reference).The more general choice of a strictly negative c < 0, strictly positive and equal C 2 = D 2 > 0, and C 1 ≥ D 1 (and, moreover, setting D 1 = −C 1 to remove an overall constant) gives the largest set of models introduced in [16] (i.e., in Eq. ( 20) and ( 21) of this reference).The parameters a, b, c and d of [16] can be read off from the above as where the left hand side of the above refers to the notation of [16] (and we stress, in particular, that c on the left hand side of the third equality is introduced in [16] and is not the same as the c of the present work).We see that the conditions Eqs. ( 73) and (75) imply that a < 0 and d > 0, and the present work provides a stability proof that was not given in [16].Note further that the condition on c given in [16] (i.e.c ≤ 0) can, in fact, be relaxed and still yields stable motion.

C. A polynomial subclass of the Liouville model: towards more general conditions for stable motion
The Liouville model at the heart of this work, and defined in Eqs. ( 12) and (15), possess an interesting subclass for which the potential V LV (x, y) is polynomial in x and y.This subclass is generated by choosing with N ∈ N (and still C n ∈ R), i.e., the same polynomial form for f and g.In App.(D), we provide a proof that the resulting potential V LV (x, y), given by is indeed polynomial in x and y.To our knowledge, it was not noticed before that this particular choice leads to polynomials, and hence, for both the PP and PG cases, to integrable systems with polynomial interactions whose specific form can be found below and in the App.(D).
The explicit form of V LV (x, y), in the PG case, and for the simplest non trivial case of interest here, is obtained from the expression in App.(D), where one trades the constants C 2 and C 3 for quadratic terms and chooses the constant C 1 to eliminate an irrelevant constant term.One obtains, ordering the interactions by powers of x and y, Except for the degenerate choice of ω x = ω y (for which the N = 3 case reduces to the N = 2 case), the N = 3 case does not lead to fully stable motion, as conditions (i) in Eqs.(37) and ( 38) and (ii) in Eqs.(39) and (40) of Subsec.(III A) are not fulfilled 5 .Indeed, either one chooses C 3 positive and for large negative v 2 , g(v) diverges to −∞, or one chooses C 3 negative and for large positive u 2 , f (u) diverges to −∞.An example of two near-by sets of initial conditions -one stable and one unstable -is visualized in Fig. (4).
In contrast, the N = 4 case leads to bounded motion since, for C 4 > 0, conditions (i) in Eqs. ( 37) and ( 38) and (ii) in Eqs. ( 39) and ( 40 6), we show an example of the motion for which the ghost has a tachyonic term, i.e., for ω 2 x = 5, but ω 2 y = −5.Finally, in Fig. ( 7), we show an example of the motion with a significant frequency hierarchy, i.e., for ω 2 x = 1 and ω 2 y = 100.In all cases, C 4 = 1.The "potential" in Eq. ( 42) that stabilizes the motion is bounded from below and grows without bound at large x and y for C 4 > 0. It is given here by where we rewrote Eq. (D12) using c > 0. Note that both the above U LV and V LV (x, y) are symmetric with respect to separate reflections x → −x , and y → −y , as they depend only on x and y via their squares x 2 and y 2 .The potential V LV (x, y) has up to power six selfinteractions and PG cross-couplings between the usual DoF x and the ghost y.These crosscoupling terms vanish when one of the coordinates is zero.However, elsewhere these cross-couplings cannot be neglected in comparison with the self-interactions.Without cross-couplings the positive energy degree of freedom x 5 Even though our proof in Subsec.(III A) only provides sufficient not necessary conditions for the boundedness of motion, on physical grounds, we expect that integrals of motion unbounded from below and from above cannot bound the motion.
has a self-interaction potential while the the ghost self-interaction potential is (84) Without interactions one could flip the sign in front of the ghost's Hamiltonian, so that uncoupled ghost effectively experiences inverted self-interaction potential (−V G (y)).Note that (−V G (y)) can be obtained from (V P (x)) by changing x to y everywhere, including indices.
The next non-trivial polynomial and stable PG theory is obtained for N = 6 (and positive C 6 ) and the corresponding potential is given in App.(D).
More generally, all polynomial theories with even N (and positive C N ) are stable according to the above proof.This becomes apparent when looking at large In writing the second expression, we arrange minus signs such as to account for −v 2 > 0.
For odd N , the proof conditions cannot be fulfilled since (for any choice of C N ) either f (u) or g(v) is unbounded below.While islands of stable initial conditions are possible (cf.left-hand panel in Fig. (4)), we expect that there are always some initial conditions for which divergent behaviour remains (cf.right-hand panel in Fig. (4)).
So, how do the proof conditions relate to the behavior at large x 2 ≫ c and x 2 ≫ |C n | (as well as large y 2 ≫ c and y 2 ≫ |C n |)? Taking both limits, Together with the proof conditions, i.e., even N and C N > 0, this implies: • The decoupled potential V P (x) (see Eqs. ( 7) and ( 83)) is bounded below and unbounded above and thus serves as a stable potential for the decoupled motion of the x mode.
• Similarly, V G (y) (see Eqs. ( 7) and ( 84)) is bounded above and unbounded below and thus serves as a stable potential for the decoupled motion of the y mode.
t 0 100 200 300 We show trajectories (with initial conditions tuned to the boundary of an apparent island of stability) for N = 3, i.e., for ω 2 x = 1, ω 2 y = −1, C4 = 0, and c = 1, in Eq. ( 79).We superimpose a contour plot of the potential (cf.left-hand panel in Fig. (3)) where increasingly light shades indicate increasingly large values of the potential.In the left-hand panel, we show a bounded trajectory.In the right-hand panel, we show an adjacent unstable trajectory.To be specific, and in case the reader wants to reproduce our results, we have chosen x(t = 0) = 1, x ′ (t = 0) = 1, y(t = 0) = A, and y ′ (t = 0) = A, for which the difference between bounded and unbounded motions is obtained by choosing respectively A ≈ 0.468 in the left-hand panel and A ≈ 0.469 in the right-hand panel.These values hence also locate the boundary between the two kinds of motions considered here.
• At sufficiently large 6 x 2 and y 2 , the interaction potential V (x, y) is bounded by the decoupled potential, i.e., • As we will see, the special case of |x| = |y|, for which the full and the decoupled potential agree at sufficiently large x 2 and y 2 , seems to require an additional criterion.
To elucidate the latter subtlety further, we take a look at another class of integrable polynomial potentials.
D. An unstable polynomial subclass: refining more general conditions for stable motion To elucidate the possibility of stable motion in the Liouville class further, we can compare it to the motion in another class of integrable two-particle Hamiltonians.As explained in Sec.(II) real-valued PG Hamiltonians can sometimes also be obtained from complex-valued PP Hamiltonians.This turns out to be the case for another class of integrable two-particle Hamiltonians, closely related to the Liouville model, cf.App.(A).The details of this integrable model are discussed in App.(E) and, following [49], we refer to it as class 3 from here on.
In fact, just as the integrable Liouville model, the integrable class 3 contains a polynomial subclass which is analogously obtained, cf.App.(E).Introducing frequencies, removing a constant shift (as in Eq. ( 79)), and specifying to equal frequencies (i.e., ω x = ω y = ω), the first few terms can be written as The general expression with unequal frequencies is given in App.(E).For C 4 = 0, the potential V (class 3) is (depending on the other parameter choices) either bounded from below or bounded from above.More generally, just as for the Liouville polynomial subclass, each odd N case is thus unstable because one of the decoupled potentials is necessarily unstable.
For C 4 > 0, cf. also Fig. (8), the potential V fulfills all the conditions stated around Eq. ( 87): Both of the decoupled potentials are stable and the full potential is bounded by the decoupled potential.Nevertheless, we find that this case exhibits runaways with polynomial growth rate and is thus not stable.We thus conclude that the conditions stated around Eq. (87) are not sufficient to guarantee the absence of runaway behaviour.
We first look at the case (ω = 0), for which we observe that the runaway behaviour occurs either along the x−y = 0 or along the x+y = 0 direction.In fact, the potential in Fig. (8) exhibits sets of saddle points (marked with thick lines) extending out to infinity along flat directions of the potential.As a consequence of this, the local potential around x ± y = 0 is no longer dominated by the respective decoupled potentials.
This can be seen from the Hessian matrix 7 , which vanishes at x = ±y, i.e., Hess V (4) In contrast, the Hessian matrix of the decoupled potential (at x = ±y) is non-vanishing, i.e., Hess V (class 3) (x, 0) + V (90) 7 For the case at hand, the Hessian matrix is defined by Hess V Along the saddle surfaces, the local potential V (4) (class 3) does, therefore, not resemble the decoupled potential.In fact, it does not resemble any decoupled potential.Put differently, there exist regions in which decoupled interactions no longer dominate the potential at large phase-space variables.
For V (4) (class 3) with non-vanishing but equal frequencies (and still for C 4 > 0), the issue of runaways is lifted along one of the two directions.Without loss of generality, picking the '−' sign in the last term in Eq. (88), the x = −y direction is stable but the x = y direction remains unstable.Once more, this can be understood in terms of the Hessian matrix.For x = −y, the Hessian matrix now reads This agrees with the Hessian matrix of a free theory, to be specific, the free theory with V (x, y) = ω 2 2 (x 2 − y 2 ) + 2 C 4 (x 4 − y 4 ).In contrast, for x = y, we find Hess V (4) for which the behaviour at large phase-space variables still has no free-theory equivalent.In agreement with the above, we find that, with equal frequencies, polynomial runaway behaviour no longer occurs along the x = −y but still occurs along the x = y direction.
Numerically, we observe that the runaways can be fully removed by the addition of decoupled interactions of the form x 4 − y 4 .(Numerics are necessary since the addition of these terms breaks integrability.)With this addition, the Hessian matrix has the same structure as for the Liouville potential V (4) LV (x, y).In both cases, it does still not agree with the Hessian matrix of the respective decoupled potential.However, one can identify a decoupled potential for which the Hessian matrices (at sufficiently large x = ±y) agree.We conclude that it is sufficient (at least for the specific case at hand) that the Hessian matrix at large |x| = |y| agrees with some decoupled potential, not necessarily with the decoupled potential of the respective model.
In passing, we mention that similar polynomial runaway behaviour also occurs in an integrable model previously obtained (from a supersymmetric construction) by Smilga et al. [21,51].Its Hamiltonian is given by where λ and ω are constants.It follows as a degenerate limit of Eq. (88), i.e., when picking the + sign in the last term of Eq. ( 88) and taking C 4 → 0 as well as c → ∞ such that C 4 c = −λ/4ω.(Alternatively, Eq. ( 93) can also be obtained as a specific case 8 of class 8 in App.(A).)The motion of this theory does not stay bounded (even though the ghost there is argued to be "benign" due to the slow (polynomial) growth-rate of the runaway solutions found numerically [21,51]).
We conclude that all of the discussed cases can be understood if the following criterion is added to the list in Subsec.(III C): • For |x| = |y| (or more generally, for all (x, y) for which V (x, y) = V P (x) + V G (y)), the Hessian matrix at large phase-space variables approaches the Hessian matrix of some free theory. 8To be specific, the Hamiltonian in Eq. ( 93) can be obtained from the integrable class 8 in Eq. (A15) by choosing g = 0 and f ′′ = 1 2 ω 2 + λ ω (x + iy) 2 .This results in the potential (for a PP theory) given by V = (x 2 + y 2 )( 12 ω 2 + λ 4ω (x + iy) 2 ).We then obtain the theory in Eq. ( 93) by applying a canonical transformation C − y .
Having understood these more subtle cases, we conclude: By translating the proof conditions to conditions on the potential V (x, y), we understand that stable motion seems to be possible because stable decoupled potentials V P (x) and V G (y) dominate all interaction terms at sufficiently large x 2 and y 2 .
In the following Sec.(IV), we will use this insight to formulate general conditions for stability.Before doing so, we comment on the relation to local Lyapunov stability.

E. Lagrange vs Lyapunov Stability
In the previous subsections, we were concerned with global boundedness of motion which corresponds to the stability in the sense of Lagrange.However, one can be also interested in local stability in the sense of Lyapunov, see, e.g., [37].In particular, this is crucial to specify the nature of equilibrium points which are fixed points of the system of equations of motion and can be considered as vacua.The two notions of stability are not equivalent.Indeed, a Lagrange stable system may not have locally stable solutions, while a system with Lyapunov stable equilibrium points may have other solutions with runaways invalidating stability in the sense of Lagrange.Here we will show that the Liouville class of systems in Eq. ( 11), describing ghosts interacting with usual degrees of freedom (σ = −1), contains plenty of theories with Lyapunov stable equilibrium points.More specifically, we focus here first on systems in this class with functions f (u) and g(v) which are respectively analytic functions of u 2 and ṽ2 around the origin (x, y) = (0, 0).These systems contain the polynomial case Eqs. ( 76) and (77) as a subcase.
For such theories, the origin (x, y) = (0, 0) is always an equilibrium point.Indeed, at the origin one has u 2 = 0 and ṽ2 = c, while because all partial derivatives around the origin and obtain corresponding frequencies Note that ωx and ωy are just defined expanding around the origin and are hence conceptually different from unbarred frequencies defined previously, hence the different notation.Of course they could be the same depending on the potential.Naive stability by the linearized approximation requires that both ω2 x and ω2 y are positive.In particular, adding the two equations above, we see that stability requires that Let us now confirm that the conditions ω2 x > 0 and ω2 y > 0 are sufficient for the Lyapunov stability of the equilibrium point in the origin.Consider the first integral J LV given by Eq. ( 41).Expanding J LV around the origin and around vanishing momenta one obtains where ellipsis stands for higher order terms.Clearly, for ω2 x > 0 and ω2 y > 0 the first integral has a minimum at the equilibrium point in the origin.Thus, J LV −(g (c)+f (0)) is a Lyapunov function.Indeed, it is conserved, vanishing at the equilibrium solution and is positive definite in an open neighborhood around it (see [37]).Thus, the origin is Lyapunov stable.Here, the conserved quantity J LV plays the analog of usual energy for establishing stability.Crucially, the "kinetic energies" are all positive definite in J LV .Moreover, U(x, y) -the "potential part" of the "energy" -has an isolated minimum at the origin.Thus, exchanging the usual potential energy with U(x, y), one meets the conditions of the classical Lagrange theorem stating that the minimum of the potential energy corresponds to a stable equilibrium configuration.In contrast, it is interesting to note that this stable equilibrium point is a saddle point for the potential V LV (x, y) as It also appears above that ω2 x > 0 and ω2 y > 0 do guarantee Lyapunov stability for the usual PP systems due to the Lagrange theorem.The Lyapunov stability of the origin, implies, by continuity, the same stability for solutions in a small neighborhood.It is important to repeat that the conditions on f (u) and g(v) required for the Lyapunov stability of the origin (positivity of (98) and (99)) can be easily violated in systems which are globally stable according to Lagrange.And vice versa: a Lagrange unstable system may have a Lyapunov stable origin.This can be well illustrated by the case of V (3) LV PG system given by Eq. (78), see also Eq. (D4).
As we discussed before, this system is not stable from the point of view of Lagrange, however the origin is Lyapunov stable provided ω 2 x and ω 2 y chosen in Eq. (D6) are positive.On the other hand, as we discussed before, V LV from Eq. ( 79) is Lagrange stable for C 4 > 0 and c > 0, but can have a Lyapunov unstable origin if one of ω 2 x , ω 2 y chosen in Eq. (D6) is negative.To finish, let us here mention that the Lyapunov stability of the origin equilibrium point for a distinct class of interacting ghosty systems has previously been demonstrated in [28] using different methods.Moreover, in [16] the Lyapunov stability of the origin has been proven for a particular PG system from the same Liouville class in Eq. ( 11) using methods similar to those employed in the current paper.

All Vacua and their Lyapunov Stability for the Polynomial Subclass of the Liouville Model
The PG interacting system with polynomial potential in Eqs. ( 79) and (D7) belongs to the class discussed above and deserves to be analyzed in more detail.In Fig. (5), Fig. (6) and Fig. (7) we have already presented numerical solutions of this system for various parameters.In this subsubsection, we find all equilibrium points, or fixed points of the evolution (x e , y e ), i.e., We then determine whether these vacua are Lyapunov stable.Assuming non vanishing C 4 and c, the above two equations are equivalent to the vanishing of , whose solutions will be discussed below.We then shall demonstrate how PG interactions which are supposed to destabilize the system, instead introduce new stable vacua.Interestingly, the appearance of these new vacua leads to spontaneous breaking of the reflection symmetry in Eq. (81).To establish Lyapunov stability, we will use the result of App.(B): the extrema of V LV are also extrema of U (4) LV , see Eq. (B4).Then, similarly to the discussion above on the stability of the origin, the minimum of U (4) LV implies the Lyapunov stability of the vacuum, as J LV from Eq. ( 41) can be taken again as a Lyapunov function.
The first equilibrium point is the origin (x e , y e ) = (0, 0), which was discussed at the beginning of this section, and the Lyapunov stability just requires that ω 2 x > 0 and ω 2 y > 0 which we shall also assume to hold below.Note that this vacuum respects the reflection symmetry in Eq. (81) for both coordinates and is such that (see Eq. ( 80)) LV (x e , y e ) = 0 . (104) One can find a second type of equilibrium points where only one of the coordinates vanishes but not the other one.Indeed, x e = 0, or y e = 0, is obviously a solution to Eq. (103) (see the two above equations).When either x e = 0, or y e = 0, the stability along the direction where the equilibrium coordinate does not vanish can be analyzed looking at the potential in Eq. (84) or Eq. ( 83).It is seen there that, for sufficiently small C 4 , either the coefficient in front of x 4 in V P (x) or in front of y 4 in (−V G (y)) can become negative, indicating that the origin is not the only minimum of V P (x) or of (−V G (y)).In that case, the vacuum spontaneously breaks the reflection symmetry in Eq. ( 81) with respect to this coordinate.Let us assume for concreteness that ω x < ω y , so that we take y e = 0 as then implied by the positivity of the coefficient in front of y 4 in Eq. ( 84).The opposite case is obtained from all formulas with the replacement x ↔ y.Looking for nonvanishing solutions of ∂ x V P (x) = 0, yields the following quartic equation which has real non-vanishing solutions x 2 for The last condition is easy to fulfill for a sufficiently small C 4 , i.e., for (see also Footnote 11) The solution of Eq. ( 105) is (108) where we denoted10 To obtain a real x b(±) , one has to require that the first term in the right hand side of Eq. ( 108) is positive 11which translates to TABLE I. We list all possible minima of U (4) LV (x, y), from Eq. ( 80) where we assumed that c > 0, C4 > 0 and Ω 2 x/y > 0. In the last column, we list conditions for the existence of such minima.These minima are Lyapunov stable equilibrium configurations or stable vacua of the theory.The origin respects the reflection symmetry in Eq. ( 81).The "decoupled-x" minima are given by Eq. ( 108) and break the reflection symmetry x → −x.The "decoupled-y" minima are given by Eq. ( 108), with exchange x ←→ y, and break the reflection symmetry y → −y.The "coupled" minima are given by Eqs. ( 112) and ( 113) and break both of the reflection symmetries in Eq. ( 81).FIG. 9.In the left-hand panel, using the rescaled frequencies in Eq. ( 109), we visualize the conditions for the existence of different stable vacua listed in Tab.(I).In the white region, the origin is the only Lyapunov-stable vacuum.In the yellow stripe above (below), the stable origin is supplemented by the two Lyapunov-stable "decoupled-x" ("decoupled-y") vacua.Further up and/or to the right, in the orange region, the stable origin is supplemented by the four Lyapunov-stable coupled vacua.In the other three panels, we specify to C4 = 1 as well as c = 1 and show the three indicated examples for the respective possible cases.We super-impose the latter on contour plots of log(|U x , Ω 2 y ) = (0.5, 7.5) (i.e., the point labeled "I" in the orange region of the left-hand panel), (Ω 2 x , Ω 2 y ) = (0.5, 3.75) (i.e., the point labeled "II" in the yellow region of the left-hand panel) and (Ω 2 x , Ω 2 y ) = (0.5, 0.5) (i.e., the point labeled "III" in the white region of the left-hand panel).Larger (smaller) dots indicate minima (other real-valued extrema) of U (4) LV (x, y) and the colours white, yellow, and orange indicate origin, decoupled extrema, and coupled extrema, respectively.Even though, for illustrative purposes, we have not chosen the parameters to preserve the ratio of Ω's but rather to fit the minima at the same region of (x, y) space, it is clear that, for fixed values of c and ωy, ωx, the vacuum configurations will change from left to right with increasing strength of C4.
Note that this inequality in Eq. ( 110) is weaker than Eq.(107), so that, when the latter holds, both extrema x 2 b(±) are real and in total there are four solutions.Clearly, x 2 b(−) corresponds to a local maximum of V P (x) while x 2 b(+) corresponds to a local minimum and hence a stable direction along the x ̸ = 0 direction.Now we can investigate when such points ±x b(+) , 0 are minima of U (4) LV given by Eq. ( 80).First one notices, that at these points LV ±x b(+) , 0 = 0, while the other second derivatives are LV .The parameters are chosen as c = 1, ωx = 1, and ωy = 1.7, while C4 is 1.5 times smaller than the lower bound in Eq. ( 114) so that C4 ≈ 0.103.This case has different numerical values of constants, but its structure of vacua corresponds to case I in Fig. (9).
where we denoted One can deduce that both second derivatives are positive provided (111) Note that for the assumed ω y > ω x the lower bound above is always smaller than the upper bound.Thus there is always a space for C 4 between them 12 .Moreover, the upper bound is just what is required for the very existence of the equilibrium point, see Eq. ( 107).Hence, we have just shown that, for C 4 satisfying the bound above, the stable points ±x b(+) , 0 correspond to local minima of U (4) LV (see Fig. (11)) and these vacua are Lyapunov stable, as they actually also are turning off interactions with the ghost y.For the opposite hierarchy between frequencies, analogous equilibrium points are 0, ±y b(+) and the whole analysis and formulas are applicable after exchange x ←→ y.Interestingly, for 12 Note that one can reinterpret this inequality as x for rescaled frequencies in Eq. ( 109), see Fig. (9).
C 4 smaller than the lower bound in Eq. ( 111), the equilibrium points ±x b(+) , 0 are not local minima of LV so that one cannot warrant stability.This is a strong hint that these equilibrium configurations then become unstable.On the other hand, one can check that ±x b(−) , 0 never become minima of U (4) LV showing that these equilibrium configurations remain always unstable.Of course the same logic is applicable after exchange x ←→ y.
Finally, the last category of equilibrium points appears for such small coupling constant C 4 that the PG-interactions in Eq. ( 82) introduce new nontrivial equilibrium points.These vacua spontaneously break both reflections in Eq. ( 81) and are solutions of the system in Eq. ( 103) By subtracting these equations from each other one obtains so that now one can substitute y 2 − x 2 back into both equations to obtain, using the same definitions as in Eq. (109), and similarly The right hand side of both equations are strictly positive provided that 13 (as easily seen using Eq. ( 109)) which guarantees that these nontrivial equilibrium points exist.Further, this inequality also implies that LV > 0 along with the determinant of the Hessian Hess(U (4) LV ) > 0, see Eq. (B13).Thus, all the four equilibrium points 14 (x I , y I ) , (−x I , y I ), (x I , −y I ), (−x I , −y I ) are Lyapunov 13 These inequalities trivially transform to inequalities on rescaled frequencies, e.g., Ω 2 x − Ω 2 y 2 > 4Ω 2 y , see Fig. (9). 14Here we assumed that the square root is taken as a principal square root returning non-negative numbers.
LV (x, y) around vacua, for C4 located in the interval (111), while on the right, we plot, in the logarithmic scales, the same "uplifted" effective "potential": log 1 + U stable, as the "potential" U LV has local minima there.One can check that Thus, comparing the above to Eq. ( 104), we see that the equilibium points above are in a sense true vacua as U LV (x I , y I ) < U LV (x e , y e ), see also Fig. (10).
We would like to stress that we found all possible equilibrium points, i.e. solutions to Eq. (103).Moreover, comparing, Eq. (111) with Eq. (114) one infers that the existence of stable equilibrium points (±x I , ±y I ) and ±x b(+) , 0 (or 0, ±y b(+) ) is mutually exclusive.Hence, PG-interactions in Eq. (82) can actually introduce new Lyapunov-stable equilibrium points which would be impossible if the ghost y had not interacted with the positive-energy DoF x.Thus, at least for discrete dynamical degrees of freedom, interactions with ghost can create new stable vacua.Finally, we would like to mention that it is the strength of the coupling constant C 4 which controls the symmetry breaking pattern for fixed c and ω x/y .The larger is the constant the more symmetry is there, see Fig. (9).

IV. NUMERICAL INVESTIGATION OF NONINTEGRABLE MODELS
We now relax the technical assumption of integrability.In Sections III C and III D, we have rephrased the conditions for stability of the proof of Subsec.(III A) (which are stated in terms of free functions, specific to the integrable model) to conditions in terms of the shape of the potential V (x, y).
This motivates us to investigate whether similar conditions lead more generally to stable motion and, in particular, when there is no integrability.
We find numerically that: (ℵ), a sufficient condition for the stable motion of a Hamiltonian point-particle system with coupled positive-and negative-sign kinetic terms for the respective coordinates x and y is that its full potential V (x, y), at large coordinates x and y, is sufficiently dominated everywhere by a potential of the form ṼP (x) + ṼG (y) where both ṼP and − ṼG are yielding by themselves stable motions along x and y.I.e.we need the interactions to be sufficiently suppressed with respect to a stable "decoupled" part at large values of the dynamical variables.
Note that these conditions are reminiscent of the Kolmogorov-Arnold-Moser (KAM) theorem (for a pedagogical discussion see [39,Appendix 8]), as the theory with the Hamiltonian just given by ṼP (x)+ ṼG (y), which dominates at large distance, would be integrable (as it has two separately conserved quantities: the respective energies of x and y), our statement is just that, it is enough for the interactions between x and y to be sufficiently subdomi-nant at large x and y to ensure stability.The conditions stated in Sections III C and III D are an explicit expression of the above statement for the specific integrable polynomial subclasses at hand.Beyond integrability, we still find that, whenever the above statement (ℵ) holds, all (at least, all numerically investigated) initial conditions avoid runaways (at least, up to the finite numerical evolution time).
Naturally, the lack of integrability also means that the analytical proof strategy in Sec.(III) will no longer apply.Hence, we need to resort to numerical evidence.The latter comes with the difficulty that -in absence of rigorous analytical proof -runaways may only be encountered at very late evolution times and/or only for a small subset of the initial conditions.In Subsec.(IV A), we introduce a specific numerical setup to alleviate these issues, at least to some degree.In Subsec.(IV B), we demonstrate that the numerics are able to confidently distinguish between stable and unstable potentials in the presented subclass of integrable and polynomial Liouville models.With this benchmark at hand, we apply the same numerical analysis to several nonintegrable models in Subsec.(IV C), i.e., to two classes of non-polynomial and to one class of polynomial potentials.This provides evidence that the above statement (ℵ) holds true, even in the absence of integrability.
We emphasize that all of our results suggest that the statement (ℵ) only needs to hold outside of a finite region in phase space.As a special case, it thus encompasses all interaction potentials which are localised in a finite phase-space region and supplemented by standard quadratic terms V P = 1 2 ω 2 x x 2 and V G = − 1 2 ω 2 y y 2 with ω x , ω y ∈ R + .Two explicit examples of this are presented in Subsec.(IV C).Note that this contrasts with the usual conditions used in the KAM theorem which uses perturbations of an integrable system which are everywhere small with respect to the integrable part, and not only at large values of the dynamical variables.
We also highlight that, if the above holds also for nonintegrable systems, it would imply that the addition of sufficiently strong self-interactions can stabilize any PG system.Indeed, all of our numerics support this conclusion: An explicit example is discussed in Subsec.(IV C).
At the same time, we caution that there is a -potentially crucial -distinction between integrable and nonintegrable models.Integrability ensures that whenever the motion returns to the initial position (x 0 , y 0 ), also the momenta must return to their initial values.For nonintegrable systems, this need not be the case.Nonintegrable systems may therefore be distinct in that they can potentially extract a pair of positive and negative kinetic energy from the finite region in phase space in which the above conditions are violated.This may lead to a slow growth of the explored phase-space region.While we do not observe any evidence for such slow growth in our sim-ulations, we cannot exclude it.In order to avoid such a growth, one could for example require that the Hamiltonian in the finite region be sufficiently close to integrable so that when the system returns to the initial position (x 0 , y 0 ), the momenta return to values sufficiently close to their initial values.Intriguingly, this does not a priori forbid arbitrarily strong interactions between a ghost and a positive energy degree of freedom in the finite region.Whether it is possible to make this kind of statement more quantitative and precise is an interesting question which could be investigated in the future.

A. Numerical Setup
For our numerical experiments, we draw random initial conditions x(t i ), y(t i ), x ′ (t i ), y ′ (t i ) at t i = 0 from a normal distribution with standard deviation SD.We then evolve the respective initial value problem with the Runge-Kutta method 15 up to either t f = t max or until the numerical routine detects a runaway (see conditions below), in which case t f < t max .We make sure that t max is sufficiently large such as to differentiate between cases with and without runaways.We then repeat the experiment N evol times and determine the mean evolution time f .Assuming that the criteria below confidently detect a runaway, a value of t mean = t max signals the absence of any runaway solutions in the numerical experiment.
In the above numerical experiment, we distinguish between two criteria to detect runaways: (A) We use stiffness detection of our numerical evolution scheme as a proxy for catastrophic (i.e., faster than exponential) runaways.Stiff numerical evolution typically occurs whenever very different scales enter the evolution.In practice, given that we evolve with a Runge-Kutta method, stiffness is detected by the dominant eigenvalue of the Jacobian, see [53,54] and references therein.
(B) We find that polynomial growth is harder to detect since the eigenvalues of the Jacobian can remain bounded throughout the runaway.In particular, this also includes asymptotically flat potentials in which trajectories can escape to infinity once they escape the potential in a confined phase-space region.To detect such cases, we perform the same type of experiment but we add an additional abort criterion that detects if x(t) > C or y(t) > C outgrow a specified bound C.This criterion thereby detects whenever the evolution escapes a specified We plot the ratio of the average stable evolution time tmean and the maximal chosen evolution time of tmax = 10 2 as a function of the polynomial degree N in the polynomial subclass of integrable Liouville models, cf.Eq. ( 78).For the coefficients, we choose CN = 1 and Cn<N = 0 such that all even cases are stable according to Subsec.(III C).Each point in the plot represents an average of the evolution of N evol = 100 random initial conditions drawn from a normal distribution with standard deviation SD = 1.
phase-space region.When we pick C ≫ SD, we expect that only rare cases of stable motion are falsely detected as unstable.
We will highlight below in which of the cases it makes a difference whether criterion (B) is applied in addition to criterion (A).

B. Benchmark: The Liouville Model Revisited
Here, we use the previously discussed integrable Liouville model, cf.Sec.(III), to demonstrate that the numerical experiment meaningfully distinguishes between the presence and the absence of runaways.
We first focus on the polynomial subclass, cf.Eqs. ( 78) and (85), for which the resulting potential V (x, y) reduces to a polynomial of degree N .As discussed in Subsec.(III C), the proof conditions imply stability at even degree N (assuming also that C N > 0) and, moreover, align with conditions stated in in terms of the potential.
For each polynomial degree 2 ⩽ N ⩽ 9, we perform our numerical experiment N evol = 10 2 times.We choose a maximal evolution time t max = 10 2 and draw random initial conditions from a normal distribution with standard deviation SD = 1.
The result is shown in Fig. (12) and confirms stability (instability) at even (odd) degree N .It turns out to be irrelevant whether or not we apply criterion (B).We take this as an indication that the runaways occur with faster than exponential growth.Going beyond the polynomial subclass, we also perform a numerical experiment with It turns out that, for either η ⩽ 2 (irrespective of ζ) or ζ ⩽ 2 (irrespective of η), the potential exhibits flat directions.This seems to allow for non-catastrophic runaway behaviour, the detection of which requires criterion (B).Without criterion (B), we do not detect any runaway behaviour for η ⩽ 2 or ζ ⩽ 2 which indicates that the runaways occur with polynomial growth rate.We vary η at fixed ζ = 2 and vice versa: The results are shown in Fig. (13) and reproduce the conditions required in the proof of Subsec.(III A).We note that close to the critical case, the potential becomes almost flat and the runaways are, therefore, increasingly hard to confidently detect numerically.

C. Numerical Evidence in Absence of Integrability
Going beyond integrability, we specify the following two classes of interaction potentials, i.e., both supplemented with standard quadratic terms such that (For simplicity, we only consider equal frequencies for x and y.)In both cases, the exponent δ is chosen We plot the ratio of the average stable evolution time tmean and the maximal chosen evolution time of tmax = 100 as a function of the power-law exponent δ for the three classes of potentials V I,δ (upper panel, cf.Eq. ( 118)), V II,δ (middle panel, cf.Eq. ( 119)), and V III,δ (lower panel, cf.Eq. ( 122)).Each point in the plot represents an average of the evolution of N evol = 100 random initial conditions drawn from a normal distribution with standard deviation SD = 1.
such that, as expected from the statement (ℵ), δ ≶ 0 determine the presence/absence of runaways.This is also what we find numerically (see below).
The above two example potentials involve nonpolynomial interactions (except for special choices of δ).Our third choice of example potential starts with a strictly polynomial interaction, i.e., Supplemented only with standard kinetic terms, this interaction is unstable.However, we can generalize the decoupled potentials to in which case, we find that for the full potential, i.e., the exponent δ ≶ 0, once more, seems to determine the presence/absence of runaways (at least within out finite evolution time).
We show the result of the respective numerical experiments in Fig. (14) for V = V I,δ (upper panel), for V = V II,δ (middle panel), and for V = V III,δ (lower panel).In all of these cases, the inclusion of criterion (B) does not affect the results.The latter suggests that, for δ < 0, the runaways are catastrophic, i.e., occur with exponentially fast growth rate.
Model III demonstrates in what sense added selfinteractions can apparently stabilize the motion of nonintegrable systems or, at least, dramatically increase the longevity of stable motion.In the following, we give further detail on the apparent absence of runaways and the respective characterization of the motion as longlived.
We start by noting that the numerical experiment summarized in the lower panel of Fig. (14) shows that for δ < 0, a significant fraction of initial conditions leads to a catastrophic runaway before t max = 100.This is in stark contrast to all cases with δ > 0, for which we observe no sign of instabilities at all.The specific case of δ = 1 is further exemplified in Fig. (15) where we evolve two distinct and random initial conditions up to t max = 4 × 10 4 .We generically observe two different types of evolution, depending on the specified initial conditions.
On the one hand, initial conditions which are sufficiently small, lead to motion which remains within the central region, throughout the entire evolution.For these cases, the nonintegrable character is always apparent.An example is shown in the upper panels of Fig. (15).
On the other hand, initial conditions which are sufficiently large, lead to motion which is dominated by the decoupled -and thus approximately integrable -selfinteractions.For these cases, the motion is close to integrable.In general we find that, the larger the initial conditions, the harder it is to reveal the nonintegrable character of the motion.We refrain from explicitly showing one of the cases for which we cannot resolve the nonintegrable character up to t max = 4 × 10 4 .
We also observe that, for some intermediate cases of initial conditions, the motion can transition between both types of behaviour.An example of this is shown in the lower panels of Fig. (15).At apparently random evolution times, the motion transitions back and forth between the nonintegrable character in the central region and the approximately integrable character at larger phase-space values.We show the evolution of two exemplary cases of random initial conditions (drawn from a Gaussian distribution with SD = 1 but selected with regards to the discussion in the main text) for model III with δ = 1.In the upper panels, the motion remains in the local region in which non-integrable interactions dominate: hence the motion is apparently chaotic on all timescales.In the lower panels, the motion transitions to larger phase-space variables for which the decoupled -and thus approximately integrable -potential dominates: hence the motion exhibits seemingly non-chaotic phases (cf.lower left-hand panel).At very long timescales, both initial conditions lead to chaotic motion.
We emphasize once more that, in all of these cases, the motion remains bounded in a finite phase-space region, at least up to the investigated evolution time.

V. DISCUSSION
Two point particles with opposite-sign kinetic energies can interact and yet exhibit stable motion.We prove this statement for a large class of integrable models, for which we show boundedness of motion, irrespective of the initial conditions (Lagrange stability).We gather substantial numerical evidence that nonintegrable models can be stable as well or at least longlived: within the necessarily finite numerical evolution time, we detect bounded motion and no indication of growth modes.We observe that our integrable and nonintegrable results seem to be connected by physical conditions on the interaction potential, which can be summarized as follows (as expressed in (ℵ) of Sec.(IV)): suppose a coupled Hamiltonian point-particle system with opposite-sign kinetic terms can be decoupled, in such a way that uncoupled ghosts and positive energy degrees of freedom are both separately stable due to their own self-interactions.Then, the coupled system avoids runaways, or at least exhibits very long-lived motion, if the full potential at sufficiently large coordinates is sufficiently dominated everywhere by its decoupled counterpart.
Our analytical results are contained in Sections II and III and concern integrable two-particle Hamiltonians with an additional first integral which is quadratic in the momenta.(All the corresponding known cases with positive energies are collected in App.(A).)The analytical results are organized as follows.
First, we show how one can relate ghostly and nonghostly theories by way of complex (univalent) canonical transformations and use this to build an integrable ghostly theory from a one with no ghost, i.e. to "ghostify" a positive kinetic energy theory while maintaining its integrability.We note that such ghostifications are not a mere recasting of variables as we do not complexify the initial conditions accordingly.Rather, real initial conditions on one side of the pair of theories correspond to purely complex initial conditions on the other side, and vice versa, while we work only with real variables when considering a given theory in isolation.Nevertheless, the ghostification procedure preserves integrability and thus serves as a generic construction principle for integrable models with ghostly interactions.Depending on the specific case, the Hamiltonian on either side of the pair of ghostly and non-ghostly theories can be complex-valued (and thus of little physical interest) while the other side may or may not be real-valued.With regards to the systematic construction of integrable ghostly Hamiltonians, two cases are of particular interest: First, the case in which both sides of the pair are real-valued.Second, the case in which a real-valued ghostly Hamiltonian can be obtained from a complex-valued non-ghostly Hamiltonian.We discuss examples of both these cases.
Second, we specify to a class of models with two free functions of (specific combinations) of the position variables, first obtained by Liouville.We delineate conditions on the free functions of the Liouville class and -assuming that said conditions hold -provide a rigorous proof that the resulting motion is bounded in a finite region in phase space, irrespective of the initial conditions.As in [16], the key to the proof is the additional integral of motion which, in combination with the Hamiltonian, can be used to bound the absolute value of all phasespace variables.This integral of motion is constructed from a positive definite quadratic form of the canonical momenta, similar to the usual kinetic energy, and a "potential" part which is bounded from below and is unbounded from above.The conditions enunciated on the free functions ensure a sufficiently fast growth rate of the "potential" at large phase-space variables.These properties are in direct analogy with the usual positive definite conserved energy usually employed to establish stability.Note that for this question it is irrelevant which of the two conserved quantities drives the evolution.We thus expect that other integrable models can be analysed in a similar manner.
Third, we systematically identify an exhaustive polynomial subclass of the Liouville model.This polynomial subclass is real-valued both for the ghostly and the nonghostly case.To the best of our knowledge, it has not been previously identified in the literature.As for stable motion with ghostly interactions, this polynomial subclass contains, in particular, a tower of polynomial potentials for which the proof conditions hold and for which we have thus established boundedness of motion in the presence of a ghost.The lowest non-trivial order features a sextic interaction term and quadratic terms which can be chosen at will.This provides a strikingly simple model of ghostly interactions which is polynomial, integrable, and globally stable.Moreover, this system exhibits spontaneous symmetry breaking entirely due to the interactions with the ghost.Remarkably, the very same interactions which are often thought to catastrophically destabilize the system, instead introduce novel (Lyapunov) stable vacua.
We highlight that all the above three parts of our analytical work are systematically presented and can, presumably, be extended to other integrable systems.Our analysis is therefore far from exhaustive: Even among the known set of all integrable two-particle systems with a constant of motion, which is quadratic in the momenta (see App. (A)), it is thus possible that further stable ghostly systems can be obtained in the same way.
Our numerical results are contained in Sec.(IV) and concern the extension to nonintegrable models.The analytical proof for boundedness of motion in the integrable Liouville class is quite striking.It is, however, equally striking that we find no indication that integrability is necessary.All of the investigated nonintegrable models seem to behave like the integrable ones.If the interactions at large phase-space variables do not obey the observed conditions (ℵ), we quickly detect runaway solutions which either lead to a breakdown of the evolution in finite time or to a less severe runaway.In contrast, whenever the interactions at large phase-space variables obey the observed conditions (ℵ), we find no indication for the onset of runaways.
The present work calls for various extensions in obvious directions.First, it is clear that our results should allow for an extension to the multiparticle case, i.e., to one or several ghost DoFs interacting with several positive-kinetic-energy DoFs.Second, an obvious question is that of quantization.A similar issue has been studied for some specific analogous PG system where stability could only be studied numerically and/or where the ghost still caused an instability, albeit an arguably "benign", i.e., a slowly-growing one (see e.g.[29]).In this context the simple sextic polynomial potential which we have established to be integrable and stable presents a starting point to investigate whether stability persists, as we actually expect, at the quantum level.Last, we are planning to extend our work to field theories.In particular, one can construct a theory of two interacting scalar fields ϕ and ψ, one with positive and one with negative kinetic energy, for which the interaction potential is equivalent to the one investigated here.Once more, the simple sextic polynomial potential found here, appears to be especially amenable for such a study.Finally, an extension of our analysis to gravitational interactions remains an important question for future work.
In comparison to [49], we have corrected a typo in the integral of motion I of class 3 (given here in Eq. (A10)) which differs from [49, Eq. (3.2.13)] but agrees with [55, Eq. ( 41)] (given there in terms of A = c/2, F 1 = 2f , and F 2 = 2g).Invariants which include linear and quadratic terms in the momenta are prohibited by time reflection symmetry, cf.[49].
The above functional classes are known to contain several polynomial subclasses, cf.[49] for a collection of those.In Appendices D and E, we identify the polynomial subclasses for class 1 (in the main text also referred to as Liouville class) and class 3. To the best of our knowledge these polynomial subclasses have not been previously identified.This equation can only be satisfied for arbitrary p x and p y provided the following conditions hold As we assume that K ik (x, y) is a smooth function of x and y, the first two equations of this system in Eq. (B4) imply that extrema of V (x, y) are also extrema of U (x, y).In particular, this means that equilibrium points of the system in Eq. (B1) are always located at the extrema of U .The last equality of the system in Eq. (B4) can be written as Again, this equation should hold for arbitrary p x and p y , which results in the following conditions ∂ y K xx = −2σ∂ x K xy , and ∂ x K yy = −2σ∂ y K xy , (B7) where we used that σ 2 = 1 to go from Eq. (B5) to Eq. (B7).Differentiating the first of the conditions in Eq. (B7) with x and second with y and using Eq.(B6) we obtain The integral of motion I is useful to establish stability when the kinetic part of I is positive definite.This can be achieved when the terms linear in coordinates are absent i.e. when κ x = κ y = 0. Let us compare the expression in Eq. (B12) with the class 1 integral of motion J LV from Eq. ( 41) with σ = −1 where we have κ 0 = κ x = κ y = 0 and the non-vanishing coefficients are , and κ xy = 1 .
After the form of K ik is fixed, one can solve the system of first-order PDEs in Eq. (B4) to find V (x, y) and U (x, y).This task is beyond the scope of the paper.The corresponding solutions of this system are reviewed in e.g [49] see also App. (A).
The main point of this appendix was to demonstrate that equilibrium configurations (∂ x V = ∂ y V = 0) of integrable system in Eq. (B1) are all located at extrema of U (x, y).However, the opposite is not true in general.Indeed, as the structure of K ik in Eq. (B12) clearly shows, the matrix of the system of the first two equations in Eq. (B4) may not be invertible.In particular, one can choose κ xy = κ y = κ yy so that K yy = 0 and one cannot invert Eq. (B4) to uniquely express ∂ i V through ∂ k U .This happens, e.g., for the models in Eqs.(A11) and (A15).
Finally let us discuss a useful relation between Hessian determinant of the potential U , and that of the potential V at an equilibrium point.The vanishing of the derivatives of V at such points implies, using Eq.(B4), that the derivatives of K ik will make no contribution to second derivatives of U there.Using then the first equation in (B4) to compute ∂ y ∂ x U and the second equation in (B4) to compute ∂ x ∂ y U , which facilitates the computation at equilibrium points, one finds after some algebra, that one has there If the the first integral is useful to establish stability, then K ik is positive definite and not degenerate, so that det K > 0. In that case, for PP case σ = 1 so that minimum of U corresponds to the minimum of V .However, for PG case σ = −1 and Hess U > 0 corresponds to Hess V < 0. Hence the minimum of U can only be achieved at saddle points of V .Thus, stable equilibrium points in systems where a ghost y interacts with a usual DoF x are located in saddle points of the potential V .
where as above, H LV,GP is a function of c but a functional of Φ : u 2 → Φ(u 2 ) and Ψ : u 2 → Ψ(u 2 ).So we have that showing that the PG and GP theories just differ by a global sign in the Hamiltonian, a proper redefinition of the function Φ and Ψ (or f and g), and a change of sign of c; hence these theories are equivalent.In particular the conditions enunciated in (i) and (ii) on the functions f an g of the proof of Subsec.(III A) yield a stable theory PG Eq. (C12) for c < 0 and a stable theory GP Eq. (C12) for c > 0.
) This rewriting captures the first two non-trivial cases, i.e., N = 3 and N = 4, for which we visualize examples of the potential in Fig. (3) and of the resulting phasespace motion in Fig. (4) (for N = 3) and Figs. 5 to 7 (for N = 4).

(class 3 )FIG. 5 .FIG. 6 .
FIG.5.We show the trajectory (for random initial conditions drawn from a Gaussian distributions with standard deviation SD = 2) for N = 4, i.e., for ω 2x = 1, ω 2 y = 1, C4 = 1, and c = 1 in Eq. (79).In the left-hand panel, we superimpose a contour plot of the potential (cf.right-hand panel in Fig.(3)) where increasingly light shades indicate increasingly large (positive) values of the potential.In the right-hand panel, we show the respective evolution of the momenta.

FIG. 10 .
FIG.10.For values of C4 satisfying Eq. (114), we show a logarithmic plot of the "uplifted" effective "potential" U (4) LV (x, y), i.e., the plot of log 1 + U (4) LV (x, y) − U (4) LV (xI , yI ) , such that we take the logarithm of a strictly positive quantity.The red points indicate the four global vacua (xI , yI ) given by Eqs.(112) and (113) with all combinations of signs for the square roots.The black point indicates the vacuum at the origin which, however, is only a local minimum of U (4)

FIG. 11 .
FIG. 11.On the left, we show the effective "potential" U

( 4 )
LV (x, y) − U (4) LV x b(+) , 0 , on larger scales.The red points indicate the two global vacua ±x b(+) , 0 from Eq. (108).The black point indicates the trivial vacuum at the origin, which however is only a local minimum of U (4) LV .The parameters are set to c = 1, ωx = 1, and ωy = 1.7, as in Fig. (10), while C4 is higher than that of Fig.(10) and is in the middle of the interval (111): C4 ≈ 0.221.This case has different numerical values of constants, but its structure of vacua corresponds to case II in Fig.(9).
FIG. 12.We plot the ratio of the average stable evolution time tmean and the maximal chosen evolution time of tmax = 10 2 as a function of the polynomial degree N in the polynomial subclass of integrable Liouville models, cf.Eq. (78).For the coefficients, we choose CN = 1 and Cn<N = 0 such that all even cases are stable according to Subsec.(III C).Each point in the plot represents an average of the evolution of N evol = 100 random initial conditions drawn from a normal distribution with standard deviation SD = 1.

FIG. 13 .
FIG.13.We plot the ratio of the average stable evolution time tmean and the maximal chosen evolution time of tmax = 10 3 as a function of (i) the exponent η (at fixed ζ = 2) in the left-hand panel, and for (ii) the exponent ζ (at fixed η = 2) in the right-hand panel, cf.Eq. (115).Each point in the plot represents an average of the evolution of N evol = 100 random initial conditions drawn from a normal distribution with standard deviation SD = 1.The numerical evolution is terminated if stiffness is detected or if x(t) > C or y(t) > C, with C = 10 3 .
FIG. 14.We plot the ratio of the average stable evolution time tmean and the maximal chosen evolution time of tmax = 100 as a function of the power-law exponent δ for the three classes of potentials V I,δ (upper panel, cf.Eq. (118)), V II,δ (middle panel, cf.Eq. (119)), and V III,δ (lower panel, cf.Eq. (122)).Each point in the plot represents an average of the evolution of N evol = 100 random initial conditions drawn from a normal distribution with standard deviation SD = 1.
FIG.15.We show the evolution of two exemplary cases of random initial conditions (drawn from a Gaussian distribution with SD = 1 but selected with regards to the discussion in the main text) for model III with δ = 1.In the upper panels, the motion remains in the local region in which non-integrable interactions dominate: hence the motion is apparently chaotic on all timescales.In the lower panels, the motion transitions to larger phase-space variables for which the decoupled -and thus approximately integrable -potential dominates: hence the motion exhibits seemingly non-chaotic phases (cf.lower left-hand panel).At very long timescales, both initial conditions lead to chaotic motion.
Appendix B: Integrability and VacuaFor self-completeness of the paper, as well as to obtain some intermediate results useful for the discussion of Subsec.(III E), let us briefly review when a system with Hamiltonian integral that is quadratic in momentaI = K ik (x, y) p i p k + U (x, y) ,(B2)where i, k run over {x, y} and σ = ±1.The vanishing of the Poisson bracket {H, I} PB implies thatp x [∂ x U − 2 (K xx ∂ x V + K xy ∂ y V )] + (B3) + p y [σ∂ y U − 2 (K yx ∂ x V + K yy ∂ y V )] + + p i p k p x ∂ x K ik + σp y ∂ y K ik = 0 .