Perpetual emulation threshold of PT-symmetric Hamiltonians

We describe a technique to emulate the dynamics of two-level -symmetric spin Hamiltonians, replete with gain and loss, using the unitary dynamics of a larger quantum system. The two-level system in question is embedded in a subspace of a four-level Hamiltonian, with the exterior levels acting as reservoirs. The emulation time is normally finite, limited by the depletion of the reservoirs. We show that it is possible to emulate the desired behaviour of the -symmetric Hamiltonian without depleting the reservoir levels, by including an additional coupling between them. This extends the emulation time indefinitely, when in the unbroken symmetry phase of the non-unitary dynamics. We propose a realistic experimental implementation using dynamically decoupled magnetic sublevels of ultracold atoms.


Introduction
Quantum theory is our most successful description of nature. It describes the dynamics of closed systems by means of Hermitian Hamiltonians. The Hermitian constraint leads to unitary evolution and a real, and hence measurable, eigenvalue spectrum. Remarkably, this constraint is stronger than necessary, and a subclass of non-Hermitian Hamiltonians also affords real-valued eigenvalues as long as the Hamiltonian operator commutes with the joint parity-time PT operator. These models emerged from perturbative approaches in quantum field theory [1,2] and can have all-real eigenvalues in some parameter regimes [3,4]. In such models, a free parameter of the Hamiltonian drives a phase transition between two regimes where the PT symmetry is either unbroken or broken, leading to real and complex eigenvalues respectively [5].
The language of PT -symmetric Hamiltonians can offer a simplified, alternative formulation of the Lindblad equation description of open systems [6], and was used to study numerous interesting and disparate phenomena ranging from photosynthesis [7] to the slowing down of decoherence in the vicinity of the critical point [8]. Generally, it is possible to derive a non-Hermitian Hamiltonian from a given Lindblad equation, that is a valid description of the systems' dynamics until the probability of a quantum jump becomes too large.
A number of experiments in the classical domain have demonstrated PT -symmetric systems, e.g., using coupled optical waveguides [9,10,11], and Floquet systems [12]. The first demonstration in the quantum domain used Sinai billiards [13,14]. The physics of these systems becomes even richer close to the critical point between the broken and unbroken regime; the phase associated with encircling the critical point is similar to the Berry phase associated with a Dirac point in non-trivial topological materials [15,16].
In a series of papers [17,18,19], Wunner and coworkers showed that it is possible in principle to reproduce the dynamics of a two-level PT -symmetric Hamiltonian with a four-level Hermitian system. In this configuration two of the levels map to the PTsymmetric levels and the remaining two act as probability amplitude source and sink. The embedded PT -subspace is coupled to the rest of the system using time-dependent transition matrix elements that emulate the non-unitary particle flow between the subspaces; we refer to this scheme as coherent emulation. Here we extend these works, initially by adopting a different calculational approach that permits analytic closed-form solutions in these four-level systems. This approach makes it explicit that the emulation of the two-level PT -symmetric dynamics in these works is always limited in time, and that it breaks down due to a depletion of the probability amplitude in the source level.
We then consider an extended scheme that includes an additional coupling between source and sink levels, to counter the depletion of the source level. We show that this amplitude recycling increases the emulation time for any given set of parameters. Moreover, when the coupling strength of the amplitude-recycling field exceeds a critical value, the duration of the emulation is extended indefinitely. Above this threshold the system dynamics dramatically change such that a genuinely periodic behaviour occurs. Implementation of this scheme requires the solution of a pair of coupled, nonlinear differential equations to pre-compute the time-dependent Hamiltonian required to emulate PT -symmetric behaviour in the target subspace. Our protocol therefore consists of a classical computation device to design a control sequence, different for all initial conditions, and a quantum system to implement the sequence. We use these results as the basis for the discussion of a realisation of the PT -symmetric dynamics in pseudospins corresponding to internal states of the groundstate manifold of alkali atoms cooled to degeneracy.

Two-level PT -symmetric Hamiltonian
We focus on an extensively studied minimal model [20,21,22,23,10,13], that captures all the relevant information of such systems, and aim to describe the full dynamics governed by h for all times and Γ, Λ > 0. The parity operator for this model is P = σ x , where σ i , i = x, y, z are the Pauli operators, and time reversal gives T hT −1 = h * . For this Hamiltonian, PT -symmetry implies Ph − h * P = 0, which is manifestly obeyed by Eq. (1). The eigenvalues of h are This sets a natural energy scale for the system and allows us to use dimensionless time τ = αt; coupling λ = 2Λ/α; and gain γ = 2Γ/α. We assume here that we stay in the PT -unbroken phase such that α is real, i.e. Γ < Λ. As with Hermitian Hamiltonians we define a time evolution operator U(τ ) that propagates the state ψ(τ ) = (ψ 1 (τ ), ψ 2 (τ )) T to time τ : where I is the 2 × 2 unit matrix. We consider initial conditions with real-valued amplitudes throughout ψ(0) = (cos(θ/2), sin(θ/2)). Several aspects of this solution are of particular interest here. First, we consider the norm n(τ ) ≡ |ψ 1 (τ )| 2 +|ψ 2 (τ )| 2 and the atomic magnetisation w(τ ) ≡ |ψ 1 (τ )| 2 −|ψ 2 (τ )| 2 . Although our technique is valid for all initial states, for brevity of exposition we base our subsequent discussion on the symmetric case when θ = π/2. In this case, we obtain n(τ ) = λ 2 − γ 2 cos(τ ); (3) For γ → 0, the norm becomes λ 2 → 1, and the magnetisation disappears. In this limit, there is no evolution since the initial state is an eigenstate of h. For θ = π/2, oscillations in the magnetisation only occur when h is non-Hermitian, and the oscillation frequency is given by the splitting α. Also of interest is the relative phase between wave-function components which, again for θ = π/2, is defined by tan (Φ(τ )) = γλ(1 − cos(τ )).
Perpetual emulation threshold of PT -symmetric Hamiltonians  This phase difference is plotted in Fig. 1.
Two important timescales emerge from these equations. The first is the period of the PT -symmetric oscillations in Eq. (4), which, with our choice of scaling for τ , is simply 2π. The second corresponds to important changes that take place at shorter times. These are apparent in the phase difference in Eq. (5), plotted in Fig. 1. For small γ/λ, we may approximate Φ(τ ) ≈ 2γ/λ sin 2 (τ /2) and the dynamics unfold with the same period as the magnetisation of Eq. (4). However, for increasing γ, the nonlinearity of the arctan function becomes important, and gives rise to sharp changes in the relative phase.

Emulation with a four-level system
Given that our understanding of Nature is in terms of Hermitian dynamics, we seek to create designer subspaces of larger systems that evolve according to a desired PTsymmetric Hamiltonian. We simulate the non-unitary dynamics of Eq. (1) using a four-level system with wavefunction φ(τ ) = (φ 0 (τ ), φ 1 (τ ), φ 2 (τ ), φ 3 (τ )) T . We encode the PT dynamics in the two central levels φ 1 , and φ 2 , and use φ 0 and φ 3 as the sink and source levels. Our 'emulator' will thus be a system with dimensionless Hamiltonian Perpetual emulation threshold of PT -symmetric Hamiltonians These couplings are sufficient to completely emulate the time-dynamics of PTsymmetric Hamiltonians as in Eq. (1). The dimensionless detunings δ 0 (τ ), δ 3 (τ ) and the couplings Ω 01 (τ ), Ω 23 (τ ) are time-dependent functions that we choose such that the behaviour of the two levels in the central subspace matches that of the PT -symmetric system, i.e. such that φ 1 (τ ) = ψ 1 (τ ) and φ 2 (τ ) = ψ 2 (τ ). The static coupling Ω 03 between the source level φ 3 and sink level φ 0 replenishes the source population by transferring particles from the sink level. We initialise the system such that the population of the central subspace always starts as ψ 1 (0) 2 + ψ 2 (0) 2 = 1; this implies a rescaling of φ 1 , φ 2 by an appropriate factor so that it is normalised to unity at t = 0. We determine the dynamic detunings and couplings using the following procedure. We first split the wavefunction components into real and imaginary parts, φ i = φ R i + iφ I i , and use the Schrödinger equation for our four-level system to get the equations: where we have replaced φ 1 and φ 2 with their target wavefunctions ψ 1 and ψ 2 , and wherė φ i = dφ i /dτ . We then write the second derivatives of ψ(τ ) with respect to time, for the real part and,ψ the imaginary part respectively. We use the equations forψ 1,2 in Eqs. (7), together with Eqs. (8), (9) to solve for φ 0,3 , andφ 0,3 , which we then eliminate from the equations forφ 0,3 in Eqs. (7). From this set of four equations we obtain explicit expressions for δ 0 (τ ) and δ 3 (τ ) in terms of Ω 01 (τ ), Ω 23 (τ ) plus known quantities, together with a pair of coupled first-order differential equationṡ where f i ≡ f i (τ ) and g i ≡ g i (τ ) are functions defined solely in terms of φ 1 , φ 2 and their derivatives. ‡ The corresponding physical Hamiltonian is H ′ (t) = 1 2 αH(αt).
Perpetual emulation threshold of PT -symmetric Hamiltonians PSfrag replacements The two time-dependent couplings Ω 01 (τ ) (black) and Ω 23 (τ ) (blue) required for emulation, as a function of dimensionless time τ /π for parameters γ = 0.5λ, The corresponding sink and source populations p 0 = |φ 0 | 2 (maroon) and p 3 = |φ 3 | 2 (green). At time τ = 41.09π (indicated with the vertical dashed line), the coupling Ω 23 diverges as the source level is depleted. This breakdown is a common feature of all solutions without amplitude recycling.

No amplitude recycling
While Eqs. (10) must in general be solved numerically, in the absence of amplitude recycling, Ω 03 = 0, they decouple and admit analytic solutions. Here we focus on initial conditions with equal amplitude in each of the target PT symmetric states, giving and where Ω 01 (τ ) = Ω 23 (−τ ) (so that when time is reversed the source acts as a sink and vice versa), and δ 0 (τ ) takes the upper sign and δ 3 (τ ) the lower. The dynamic couplings are plotted for representative parameters in Fig. 2a, making clear that the coupling Perpetual emulation threshold of PT -symmetric Hamiltonians  The breakdown time with amplitude recycling for a range of amplitude recycling strengths Ω 03 and for fixed γ = 0.5λ. For large initial coupling strength Ω init , the amplitude recycling field has little effect. As Ω init is lowered, however, we see a dramatic divergence of the breakdown time which corresponds to the onset of periodic behaviour. Parameters not explicitly mentioned are the same as in Fig. 2 Ω 23 diverges. Physically, this break down arises due to the depletion of source level, as shown in Fig. 2b. As the population of this level decreases, a stronger coupling Ω 23 is required to maintain the constant probability flux at rate γ into the system, and this necessarily stops when the level is depleted.
Let us denote as τ * the time at which Ω 23 diverges, the maximum time we can expect our emulation to run. We obtain τ * as the smallest time at which the denominator in Eq. (11) vanishes, and plot it in Fig. 3a as a function of the starting coupling Ω init for different values of γ. For γ in the range 0.2λ γ 0.8λ, changing γ does not alter the τ * very much. However, increasing γ above 0.8λ leads to a significant drop in maximum emulation time; we observe similar behaviour for decreasing γ 0.2. We can obtain a simple approximation to the breakdown time by setting the oscillating terms in the denominator of Eq. (11) to zero. This gives valid for Ω init ≪ 1. In the opposite regime, Ω init ≫ 1, we obtain an alternative approximation by expanding the denominator to first order in time, and setting the result equal to zero. Both these expressions are good approximations in their respective regions of validity as shown in Fig. 3a.
Perpetual emulation threshold of PT -symmetric Hamiltonians

8
From the limits of Eqs. (13), (14), as well as from the full results in Fig. 3a, we see that the dominant behaviour of the breakdown time τ * is an approximate scaling with Ω −2 init . Thus, in principle, we can always arrange our initial coupling strength to enable the emulation to cover any time interval of interest. However, making Ω init arbitrarily small is not without cost: the initial values of the source and sink populations required by the emulation are Since the total probability is not contained in the PT -symmetric subspace, we define the PT -symmetric fraction in terms of the populations p i = |φ i | 2 as This expresses the minimum probability of finding the system in the PT -symmetric subspace. This definition is motivated by an experimental resource limit since detecting the probability distribution in all states becomes increasingly cumbersome as PT -symmetric fraction gets smaller. The PT -symmetric fraction is equal to r = min In the limit Ω init /γ ≪ 1, this gives r ∼ Ω 2 init /γ 2 . Thus, decreasing Ω init to extend the emulation time necessarily leads to a corresponding decrease in the PT -symmetric fraction.

Amplitude recycling and the perpetual emulation threshold
We now consider the effects of the amplitude recycling coupling field. In this case, Eqs. (10) for Ω 01 and Ω 23 do not decouple and we obtain our results through numerical integration. Figure 3b shows the breakdown time τ * of the emulation as a function of the initial coupling Ω init for several values of the amplitude recycling strength. For large Ω init , the amplitude recycling does not significantly affect the breakdown time. However, as Ω init decreases, the breakdown time increases until it rapidly diverges at a value that increases with the amplitude-recycling strength. Figure 4 shows the populations of the sink and source levels as a function of time for several values of the amplitude recycling strength Ω 03 and for fixed initial coupling Ω init = 0.05λ. As Ω 03 is initially increased from zero, the decline of p 3 with time becomes shallower. The solution nevertheless terminates with p 3 = 0, as in the Ω 03 = 0 case.  This trend continues until the critical value Ω 03 = 0.01412λ of the amplitude recycling strength is reached and the cycle-averaged gradients of p 0 and p 3 are both zero. At this point, the behaviour of the populations is purely periodic; they never drop to zero and the coupling Ω 23 does not diverge. For yet larger values of Ω 03 the populations reverse order such that p 3 > p 0 , and, in the long-time limit, become periodic once again. The emulated PT -symmetric dynamics are correct for all times when Ω 03 exceeds its critical value.

Perpetual emulation threshold of PT -symmetric Hamiltonians
In Fig. 5 we explore the dynamics of the coupling functions and plot Ω 01 (τ ) as a function of Ω 23 (τ ). Without amplitude recycling (Fig. 5a) the trace simply diverges.
With Ω 03 equal to its critical value (Fig. 5b) the functions form a single closed orbit which illustrates the periodic nature of this solution. Then, for higher values of Ω 03 (Fig. 5c), after some initial transients, the behaviour collapses onto a closed orbit. We note that these closed orbits do not represent limit cycles because different initial values of the couplings give rise to different oscillating states. Finally, in Fig. 5d, we show results for a value of Ω 03 far in excess of its critical value. In this case, the long-time orbit develops a sharp kink feature, which is a sign of the pronounced lack of timereversal symmetry in the coupling functions in this regime.
We have shown how the addition of amplitude recycling can elicit one of two distinct responses. Below some critical coupling value, the emulation still terminates, albeit now with a larger breakdown time cf. the case without amplitude recycling. Above the critical coupling, however, the amplitude-recycling coupling field stabilises the dynamics, PSfrag replacements such that the emulation enters a purely periodic mode that can run indefinitely; a situation that only occurs when γ = 0 in the absence of amplitude recycling. A clear boundary between the two regimes exists and is shown in Fig. 6 as a function of Ω init and Ω 03 .

Experimental implementation and constraints
We can readily implement our amplitude recycling scheme to coherently emulate a twolevel PT -symmetric Hamiltonian using ultracold atoms. Let us consider a 87 Rb Bose-Einstein condensate in its electronic groundstate 5 2 S 1/2 . 87 Rb is an alkali atom with nuclear spin I = 3/2 so that the groundstate splits into two hyperfine levels, F = 1, 2.
We identify the PT -symmetric subspace with the |F, m F = |1, −1 and |1, 0 levels and the reservoir levels with |1, 1 and |2, 0 levels. The |2, 0 level can be connected to the stretched levels of the F = 1 manifold via two microwave transitions that serve as the coupling to one of the reservoirs and the amplitude-recycling field. However, binary collisions between the two hyperfine manifolds might limit the lifetime of the atomic cloud to 100 ms. For realistic rf coupling parameters this is in the "long emulation time" limit, and does not degrade the utility of this approach for simulating PT -symmetric Hamiltonians. Alternatively, the continuous dynamical decoupling techniques described in [24,25] allow this to be equally well realized all in the F = 2 hyperfine manifold. Ideally, we would like to be able to emulate the PT -symmetric system with arbitrary parameters for as long a time as desired. While this is possible in principle, both with and without amplitude recycling, there exist a number of aspects inherent to this scheme that provide constraints on what is possible in practice. The appearance of the timescale τ sharp ≪ 1, discussed in Sec. 2, represents a challenge since for γ → λ the period of the population oscillations diverges, whilst τ sharp remains fixed. This timescale is directly reflected in the detunings. Expanding Eq. (12) about one of its maxima, we obtain δ 0 (τ ) ≈ (λ + γ) [1 − γ(λ + γ)(τ − τ min ) 2 /2], which exhibits a timescale [γ(γ + λ)/2] −1/2 ∼ 2/γ = √ 2τ sharp . The corollary is that experimental control must be able to simultaneously cover both small, √ 2/γ, and large, unity, time scales. Imperfect preparation of the initial state has a similar effect in the faithful reproduction of the PT -symmetric dynamics. Only the squared modulus of the wavefunction, i.e., the level populations, needs to be taken into account when preparing the initial states. The required phases can be absorbed into the driving fields through an appropriate gauge transformation, which means that control parameters do not need to be recalculated for different initial phases. An initial error of a few percent in state preparation of one of the PT -symmetric levels induces a small discrepancy between the realised population dynamics and the target. Due to the larger absolute population of source and sink levels, an error of a few percent in their preparation will cause larger subsequent errors in the dynamics. This discrepancy, however, can be made smaller if the preparation error is common-mode, where any offset in initial population is shared by both source and sink levels. In any case, we can eliminate a state preparation error with an initial weak measurement right after state preparation takes place. The The solid lines indicate the parameter values for which we obtain a simulation break-down time τ * = 10π for the indicated values of γ. In the regions above these curves, emulations last long enough to capture a minimum of ten Rabi oscillations of the PT system. On the other hand the dashed lines indicate parameters for which r = r min = 1/10, which from Eq. (17) are the straight lines Ω init = γ/3. Parameters to the right of these lines lead to values of the PT -fraction r > r min . To obtain feasible emulations then, we require that these two lines cross for a given γ. For the parameters considered here (same as Fig. 2), this is the case for γ/λ = 0.5, 0.85, 0.875, but not for γ/λ = 0.9.
weak measurement does not affect the initial populations, but provides the necessary information for recalibrating the dynamics of the emulation. This way, the timedependence of the coupling fields and detunings is dictated by the measured initial state but the subsequent dynamics always follow that of a PT -symmetric Hamiltonian albeit with slightly different parameters than the intended ones. Similarly, weak measurements performed while the dynamics take place can also be used to correct for errors during the emulation.
Further constraints arise from the emulation scheme. The main practical constraint limiting the total emulation time without amplitude recycling is the fact that longer times require larger populations of the source and sink levels, and this reduces the PT -symmetric fraction r. Small values of r demand a large dynamic range from the experimental detection scheme. In practice, the dynamic range of the detectors sets the minimum signal-to-noise ratio necessary for the chosen emulation parameters. This situation however can be improved dramatically by amplitude recycling.
The equations in Eq. (7) used to determine the initial wave function components of the PT -symmetric subspace at time τ = 0, are independent of the amplitude recycling strength Ω 03 . N 0 and therefore r are unaltered by the addition of the amplituderecycling field and Eq. (17) applies irrespectively of Ω 03 . Thus, the emulation times can PSfrag replacements ∆δ 0 /λ be extended by the amplitude recycling without changing the PT -symmetric fraction. The minimum value of r = r min , set by experimental constraints, gives a minimum value for the ratio: through Eq. (17) and this applies independently of Ω 03 . Figure 7 illustrates how this constraint can be taken into account to determine which values of Γ are experimentally accessible for a given r min . Assuming we desire to simulate for a fixed time of 10πα −1 , say, we require a breakdown time τ * > 10π. This requirement defines an accessible region in the parameter space of Ω init and Ω 03 ; the solid lines in the Fig. 7 delineate the boundary of this region for various values of γ. Let us also assume a minimum PT -symmetric fraction of r min = 1/10. The constraint Eq. (18) then also defines a region in the same parameter space, and these boundaries are plotted with dashed lines. Clearly, for an emulation to both run long enough and have a sufficiently large r, we require that the two regions overlap. Whether this happens or not for a given γ can be seen from Fig. 7 as, for a feasible emulation, the two boundary lines will cross. This happens for γ/λ = 0.5 for example, but not for γ/λ = 0.9. Indeed, from this plot we can obtain the maximum value of γ that it is possible to emulate. If we restricted ourselves to Ω 03 /λ ≤ 2 such that the parameter range shown in Fig. 7 is all that is accessible, then we should be able to emulate the range 0 ≤ γ/λ ≤ 0.875. By way of contrast, we can use the analytical approximation of Eq. (13) to obtain the accessible range without amplitude recycling, and, for the same parameters, this turns out to be 0 ≤ γ/λ ≤ 0.257. Use of the amplitude recycling therefore significantly enlarges the range of possible γ values that can be explored.
One further constraint comes from considering the change in size of the dynamic couplings and detunings. This is not such a problem for the couplings Ω 01 and Ω 23 as these are made small by construction to extend the simulation times. On the other hand, it is desirable for the detunings to also remain small but there is no a priori reason why this should be the case. Without amplitude recycling, Eq. (12) shows that the detunings are always bounded as −λ − γ ≤ δ 0,3 (τ ) ≤ −λ + γ. These bounds change however, when we turn on the amplitude recycling. Fig. 8 shows the size of the range of the detuning ∆δ 0 ≡ |max(δ 0 ) − min(δ 0 )| as a function of amplitude recycling strength Ω 03 .
The overall trend is that the range of values that the detunings take increases as we increase the amplitude-recycling strength. The same quantity for detuning δ 3 shows little variation and is always order 2γ. Thus, to avoid requiring large detunings, we need to operate with as small Ω 03 as possible, in which case we obtain ∆δ 0 ∼ 2γ. We note that although we show the range of δ 0 here, its increase is largely due to increase in magnitude of the minimum value, |min δ 0 |.

Conclusions and outlook
Coherent emulation of PT -symmetric or other non-Hermitian Hamiltonians can be used for gaining an understanding of these systems. In our scheme, a single unitary system is partitioned into two subspaces coupled together with coherent couplings. Even though the couplings are coherent, we can make them resemble particle gain and loss by appropriate manipulation of the dynamical equations describing them. The subspace of interest is therefore made to mimic the dynamical evolution of an open system coupled to two reservoirs. The coherent couplings between the embedded system and the 'environment' can be controlled with standard experimental techniques [26]. This way, the effective decoherence and dynamic timescales can be adjusted arbitrarily to suit the experimental constraints. Importantly, the choice of emulation parameters also dictates the properties of the reservoirs themselves.
Our results show how a closed, unitary system with no intrinsic gain or loss can emulate the dynamics of PT -symmetric Hamiltonians. In the future it will be interesting to extend these results to a wider class of non-Hermitian Hamiltonians using arbitrarily connected spin Hamiltonians as emulators.