Integrability and trajectory confinement in $\mathcal{PT}$-symmetric waveguide arrays

We consider $\mathcal{PT}$-symmetric ring-like arrays of optical waveguides with purely nonlinear gain and loss. Regardless of the value of the gain-loss coefficient, these systems are protected from spontaneous $\mathcal{PT}$-symmetry breaking. If the nonhermitian part of the array matrix has cross-compensating structure, the total power in such a system remains bounded -- or even constant -- at all times. We identify two-, three-, and four-waveguide arrays with cross-compensatory nonlinear gain and loss that constitute completely integrable Hamiltonian systems.


I. INTRODUCTION
The concept of PT symmetry, originally introduced in nonhermitian quantum mechanics [1,2], has led to significant developments in photonics, plasmonics, quantum optics of atomic gases, metamaterials, Bose-Einstein condensates, electronic circuitry, and acoustics [2,3]. The PT -symmetric equations model physical structures with a built-in balance between gain and loss.
In a linear nonhermitian system, raising the gain-loss coefficient above a critical level causes trajectories to escape to infinity. The setting in of this blow-up instability can have a detrimental effect on the structurein optics, for example, an escaping trajectory implies an uncontrollable power hike.
Taking into account nonlinear effects and adding nonlinear corrections to the equations may arrest the blow-up through the emergence of conserved quantities confining trajectories to a finite part of the phase space [4]. In this paper, we study a class of PT -symmetric systems with the nonhermiticity induced entirely by the nonlinear terms. On one hand, these systems provide access to the full set of behaviours afforded by the presence of gain and loss. On the other hand, they exhibit remarkable regularity properties such as the existence of the Hamiltonian structure and trajectory confinement.
The Hamiltonian structure in a dynamical system imposes a deep symmetry between two sets of coordinates parameterising its phase space. In the presence of additional first integrals, the Hamiltonian structure establishes an even higher degree of regularity: the Liouville integrability. We identify two-, three-and fourcomponent integrable systems with nonhermitian nonlinearities that are free from the blow-up behaviour.
The two-component complex systems that we study have the form of the nonlinear Schrödinger dimer: a discrete Schrödinger equation, defined on only two sites [5][6][7][8][9][10]. The nonhermitian dimer serves as an archetypal model for a pair of optical waveguides or a pair of micro-ring resonators with gain and loss, coupled by their evanescent fields [11,12]. It also arises in the study of Bose-Einstein condensates [13,14], plasmonics [15], spintronics [16], electronic circuitry [17] and several other contexts.
The literature suggests several recipes for preventing the blow-up in dimers, including linear vs nonlinear gainloss competition [18,19] and nonlinear gain-loss saturation [12,20]. This paper explores cross-stimulationan alternative mechanism that, in addition to ensuring nonsingular evolution, conserves the norm |u| 2 + |v| 2 (interpreted as the total power of light in the optical context). We present two cross-stimulated PT -symmetric dimers that describe completely integrable Hamiltonian systems.
The cross-stimulation is a special type of a more general notion of cross-compensation of gain and loss in a multichannel structure. To illustrate this concept, we invoke a three-and four-site discrete Schrödinger equation -the PT -symmetric trimer and quadrimer, respectively [21]. In the optical domain, the trimer and quadrimer model an optical necklace -an array of three or four coupled waveguides or resonators. We produce examples of a completely integrable PT -symmetric trimer and quadrimer, with solutions free from the blow-up behaviour. Similar to the dimers in the first part of this study, the nonherimiticity of these necklaces is entirely due to the nonlinear terms and has a cross-compensatory character.
The paper is organised into five sections. We start with a linear hermitian dimer with the cross-stimulating cubic gain and loss (section II). The subsequent section (section III) deals with a slightly more complex crossstimulating system whose hermitian part is cubic itself. We uncover the hidden Hamiltonian structure of these systems, determine their integrals of motion and construct analytic solutions. An integrable PT -symmetric trimer and quadrimer are identified in section IV. Section V summarises results of this study.

A. Gain and loss cross-stimulation
As the gain-loss coefficient is varied, the topological structure of the phase portrait of a dimer with linear gain and loss undergoes a spontaneous change. Consider, for example, the so-called standard PT -symmetric dimer, a model that arises in a wide range of physical contexts [6,[8][9][10]: As γ is raised above the critical value γ c = 1, the fixed point at u = v = 0 loses its stability and the PTsymmetry is said to become spontaneously broken. In the symmetry-broken phase (γ > 1) small initial conditions give rise to exponentially growing solutions. Since the change of behaviour concerns solutions of the linearised equations, we refer to this bifurcation as the linearised symmetry breaking. In a generic Schrödinger dimer, trajectories resulting from the exponentially growing solutions of the linearised equations may escape to infinity. For example, in the standard dimer (1) all growing linearised solutions give rise to escaping trajectories [9]. (Note that unbounded solutions may occur in the symmetric phase too -they just need to evolve out of large initial data [8,9].) To preclude the blow-up of small initial data, we consider a system with a hermitian linearised matrix. Furthermore, our dimer is assumed to be cross-stimulated. This means that the u-channel gains energy at the rate proportional to the power carried by its v-neighbour while the v amplitude loses energy at a rate proportional to the power carried by u: As a result of the cross-stimulation, the total power P = |u| 2 + |v| 2 is conserved. This keeps all trajectoriesboth with small and large initial data -in a finite part of the phase space.
Although the cross-stimulation may come across as a purely mathematical construct, the system (2) originates in a well-established physical context. It describes the spin-torque oscillator -an isotropic ferromagnet in an external magnetic field with polarised spin current driven through it [16]. The components of the total spin vector S = {X, Y, Z} in the free layer of the oscillator are expressible through the complex amplitudes of the dimer's channels: The length of the vector, R = √ X 2 + Y 2 + Z 2 , is proportional to the total power carried by the dimer, which is fixed by the initial condition. The direction of the vector is determined by a system of three equations [16]:Ẋ where X, Y and Z are considered to be functions of τ = 2t, and the overdot stands for the derivative w.r.t. τ . Another area of applications of the cross-stimulated dimers is the classical limit of quantum theory. A system equivalent to (2) governs the normalised state vector of the nonhermitian two-level atom [13]. The corresponding Bloch-sphere dynamics obey equations (4).

B. Hamiltonian structure and integrability
Upon defining a set of polar coordinates by equations (4) simplify tȯ It is not difficult to see that the system (6)- (7) is Hamiltonian. We choose the integral R 2 − γ −2 as the Hamilton function: Here η and θ are the canonical coordinates and p η = 2γ −1 Z is the canonical momentum conjugate to η. The equationsη reproduce equations (6), anḋ amounts to (7). The last Hamilton equation, serves as a definition of the momentum p θ : The existence of the Hamiltonian structure and an additional first integral -the azimuthal angle θ -establishes the complete integrability of the dimer (2).
FIG. 1. Trajectories of the system (4) on the surface of the sphere of the radius R > γ −1 (a); R = γ −1 (b) and R < γ −1 (c). (In this plot, γ = 1 2 .) Note that the largest sphere (panel a) has been scaled down while the smallest one (panel c) has been scaled up. In each ball, tinted is the equatorial plane.

C. Solutions
The trajectories of the system (4) are circular arcs resulting from the section of the sphere X 2 + Y 2 + Z 2 = R 2 by vertical half-planes θ = const. These were classified in an earlier study [13] in terms of the fixed points of the system. Changing to the canonical variables (5) allows us to express solutions of the dimer (2) explicitly.
Using (8) and the first equation in (6), we obtaiṅ where K = e η and ν 2 = γ 2 R 2 − 1. This is a conserved quantity of the linear equation The sphere with radius R > γ −1 encloses a section of the vertical axis X = 0, Y = −γ −1 . Accordingly, the sphere is crossed by vertical half-planes with all possible azimuthal angles, from θ = 0 to θ = 2π. The resulting trajectories are circular arcs that connect an attracting fixed point (a stable node) at to a repelling point (unstable node) that has the same X, Y and opposite Z (Fig 1(a)). The corresponding solutions of (10) are given by where the relation between the constants of integration is found by substituting (11) into (9): For each value of spherical radius R the equations (4) represent a dynamical system on the (X, Y )-plane, with a conserved quantity θ(X, Y ). We note an ostensible paradox, where the existence of the conservation law should be precluded by the presence of the attractor and repeller in the phase space. The paradox is resolved, however, upon observing that θ is undefined at either fixed point.
As R is decreased through γ −1 , the attracting and repelling points approach each other, collide and then diverge along the equator of the sphere. The sphere with R < γ −1 has two elliptic fixed points, at Each point is surrounded by a family of circular orbits (Fig 1(c)). The corresponding solutions of (10) are periodic: Here, The vertical axis X = 0, Y = −γ −1 lies outside the sphere with R < γ −1 . As a result, the sphere is only crossed by the half-planes with θ in the interval defined by equation (12): arccos(γR) < θ < π − arccos(γR).
In the borderline case R = γ −1 , the system (4) has a single fixed point, at X = Z = 0 and Y = −γ −1 . Trajectories are homoclinic: they emerge out of the fixed point, wrap around the sphere and flow into the same semistable point (Fig 1(b)). The corresponding solution of equation (9) is Once we have an explicit expression for X(τ ), Y (τ ) and Z(τ ), the corresponding dimer components can be easily reconstructed:

D. Related systems
We close this section with three remarks. Firstly, the absence of any conservative nonlinearity in the crossstimulated dimer is not a prerequisite for its integrability. A simple example of the model with a nonlinear hermitian part and cross-stimulated channels is The dimer (13) maps onto our system (2) and inherits its integrability and trajectory-confinement property. The gauge transformation relating the two systems is simply Secondly, the transformation takes a solution of the system (2) with |u| 2 + |v| 2 = 2R to a solution of the linear dimer This correspondence accounts for the reducibility of the cross-stimulated dimer to a linear equation (equation (10)).
Finally, it is appropriate to mention an earlier study [19] in which the authors constructed integrals of motion of the system The model (15) includes our dimer (2) as a particular case with α = β = 0. However, the integrals of motion of (15) do not persist as β → 0; see Ref [19].

A. The system
Having outlined the effect of cross-stimulation on a model dimer (2), we observe that this mechanism remains available to a broad class of systems in nonlinear optics. If we retain the generic Kerr nonlinearity in equations (2), we obtain another cross-stimulated dimer conserving the total power |u| 2 + |v| 2 : Equations (16) result from the model of a birefringent single-mode fibre amplifier with a saturable nonlinearity [20]: Here E 1 and E 2 are the amplitudes of the orthogonally polarised modes; R > 0 is the nonlinearity parameter andγ > 0 is the gain-loss coefficient. The components E 1 and E 2 can be obtained from solutions of (16) with |u| 2 + |v| 2 = 2R and γ =γ/R by means of the transformation (14).
The cross-stimulated dimer (16) can be written in terms of the spin variables (3): Here X, Y, Z are functions of τ = 2t, and overdots denote derivatives w.r.t τ . The spin vector S = {X, Y, Z} will provide an appropriate basis for visualising the solutions of the Kerr dimer. We further observe that the system (17) arises in the mean-field approximation of the many-body Bose-Hubbard model [13].

B. Hamiltonian structure and integrability
If we define the polar coordinates η and θ such that X = 1 1 + γ 2 +e −η cos θ, Y = − γ 1 + γ 2 +e −η sin θ, (18) equations (17a) and (17b) reduce tȯ Accordingly, is a conserved quantity, in addition to (21) To uncover the Hamiltonian formulation, we appoint η and ξ as two canonical coordinates and the first integral as the Hamiltonian: Here p η = 2γ −1 Z. The first equation in (19) and equation (17c) acquire the forṁ hence p η is the momentum canonically conjugate to η. The momentum p ξ conjugate to the coordinate ξ can be found from the Hamilton equatioṅ by integration: The existence of the canonical formulation and an additional first integral (ξ) establishes the Liouville integrability of the dimer (16).

C. Fictitious particle formalism
Transforming to the canonical variables allows us to obtain the general analytical solution of the system (17).

The symmetry identities (34), (35) and (36) yield the inequality
D. Spin trajectories from particle flight paths Consider a trajectory of the spin system (17) passing through a point on the equator (Z = 0) at time τ = 0. The corresponding fictitious particle starts its motion from rest (θ(0) = 0), with the values of integrals ξ and R defined by the initial data, θ 0 and η(0): The implicit solution τ (θ) of equation (24) is given by the integral .
The solution is valid for all θ such that the expression under the radical is positive in the interval (θ 0 , θ). Without loss of generality we may let θ 0 lie in the interval [0, 2π). If ∂U ξ /∂θ < 0 at θ = θ 0 , the θ-particle will start moving in the positive direction and we choose the positive sign in (40). The corresponding spin trajectory S(τ ) will emerge into the northern hemisphere. If ∂U ξ /∂θ > 0 at θ = θ 0 , the particle will start moving in the negative direction and we choose the negative sign in (40). In that case, the point S will move into the southern hemisphere.
The negative-time motions can be classified in a similar manner.
For a given value of γ, the character of motion is determined by the values of the parameters R and ξ. Since two or three values of θ 0 ∈ [0, 2π) can be mapped to the same U ξ , we also need to indicate the position of θ 0 relative to the minimum and maximum of U ξ (θ).
spin sphere have been numerically delineated in Ref [13].) Similar to the system (4) (Fig 1 (a)) the sphere with R > γ −1 supports a pair of latitudinal fixed points, with The fictitious particle's journeys from infinity to θ 0 and back to infinity correspond to trajectories S(τ ) emerging from the unstable focus in the southern hemisphere and spiralling into the attractor in the north (Fig 3(a)).
When ξ < ξ c or ξ > ξ d (R), the total energy R 2 of the particle is greater than U ξ at all of its local maxima. Accordingly, the particle escapes to infinity as τ → ±∞. This class of motions corresponds to the heteroclinic trajectories on the sphere flowing from the southern to the northern focus in (41) (Fig 3 (b)).
Turning to the interval ξ c < ξ < ξ d we first assume that the point θ 0 lies to the left of the local maximum θ (2) . In that case the energy of the θ-particle is insufficient to overcome the potential barrier: R 2 < U ξ (θ (2) (ξ)). The particle becomes trapped in the potential well, and its motion is periodic with period Here θ turn is the second root of the equation R 2 = U ξ (θ) to the left of the maximum of the potential in the interval [0, 2π).
To interpret the periodic motions in terms of trajectories on the sphere, we note that as R is decreased through γ −1 , a saddle-centre bifurcation brings about two new fixed points lying on the equator: For R in the interval (1 + γ 2 ) −1/2 < R < γ −1 , the point with negative X is a centre and the one with positive X is a saddle. The oscillations of the θ-particle in the potential well translate into a thicket of closed orbits surrounding the centre (Fig 3(b)).
By examining the neighbourhood of the saddle point in Fig 3(b)), one can readily reconstruct a homoclinic curve connecting the saddle to itself. This trajectory results by choosing ξ = ξ d (R). The homoclinic curve separates the family of closed orbits from the focus-to-focus flows.
The disappearance of aperiodic solutions θ(τ ) is consistent with the bifurcation occurring as R is reduced below (1 + γ 2 ) −1/2 . At this value of R, the pair of latitudinal fixed points (41) merges with the saddle on the equator forming the second centre (Fig 4(a)). The two centre points are given by equations (42); each point is surrounded by a family of closed orbits (Fig 4 (b)). The closed orbits are described by periodic motions of the particle in the potential well.

IV. INTEGRABLE PT -SYMMETRIC NECKLACES
In a nonhermitian necklace of more than two waveguides, the trapping of trajectories in a finite part of the phase space may require each site to coordinate its gain and loss rate with both of its left and right neighbours. This cross-compensation mechanism is more subtle than the cross-stimulation of two channels of a dimer. In this section, we exemplify the cross-compensation with arrays consisting of three and four elements.

A. Trimer
We identified several PT -symmetric trimers endowed with a Hamiltonian structure and possessing an integral of motion that prevents the blow-up behaviour. However, only one of these models has three first integrals in involution and defines a Liouville-integrable dynamical system.
The trimer in question is a nonhermitian extension of the closed Ablowitz-Ladik chain: The system is invariant under the product of the P and T transformations, where Similar to the dimers (2) and (16), the linearisation of equations (43) about u = v = w = 0 gives a system with a hermitian matrix and a purely real spectrum. Accordingly, the trimer (43) does not suffer the linearised PT -symmetry breaking as the value of γ is raised.
Note that since the total power P = |u| 2 + |v| 2 + |w| 2 is bounded from above by the conserved quantity Π, the trimer does not exhibit any unbounded trajectories.
It is worth noting that the absence of blow-up regimes is not related to the integrability of the model but is instead a consequence of the compensation of gain and loss rate in the neighbouring channels. This can be illustrated by the following family of nonintegrable crosscompensated trimers: iv In equations (45), f n (ρ) (n = 1, 2, 3) are positive functions with f n (0) = 1 and f n (ρ) ≤ 1 for ρ > 0. The model (45) is Hamiltonian with the Hamilton function (44) and obvious modification of the bracket. It has the first integral .
Since F n (ρ) ≥ ρ, the total power is bounded from above by I and no trajectory can escape to infinity.

B. Quadrimer
Guided by the symplectic structure and crosscompensating arrangement of the integrable PTsymmetric trimer, it is not difficult to construct a fourwaveguide necklace with similar properties: The quadrimer (46) is PT -symmetric, with the P operator defined by and T as in T u n (t) = u * n (−t). The linearisation of the quadrimer (46) about u 1 = u 2 = u 3 = u 4 = 0 gives a system with a hermitian matrix and a purely real spectrum.
The quadrimer retains a similar canonical representation to the trimer (43), with the Hamiltonian (47) and GIK bracket. Here In direct analogy to the trimer, the system conserves the total momentum as well as the product (1 + |u n | 2 ).
In addition, the quadrimer conserves a quartic quantity All four integrals H, M , Π and I mutually commute, and consequently the system (46) is completely integrable.

V. CONCLUDING REMARKS
In this paper, we explored a class of the PT -symmetric discrete Schrödinger equations with purely-nonlinear nonhermitian terms. An a priori advantage of systems of this type in physics is that the trivial solution remains stable and linearised PT -symmetry unbroken regardless of the value of the gain-loss coefficient.
We have identified two nonequivalent dimers that, in addition to preserving the linearised PT -symmetry, conserve the quantity |u| 2 + |v| 2 . (In the optical context, this means that the total power of light is conserved.) The power conservation is brought about by the crossstimulation of the two channels of the dimer. The constancy of |u| 2 + |v| 2 ensures that no trajectories of these dynamical systems escape to infinity.
Remarkably, both dimers possess a canonical structure and a first integral independent of the Hamiltonian. This establishes the Liouville integrability of the two systems. The transformation to the canonical variables allowed us to construct their general analytical solution.
The cross-compensation of the gain and loss rate in the neighbouring channels remains an efficient blow-up prevention mechanism in the case of the nonhermitian Schrödinger necklaces -the ring-shaped arrays of N waveguides [23]. We have exemplified this idea with the construction of a family of trimers (N = 3) whose trajectories are confined to a finite part of their phase space. The entire family is endowed with a canonical structure while one member of the family has 3 first integrals in involution and defines a completely integrable system. Finally, we have identified a completely integrable cross-compensated quadrimer -a PT -symmetric necklace of N = 4 waveguides.

ACKNOWLEDGMENTS
A discussion with Robert McKay is gratefully acknowledged. We thank Andrey Miroshnichenko and Tsampikos Kottos for instructive correspondence. This research was supported by the NRF of South Africa (grant No 120844).