General description of quasi-adiabatic dynamical phenomena near exceptional points
Thomas J. Milburn, Jörg Doppler, Catherine A. Holmes, Stefano Portolan, Stefan Rotter, Peter Rabl
PPreprint of September 6, 2018—not for dissemination
General description of quasi-adiabatic dynamical phenomena near exceptional points
Thomas J. Milburn, ∗ J¨org Doppler, Catherine A. Holmes, Stefano Portolan, Stefan Rotter, and Peter Rabl Institute of Atomic and Subatomic Physics, Vienna University of Technology (TU Wien), Stadionallee 2, 1020 Vienna, Austria Institute for Theoretical Physics, Vienna University of Technology (TU Wien), 1040 Vienna, Austria University of Queensland, School of Mathematics and Physics, QLD 4072, Australia University of Southampton, Department of Physics and Astronomy, Southampton SO17 1BJ, United Kingdom
The appearance of so-called exceptional points in the complex spectra of non-Hermitian systems isoften associated with phenomena that contradict our physical intuition. One example of particularinterest is the state-exchange process predicted for an adiabatic encircling of an exceptional point.In this work we analyse this and related processes for the generic system of two coupled oscillatormodes with loss or gain. We identify a characteristic system evolution consisting of periods ofquasi-stationarity interrupted by abrupt non-adiabatic transitions, and we present a qualitative andquantitative description of this switching behaviour by connecting the problem to the phenomenonof stability loss delay. This approach makes accurate predictions for the breakdown of the adiabatictheorem as well as the occurrence of chiral behavior observed previously in this context, and providesa general framework to model and understand quasi-adiabatic dynamical effects in non-Hermitiansystems.
I. INTRODUCTION
The quantum adiabatic theorem is a seminal result inthe history of quantum mechanics. Paraphrasing Born,the theorem states that for an infinitely slow paramet-ric perturbation there is no possibility of a quantumjump [1]. Many physical phenomena observed in bothquantum and classical systems can be explained by thistheorem, ranging from optical tapers [2] to robust quan-tum gates [3]. Recently, the applicability of adiabaticprinciples to non-Hermitian systems, e.g., coupled har-monic modes with gain or loss, has attracted consider-able attention. Here, the complex eigenvalue structureand the existence of so-called exceptional points (EPs)leads to new counterintuitive phenomena [4–18]. Per-haps most strikingly, adiabatically encircling an EP waspredicted to effect a state-exchange, with applications forswitching and cooling [20–22]. However, it is now knownthat the very presence of non-Hermiticity prevents a gen-eral application of the adiabatic theorem [23–25], and theinevitability of non-adiabatic transitions leads to new ef-fects, e.g., to chiral behavior [26–31].Whereas the above results point to fascinating newphysical phenomena, the complexity of the problemmostly requires one to resort to numerical studies (ascited above) or to focus on limiting cases where the sys-tem evolution is eventually dominated by a single modewith maximum gain or minimum loss. An importantstep beyond this limitation has recently been presentedin Refs. [27] and [32] in which an exactly solvable modelis considered and a connection between the appearanceof non-adiabatic transitions and the Stokes phenomenonof asymptotics [33] is thereby found. However, even forvery simple scenarios, these exact case studies are math-ematically already quite involved, and the translation of ∗ [email protected] the observed dynamics to other systems, in particularto realistic systems with imperfections and noise, is notimmediately obvious.In this work we analyse quasi-adiabatic dynamics innon-Hermitian systems near EPs with the aim to providea generalised framework for both modelling and under-standing the associated dynamical phenomena. Our ap-proach reveals that the solutions are in general composedof periods of quasi-stationary during which the solutionfollows fixed points, interrupted by abrupt non-adiabatictransitions due to the exchange of stability. However, thetime of these transitions cannot be predicted by a stan-dard stability analysis, and, intriguingly, we find thatpiece-wise adiabaticity is still a key ingredient for under-standing the evolution of the system in spite of an overallbreakdown of adiabatic principles.On a more fundamental level, our analysis shows thatthe quasi-adiabatic dynamics near an EP is a singlu-arly perturbed problem [34], meaning that, in contrastto Hermitian systems, the dynamics cannot be obtainedby perturbative corrections to the adiabatic prediction.This fact makes adiabatic principles in non-Hermitiansystems particularly interesting as well as challenging tounderstand, both from a physical and from a mathe-matical point of view. Specifically, here we connect theproblem of non-adiabatic transitions to the more gen-eral phenomenon of stability loss delay [35, 36] in dy-namical bifurcations. This concept more easily affordsintuition in complicated examples where exact solutionscannot be found and in realistic systems where noise can-not be ignored. Our results are therefore important fora variety of modern-day experiments with, e.g., waveg-uides [15, 16], coupled resonators [17, 18], semiconduc-tor microcavities [19], or electromechanical [37, 38] andoptomechanical systems [39–41], which offer sufficientlyhigh control for the observation of the predicted dynam-ical phenomena.Typeset by REVTEX a r X i v : . [ qu a n t - ph ] J un FIG. 1. (Color online.) (a) Cartoon of two coupled har-monic modes with gain or loss. (b) Example parametric pathwhere γ is fixed, ω = r sin φ ( t ), and g = γ/ r cos φ ( t ).(c) Real ( (cid:60) ) and imaginary ( (cid:61) ) parts of the spectrum, λ ∓ = ∓ (cid:112) ( ω + iγ/ + g . The curve is the trajectory of λ − forthe path defined in (b) and depicts the fully adiabatic evolu-tion. II. NON-HERMITIAN DYNAMICS ANDEXCEPTIONAL POINTSA. Model
For the following discussion we consider the genericmodel of two coupled harmonic oscillators with frequen-cies ω and ω , decay rates γ and γ , and couplingstrength g ; see Fig. 1(a). The equations of motion forthe amplitudes α and α are ddt (cid:18) α α (cid:19) = − i (cid:18) ω − iγ / gg ω − iγ / (cid:19) (cid:18) α α (cid:19) , (1)where in general ω i = ω i ( t ), γ i = γ i ( t ), and g = g ( t )are functions of time. For the following analysis it isconvenient to eliminate the common evolution with av-erage frequency Ω := ( ω + ω ) / γ + γ ) / β and β via (cid:18) α ( t ) α ( t ) (cid:19) = e − i (cid:82) t [Ω( t (cid:48) ) − i Γ( t (cid:48) ) / dt (cid:48) (cid:18) β ( t ) β ( t ) (cid:19) . (2)The remaining non-trivial dynamics in this frame is ddt (cid:18) β β (cid:19) = − i (cid:18) − ω − iγ/ gg ω + iγ/ (cid:19) (cid:18) β β (cid:19) , (3)where ω := ( ω − ω ) / γ := ( γ − γ ) /
2. Note thatwhile the global transformation (2) does not affect anyof the following results, if Γ (cid:54) = 0 then the experimen-tally observable amplitudes α , are related to β , by anexponentially large or small prefactor. Below we suppose that at least ω and g , or ω and γ canbe controlled as a function of time. This can be achieved,e.g., with optical modes propagating through waveguideswith spatially varying losses [42, 43], by applying chirpedlaser pulses to molecular systems [28], or by using two me-chanical resonators with electrically [37, 38] or optome-chanically [39–41] controlled parameters. B. Exceptional points
Let us write Eq. (3) more compactly as ˙ (cid:126)x = − i M (cid:126)x ,where (cid:126)x is the state vector and M is the dynamical ma-trix, or sometimes called in this context a non-HermitianHamiltonian [9], i.e., (cid:126)x := (cid:18) α α (cid:19) and M := (cid:18) − ω − iγ/ gg ω + iγ/ (cid:19) . (4) M has eigenvalues λ ∓ = ∓ λ = ∓ (cid:112) ( ω + iγ/ + g .The corresponding eigenvectors are (cid:126)r − = (cid:18) cos ϑ/ ϑ/ (cid:19) and (cid:126)r + = (cid:18) − sin ϑ/ ϑ/ (cid:19) (5)with ϑ such that tan ϑ = − g/ ( ω + iγ/ (cid:60) ) and imaginary ( (cid:61) ) parts of λ ± asa function of g and ω with γ fixed. The pinch points ω + iγ/ ∓ ig = 0 are EPs [4–10]. At these points theeigenvalues as well as the eigenvectors coalesce, and M becomes non-diagonalizable. Encircling an EP with aclosed path in parameter space causes the two eigenval-ues, and hence also the two eigenvectors, to swap; seeFigs. 1(b, c). Based on intuition from the quantum adia-batic theorem, it was suggested that this unique featurecould be observed in physical systems by encircling anEP over a time T such that T | λ − − λ + | is large [20–22].However, other studies contradict this result and showthat due to non-Hermiticity this picture cannot hold ingeneral [23–31]. C. Numerical examples
Before presenting a further analytic treatment ofEq. (3) we consider in Fig. 2 some typical solutions forencircling an EP with T | λ − − λ + | (cid:29)
1. For these exam-ples we choose a path in parameter space as defined inFig. 1(b). We expand the solution as (cid:126)x ( t ) = c − ( t ) (cid:126)r − ( t ) + c + ( t ) (cid:126)r + ( t ) , (6)where (cid:126)r − ( t ) and (cid:126)r + ( t ) are the instantaneous eigenvectorsof M ( t ), and we choose the initial condition c − (0) = 1and c + (0) = 0. The adiabatic prediction is c ad. − ( t ) (cid:39) e − i (cid:82) t λ − ( t (cid:48) ) dt (cid:48) and c ad.+ ( t ) (cid:28) c ad. − ( t ) [44–46].In examples (i) and (ii) in Fig. 2 we have chosenan anticlockwise and a clockwise encircling respectively, FIG. 2. (Color online.) Plots of typical numerical solutions of Eq. (3) for the path defined in Fig. 1(b) with initial eigenvectorpopulations c − (0) = 1 and c + (0) = 0. For the function φ we choose φ ( t ) = ± πt/T in examples (i, ii), and we choose φ ( t ) = − πt/T + π in example (iii). In all cases we set r = 0 . γ = 1 and T = 45, for which T | λ − − λ + | (cid:29)
1. The toprow shows the dynamical gain parameter, T (cid:61) λ − ( t ), and the total integrated gain, (cid:82) t (cid:61) λ − ( t (cid:48) ) dt (cid:48) . Note that the dynamicalgain is the gain of the adiabatic prediction but not necessarily the actual gain of the numerical solution. The middle rowshows the eigenvector populations, | c ∓ ( t ) | , along with the adiabatic prediction, | c ad . − ( t ) | . We do not plot | c ad . + ( t ) | becauseadiabatic principles imply | c ad . + ( t ) | (cid:28) | c ad . − ( t ) | . The bottom row shows a projection of the numerical solution onto the realand imaginary parts of the eigenspectrum, specifically [ | c − ( t ) | λ − ( t ) − | c + ( t ) | λ + ( t )] / [ | c − ( t ) | + | c + ( t ) | ]. The use of red andblue is to provide an indication of which population, or surface, corresponds to a gain and loss eigenvector respectively. φ ( t ) = ± πt/T . In the anticlockwise example the so-lution matches the adiabatic prediction and the corre-sponding state flips, but in the clockwise example we ob-serve a non-adiabatic transition, for which, apart from anoverall amplification, the system returns to the originalstate. This chiral behavior, first presented in Ref. [26], il-lustrates one of the key differences between the dynamicsin Hermitian and non-Hermitian systems. In the latter,the eigenvalues are complex, which causes gain or loss in c − and c + . An infinitesimally small non-adiabatic cou-pling can therefore be exponentially amplified, causingthe gain eigenvector to dominate. This mechanism intu-itively explains why the adiabatic theorem does not ingeneral hold for non-Hermitian systems.Example (iii) shows the result for a more interestingpath φ ( t ) = − πt/T + π where gain-loss behavior swapshalf-way through and the total integrated dynamical gainvanishes, (cid:82) T (cid:61) λ ( t ) dt = 0. Surprisingly, the final statematches the adiabatic prediction, | c − ( T ) | (cid:39) | c − (0) | ,even though during the interim the solution is highlynon-adiabatic. This observation cannot be explained bythe intuitive argument above because c − is non-trivially slaved to c + past the time t = T / c − to increase exponentially. Thus, considering dy-namical gain alone is insufficient to accurately predictbehaviour for quasi-adiabatic dynamics near EPs.These basic examples illustrate that the dynamics ofnon-Hermitian systems involves three characteristic ef-fects:(i) The swapping of eigenvectors due to a 4 π -periodicity about an EP, which follows from thetopology of the complex eigenvalue spectrum.(ii) The appearance of enhanced non-adiabatic transi-tions due to the presence of gain or loss.(iii) Periods of adiabatic evolution that persist signifi-cantly beyond the time of stability loss.While (i) is readily incorporated by the eigenvector de-composition (6), we will now develop a general approachto describe the non-trivial interplay between (ii) and (iii). III. DYNAMICAL ANALYSISA. Relative non-adiabatic transition amplitudes
In order to develop a general dynamical description weconsider the evolution operator, U ( t ), defined by (cid:126)x ( t ) = U ( t ) (cid:126)x (0), which contains the full dynamics independentof the initial condition. In the eigenbasis Eq. (5), U ( t ) isthe solution of˙ U = − i (cid:18) − λ ( t ) − f ( t ) f ( t ) λ ( t ) (cid:19) U , U = (cid:18) U − , − U − , + U + , − U + , + (cid:19) (7)with initial condition U (0) = , where f ( t ) = g ( t ) (cid:2) ˙ ω ( t ) + i ˙ γ ( t ) / (cid:3) − (cid:2) ω ( t ) + iγ ( t ) / (cid:3) ˙ g ( t )2 iλ ( t ) , (8)is the non-adiabatic coupling [44]. Adiabaticity usu-ally requires that the non-adiabatic coupling be muchsmaller than the distance between eigenvectors, ε ( t ) := | f ( t ) / [2 λ ( t )] | (cid:28)
1. Since ε ( t ) ∝ T − this condition isalways satisfied for an appropriate T . To set f ( t ) = 0in Eq. (7), which would imply ε ( t ) = 0, would yield thediagonal adiabatic prediction U ad. ( t ) = (cid:32) e − i (cid:82) t λ ( t (cid:48) ) dt (cid:48) e i (cid:82) t λ ( t (cid:48) ) dt (cid:48) (cid:33) . (9)However, as is evident in Fig. 2, even for arbitrarily smallyet non-vanishing ε ( t ) the actual solution is significantlynon-diagonal. This indicates that the system is singularlyperturbed by the non-adiabatic coupling, and U ( t ) cannotbe obtained as a perturbative correction to U ad. ( t ). Weshall henceforth call ε ( t ) (cid:28) quasi -adiabatic condi-tion (see Appendix A 1 for more details).In order to describe the non-adiabatic character of U ( t )for quasi-adiabatic dynamics we focus on the relativenon-adiabatic transition amplitudes [26]: R − ( t ) := U + , − ( t ) U − , − ( t ) and R + ( t ) := U − , + ( t ) U + , + ( t ) . (10)These describe the amount of non-adiabaticity in the so-lution. For example, R − ( t ) is a measure of the magnitudeof the net non-adiabatic transition from (cid:126)r − ( t ) to (cid:126)r + ( t ).If R ∓ ( t ) (cid:28) c ∓ is behaving adia-batically, while R ∓ ( t ) (cid:29) R ∓ ( t ) considered as a dynamical variable is thesolution to the Riccati equation [32, 47]˙ R ∓ = ∓ iλ ( t ) R ∓ ∓ if ( t )(1 + R ∓ ) (11)with initial condition R ∓ (0) = 0. Dynamical phenomenaassociated with quasi-adiabatically encircling EPs canthus be understood from the solutions of this equationin the limit ε ( t ) (cid:28)
1. Note that the equations of motion
FIG. 3. (Color online.) (a) Plot of (cid:61) λ ( t ) (upper panel), anda typical solution for R ≡ R − (lower panel), for the pathdefined in Fig. 1(b) with φ ( t ) = − πt/T + π . Note that (cid:61) λ ( t ) = −(cid:61) λ − ( t ), which is plotted in Fig. 2. We have cho-sen r = 0 . γ = 1 and T = 100, for which ε ( t ) (cid:39) . | R ad. ( t ) | and | R n.ad. ( t ) | respectively. The shaded area isone standard deviation about the mean of R − obtained from10 ,
000 stochastic numerical integrations of c − and c + (seeSec. IV for more details). (b) Cartoons of the global phaseportraits of the equation of motion for R near t ∗ . Arrows de-note the direction of time-evolution along an integral curve.The fixed point near the origin corresponds to R ad. ( t ), andthe fixed point far from the origin corresponds to R n.ad. ( t ). for R − and R + are related via R − ↔ /R + . In the fol-lowing, we therefore consider only R := R − without lossof generality.We remark that, assuming transients are damped, therelation R − ↔ R + has the immediate consequence thatlim t →∞ R − ( t ) R + ( t ) = 1, which agrees with Ref. [26] andprohibits simultaneous adiabatic behaviour in both c − and c + over long times. B. Fixed points and stability loss delay
The lower panel of Fig. 3(a) shows a generic solution for R during multiple quasi-adiabatic encirclements of an EP(see the caption for details). It resembles a square wave,i.e., we see fast switching between two quasi-stationaryvalues. This behavior can be understood from a sepa-ration of time-scales in Eq. (11). For short times theslowly varying parameters λ ( t ) (cid:39) λ and f ( t ) (cid:39) f can beconsidered constant and˙ R (cid:39) − iλR − if (1 + R ) . (12)On a fast time-scale set by |(cid:61) λ | − the solution thereforeapproaches one of two fixed points R ad. = − λf (cid:32) − (cid:112) λ − f λ (cid:33) (cid:39) − f λ and R n.ad. = − λf (cid:32) (cid:112) λ − f λ (cid:33) (cid:39) − λf . (13)The first fixed point, R ad . ( t ) ∝ ε ( t ) (cid:28)
1, indicates adia-batic behavior ( c − dominates) and is stable for (cid:61) λ ( t ) <
0. The second fixed point, R n . ad . ( t ) ∝ ε − ( t ) (cid:29)
1, in-dicates a non-adiabatic transition has occured ( c + domi-nates) and is stable for (cid:61) λ ( t ) >
0. These two fixed pointsare plotted in Fig. 3(a). Evidently, the periods of quasi-stationarity there exhibited correspond to following oneof these two fixed points.On a slow time-scale set by T the parameters λ ( t ) and f ( t ) may change considerably and at certain critical timesthe stability of the two fixed points swaps. For example, R ad. ( t ) becomes unstable and R n.ad. ( t ) stable when thesign of (cid:61) λ ( t ) becomes positive. Let us denote the criticaltimes by t ∗ , which are marked in Fig. 3(a). Na¨ıvely, onemight expect an immediate rapid transition between theneighbourhoods of R ad. ( t ) and R n.ad. ( t ) upon passing acritical time t ∗ , but, as is evident in Fig. 3(a), this isnot the case. Instead, we see that the solution follows,e.g., R ad. ( t ) while it is unstable for a significant amountof time; the loss in stability is delayed. Intuition forthis behaviour is obtained from the phase portraits ofEq. (11), shown in Fig. 3(b). The local phase portraitabout R ad. ( t ) goes from a spiral towards R ad. ( t ) for t < t ∗ to a spiral away from R ad. ( t ) for t > t ∗ , passing through adegenerate bifurcation at t = t ∗ when R ad. ( t ) is a centreand neither stable nor unstable. We therefore expectsome persistence in the following of R ad. ( t ) because near t ∗ it is only ‘weakly’ stable or unstable.To illustrate the existence of a significant delay be-tween the critical time t ∗ and the actual time of a non-adiabatic transition t + we consider the specific path de-fined in Fig. 1(b) with φ ( t ) = − πt/T + π and r (cid:28) γ .This is a good model for the numerical solution shownin Fig. 3(a). Then, λ ( t ) (cid:39) i √ rγe − iπt/T , f ( t ) (cid:39) iπ/ (2 T ),and ε ( t ) (cid:39) ε = π/ (4 √ rγ T ). Let us focus on the lossof stability of R ad. ( t ) at t ∗ = 3 T /
2. Assuming that thesystem is near R ad . ( t ) we can neglect the nonlinear termin Eq. (11): ˙ R √ rγ (cid:39) e − iπt/T R + ε. (14)The particular integral of this equation is found to be R ( t ) = − i E (cid:18) i ε e − iπt/T (cid:19) e i exp( − iπt/T ) / (2 ε ) , (15) where E is the exponential integral. Since ε (cid:28) E (see, e.g., 5.1.7and 5.1.51 in Ref. [48]) to obtain R ( t ) (cid:39) − εe iπt/T − iε e iπt/T + . . . − π (cid:20) cos (cid:18) πtT (cid:19)(cid:21) e i exp( − iπt/T ) / (2 ε ) . (16)The first two terms (upper line on the right) correspondto following R ad . ( t ) with higher order corrections. Thethird term (lower line) is negligible for ( t − t ∗ ) < T / t ∗ = 3 T / t − t ∗ ) > T /
2, thereby indicating a non-adiabatic transi-tion. Thus, under the ideal conditions assumed here andgiven that the solution has approached R ad. ( t ) by t = t ∗ ,the delay in the loss of stability is ( t + − t ∗ ) = T / R ad. ( t ) is stable for the entire loop aroundthe EP and therefore the solution follows the adiabaticprediction, | c + ( t ) /c − ( t ) | (cid:39) | R ad. ( t ) | . In contrast, in (ii) R ad. ( t ) is always unstable and a non-adiabatic transitionoccurs. In (iii) the solution first switches from R ad. ( t ) to R n.ad. ( t ), but then back again with a delay t + (cid:46) T / R ad. ( t ) becomes stable at t = T /
2. Note that the de-lay times exhibited in the first encircling period as seen inFigs. 2 and 3(a) differ somewhat from the value t + = T / R (0) = 0, which is not exponentiallyclose to R ad . (0), and therefore effects a transient termof the form Ae i exp( − iπt/T ) / (2 ε ) . After about one encir-cling period the system approaches the unique long-timerelaxation oscillation which is a universal signature ofquasi-adiabatically encircling EPs.We finish this section with a remark on the relationbetween the above results and the Stokes phenomenonof asymptotics , i.e., the switching-on of exponentiallysuppressed terms in asymptotic expansions [33]. InRefs. [27, 32] an exact solution for the example consid-ered in this subsection is presented (using r (cid:28) γ butnot neglecting the nonlinearity), which we review in Ap-pendix B. In this exact solution one sees that the sharp(but continuous) transition, which in Eq. (16) is repre-sented by the signum function, is precisely the Stokesphenomenon of asymptotics, leading here to a breakdownof the adiabatic theorem. In our current approach, whichwe elaborate further in the next section, this disconti-nuity is connected to the problem of stability loss de-lay. To our knowledge the connection between the Stokesphenomenon of asymptotics and stability loss delay hashitherto not been suggested, and might be worth explor-ing further. However, here we will leave such considera-tions aside and proceed with a pragmatic generalisationof these initial results to arbitrary paths in parameterspace. C. Generalised quasi-adiabatic solution
In the previous section, Sec. III B, we were able to un-derstand the generic solution exhibited in Fig. 3(a) froma separation of time-scales, which resulted in a delay inthe loss of stability of the instantaneous fixed points. Infact, slow-fast systems with dynamical bifurcations are asubject of current mathematical interest. The reader isreferred to Ref. [34] for a concise description. The reasonthat the critical times do not coincide with the observedtimes when an instantaneous fixed point loses stabilityis because our slow-fast system is singularly perturbed ;the slow system is described by an algebraic equation,while the fast system by a differential equation. Onemust therefore resort to non-standard analysis. A prin-cipal result of the non-standard analysis of slow-fast sys-tems is the existence of stability loss delay about certaindynamical bifurcations [35, 36, 49], which we observedexplicitly in Sec. III B. In the following we build uponthis to construct a generalised quasi-adiabatic solution,which, additionally, affords an estimation of delay times.We are interested in solutions that for times near criti-cal times t ∗ are in the vicinity of a fixed point. We there-fore begin by looking at the zero-crossings of (cid:61) λ ( t ), whichdetermine t ∗ . For some window [ t − , t + ] about each t ∗ ,i.e., t − < t ∗ < t + , we seek a solution R t ∗ ( t ) that follows R ad. ( t ) or R n.ad. ( t ). Since transitions between R ad. ( t )and R n.ad. ( t ) are very quick, by making a piece-wise ad-dition of segments that follow one or the other fixed pointwe arrive at the approximation for the complete solutionthus R ( t ) (cid:39) (cid:88) t ∗ [Θ( t − t − ) − Θ( t − t + )] R t ∗ ( t ) , (17)where Θ is the Heaviside step function.Let us now consider a single segment and omit thesubscript t ∗ for brevity. We may focus on the case that R ( t ) follows R ad. ( t ) without loss of generality because R n.ad. ( t ) = 1 /R ad. ( t ) and Eq. (11) is antisymmetric un-der the transformation R (cid:55)→ /R . Since we assume R ( t )to be in the vicinity of R ad. ( t ) for t ∈ [ t − , t + ] we studythe linearised equation of motion about R ad. ( t ):˙ R = − iλ ( t ) R − if ( t ) . (18)The general solution from time t = t of this equation is R ( t ) = R ( t ) e Ψ( t ) − Ψ( t ) − i (cid:90) tt dt (cid:48) f ( t (cid:48) ) e Ψ( t ) − Ψ( t (cid:48) ) , (19)where R ( t ) is the initial condition andΨ( t ) = − i (cid:90) tt ∗ λ ( t (cid:48) ) dt (cid:48) . (20)Note that, to first order about t ∗ we have λ ( t ) = λ ( t ∗ ) +˙ λ ( t ∗ )( t − t ∗ )+ O (( t − t ∗ ) ). Since (cid:61) λ ( t ∗ ) = 0 and (cid:61) ˙ λ ( t ∗ ) > (cid:60) Ψ( t ) = (cid:61) ˙ λ ( t ∗ )( t − t ∗ ) + O (( t − t ∗ ) ) is convex. We refer to this property of Ψ below. Integrating theintegral in Eq. (19) by parts N times yields R ( t ) =[ R ( t ) − R ad. ( t )] e Ψ( t ) − Ψ( t ) + R ad. ( t ) + ∆( t ) e Ψ( t ) . (21)Here we have introduced R ad. ( t ) = N − (cid:88) n =0 (cid:18) − iλ ( t ) ddt (cid:19) n R ad. ( t ) , (22)which encapsulates the following of R ad. ( t ): The n = 0term in R ad. ( t ) is simply R ad. ( t ), and the higher orderterms are corrections due to finite variations in λ ( t ) and f ( t ). However, since each term in the sum contains aderivative and therefore scales with n !, there is an optimaltruncation N = N op. beyond which the sum diverges.The precise value of N op. is problem-specific, but for mostpurposes including only the first few terms in the sumEq. (22) is sufficient.The final term in Eq. (21), ∆( t ) e Ψ( t ) , is the remain-ing part of the solution which is not included in the sumEq. (22). It therefore describes the non-trivial part ofthe dynamics that inevitably causes a departure from R ad. ( t ). Since ∆( t ) e Ψ( t ) is the remainder of an optimallytruncated sum it is negligible whenever the solution fol-lows R ad. ( t ). On the other hand, for times t ≈ t + when∆( t ) e Ψ( t ) starts to dominate, R ad. ( t ) is negligible and wemay approximate∆( t ) e Ψ( t ) (cid:39) − ie Ψ( t ) (cid:90) tt dt (cid:48) f ( t (cid:48) ) e − Ψ( t (cid:48) ) . (23)Since (cid:60) Ψ is convex and since Ψ( t ) ∝ ε − ( t ), the inte-grand in Eq. (23) is non-negligible only for times t (cid:48) ≈ t ∗ and the value of ∆ becomes quite independent of t > t ∗ and t < t ∗ . Therefore, under quite general conditions,we can approximate ∆( t ) e Ψ( t ) (cid:39) Θ( t − t ∗ )∆ e Ψ( t ) , whereΘ is the Heaviside step function, and∆ = − i (cid:90) t + t − dtf ( t ) e − Ψ( t ) . (24)The precise values of t − and t + are of little importancein the evaluation of this integral, only that they are farenough from t ∗ that the integrand is negligible at them.We thus arrive at R ( t ) (cid:39) R ad. ( t ) + [ A + Θ( t − t ∗ )∆] e Ψ( t ) , (25)where A = [ R ( t ) −R ad. ( t )] e − Ψ( t ) depends on the initialcondition.From Eq. (25), and the analogous expression for a seg-ment that follows R n.ad. ( t ), we construct our piece-wiseaddition of segments by determining the exit time t + ofa segment from the condition | R ( t + ) | = 1, i.e. when thesolution is ‘half-way between’ R ad. ( t ) and R n.ad. ( t ), andthen using this as the entry time t − for the next segment.Two effects may cause this transition. Firstly, if the solu-tion does not have enough time to approach, e.g., R ad. ( t )sufficiently closely by the critical time, then the finite dif-ference | R ( t ∗ ) − R ad. ( t ∗ ) | will be exponentially amplifiedafter t = t ∗ . This mechanism is responsible, e.g., for theinitial transitions one observes in a single encircling of anEP, where the system is initialised to R (0) = 0 (cid:54)≈ R ad. ( t ).Secondly, however, we see that even for A = 0 a destabili-sation occurs due to a dynamical mechanism representedby ∆ (cid:54) = 0, which yields the time of stability loss t + via | ∆ e Ψ( t + ) | = 1 . (26)The time t + determined in this way is independent oftransients and therefore characterises the longest timethe solution can remain stable after t ∗ . In the quasi-adiabatic limit ε ( t ) → maximal delay time t ∗ + (see Appendix A 5 for moredetails). D. Analytic examples
Here we consider three examples analytically in orderto illustrate our generalised quasi-adiabatic solution: acircular λ ( t ) as in Sec. III B; a linear λ ( t ) correspond-ing to the lowest order Taylor expansion; an elliptical λ ( t ) corresponding to the lowest order Fourier expansion.The first example serves to verify that our generalisedquasi-adiabatic solution recovers the more specific ana-lytic results in Sec. III B. The second and third examplesserve to illustrate the sensitivity of ∆ and therefore t + in Eq. (26) to the global path—stability loss delay is aglobal phenomenon. Note that a circular, elliptical, orlinear λ ( t ) does not precisely correspond to a circular,elliptical, or linear path in parameter space unless we arein, say, the limit r (cid:28)
1. We study particular paths inparameter space numerically in Sec. III E.
Circular λ ( t ) . —From Sec. III B, λ ( t ) = i √ rγe − iπt/T and (27) f ( t ) = iπ T . (28)The adiabatic fixed point with corrections is R ad. ( t ) = − εe iπt/T N − (cid:88) n =0 n ! (cid:16) iεe iπt/T (cid:17) n (29)where ε = π/ √ rγT and the optimal truncation is N ∼ (2 ε ) − . Furthermore, about t ∗ = 3 T / t ) = 12 ε ( ie − iπt/T + 1) and (30)∆ = − πe − ε . (31)Putting these expressions together in Eq. (25) recoversEq. (16), and solving for t + in Eq. (26) yields the delay time t + − t ∗ = Tπ arccos(2 ε log π ) , (32)which in the limit ε → t + − t ∗ = T /
2, inagreement with Sec. III B.
Linear λ ( t ) . —Let us now consider another importantscenario, where the line of instability is crossed in a linearsweep, λ ( t ) = λ (cid:60) + i ˙ λ (cid:61) t and (33) f ( t ) (cid:39) f ( t ∗ ) = const (34)where ˙ λ (cid:61) >
0. In this case we have Ψ( t ) = − iλ (cid:60) t + ˙ λ (cid:61) t about t ∗ = 0 and the discontinuity is∆ = − if ( t ∗ ) (cid:114) π ˙ λ (cid:61) e − λ (cid:60) / ˙ λ (cid:61) . (35)From these expressions we deduce the delay time t + − t ∗ = λ (cid:60) ˙ λ (cid:61) (cid:118)(cid:117)(cid:117)(cid:117)(cid:116) λ (cid:61) λ (cid:60) log (cid:115) ˙ λ (cid:61) π | f ( t ∗ ) | , (36)which, in the quasi-adiabatic limit becomes t + − t ∗ = λ (cid:60) / ˙ λ (cid:61) . One might na¨ıvely hypothesize that Eq. (35)describe more general paths by using λ (cid:60) = (cid:60) λ ( t ∗ ) and λ (cid:61) = (cid:61) ˙ λ ( t ∗ ). However, a comparison with the circularpath above already shows that this would only give ratherpoor quantitative results. Eq. (36) may still serve fora first estimate of the expected delay times in generalscenarios. Elliptical λ ( t ) . —As an interpolation between the twocases above we consider the lowest order Fourier expan-sion of λ ( t ) about t = t ∗ : λ ( t ) = λ (cid:60) cos( πt/T ) + i T ˙ λ (cid:61) π sin( πt/T ) and (37) f ( t ) (cid:39) f ( t ∗ ) = const . (38)With λ (cid:60) = T ˙ λ (cid:61) /π = √ rγ one recovers the circular λ ( t ),and for T → ∞ , but keeping ˙ λ (cid:61) fixed one recovers thelinear sweep of λ ( t ). For this example the discontinuityis ∆ = − iT f ( t ∗ ) e − T λ (cid:61) π I Tπ (cid:115) T ˙ λ (cid:61) π − λ (cid:60) , (39)where I is the first order modified Bessel function of thefirst kind. By taking the appropriate limits— I ( x ) ∼ x (cid:28) x ∈ R + (see, e.g., 9.6.7 in Ref. [48])for the circular λ ( t ) and I ( x ) ∼ e x / √ πx for x (cid:29) x ∈ R + (see, e.g., 9.6.30 and 9.7.1 in Ref. [48]) forthe linear λ ( t )—one recovers either Eq. (31) or Eq. (35).Therefore, Eq. (39) interpolates between the two limit-ing cases above and can be used to accurately calculatedelay times for situations where the encircling path liessomewhere in between. FIG. 4. (Colour online.) Plots of | R ( t ) | (three bottom panels)for three different paths (top row). In every plot of | R ( t ) | thesolid line is our generalised quasi-adiabatic solution, the opensquares denote the numerical solution, the dashed lines denote R ad. ( t ) and R n.ad. ( t ), and the dot-dashed lines denote criticaltimes t ∗ . For the plots of the path, the yellow-filled circlemarks the position of the EP and the dot-dashed line is thecritical line where (cid:61) λ = 0. The parameter settings chosen are:(a) r = 0 . T = 200, and g o.s. = 0 .
2; (b) r = 0 . T = 200, e = 0 .
75, and θ aa. = π/
4; (c) L = 0 . T = 200, and g o.s. =0 .
05 (see the main text for details on the parametrisation).
E. Numerical examples
Let us now demonstrate the validity of our approachnumerically for more general examples depicted in Fig. 4:(a) a displaced circular path, ω ( t ) = r sin φ ( t ) and g ( t ) = γ/ r cos φ ( t ) + g o.s. where φ ( t ) = 2 πt/T and g o.s. is a variable off-set in the coupling;(b) a tilted elliptical path, ω ( t ) = r ( t ) sin φ ( t ) and g ( t ) = γ/ r ( t ) cos φ ( t ) where φ ( t ) = 2 πt/T , r ( t ) = r (1 − e ) / { e cos[ φ ( t ) + θ aa. ] } , e is the ellipticity and θ aa. is the angle of the apoapsis;(c) oscillations along a straight path which crosses thecritical line, ω ( t ) = − L sin φ ( t ) and g ( t ) = γ/ g o.s. where φ ( t ) = 2 πt/T . FIG. 5. A plot of the departure time t + for path (c) in Fig. 4with the same parametrisation, except for g o.s. , which we vary.The solid line is t + as determined by Eq. (26) and the opensquares denote the numerically observed departure time asdetermined via | R ( t ) | = 1. Good agreement is exhibited be-tween the analytic t + s and the numerical t + s. A particularlyinteresting feature is that t + becomes infinite for g o.s. (cid:38) . g o.s. one then expects the system to remain adiabaticfor all times, as observed in Fig. 4(a). For these examples the numerically simulated solution R ( t ) is plotted in Fig. 4 and compared with the gen-eralised quasi-adiabatic solution presented in Sec. III C,with ∆ being evaluated numerically. We see that in allcases the numerical and analytic results match up per-fectly.Let us first consider path (a). For this path there aretwo dynamical bifurcations in every period, as indicated,but the departure time as determined by Eq. (26) is in-finite. Therefore, the solution never has enough time tobe significantly repelled from R ad. ( t ) before R ad. ( t ) be-comes stable once again. As a result, the system neverleaves the neighbourhood of R ad. ( t ), i.e., the solution isadiabatic. In some sense, the increased frequency of dy-namical bifurcations and the long departure time has sta-bilised the adiabatic prediction. In path (b) the oppositeis the case. Here the solution undergoes a non-adiabatictransition every period. The solution looks quite similarto that shown in Fig. 4(a), except that the non-adiabatictransitions occur earlier. One finds in this case that thedelay time is roughly 0 . T , slightly less than T /
2. Fi-nally, in path (c) we have chosen again a path with twodynamical bifurcations per period, but in this case thedeparture time is roughly 0 . T , significantly less than T /
2. Accordingly, we observe two non-adiabatic transi-tions per period.The last case can in fact be made to resemble either ofthe former two cases by tuning g o.s. , i.e., by changing thevalue of λ (cid:60) ( t ∗ ). In Fig. 5 we have plotted the departuretime as a function of g o.s. . For g o.s. < .
05 the quasi-adiabatic condition breaks down. For 0 . < g o.s. (cid:46) . g o.s. = 0 .
12 the departure time becomes in-finite, which implies that for g o.s. (cid:38) .
12 the dynamicsbecomes fully adiabatic and resembles Fig. 4(a).
IV. NOISE
Finally, it is important to address the influence ofnoise, which will be present in any experimental imple-mentation. To do so we simulated the dynamics of c ∓ in the presence of delta-correlated Gaussian noise ξ ( t )with variance (cid:104) ξ ( t ) ξ ( t (cid:48) ) (cid:105) = γN δ ( t − t (cid:48) ). The gray shadedarea in Fig. 3(a) indicates the resulting distribution ofstochastic trajectories of R ( t ) for N = 1 /
10. This re-sembles the case where the initial resonator amplitudeis a factor of 10 above the thermal noise floor. For thefirst encircling period the fixed point R n . ad . ( t ) is still ro-bust, but the delay time t + is significantly reduced. Thisagain demonstrates the extreme sensitivity of Eq. (11)upon initial conditions. However, the dynamics of R isself-correcting, and after the first encircling period it set-tles into robust periodic dynamics much resembling thecase without noise. This surprising observation can beunderstood as follows. Initially, noise causes R to loosestability early. But this means that the total populationof the system is increased, and therefore the effect of theconstant noise background is reduced. V. CONCLUSION
In summary, we have analyzed the quasi-adiabatic evo-lution of non-Hermitian systems near an EP. Our studyshows that various dynamical phenomena associated withthis process can be predicted from the analysis of thenon-adiabatic transition amplitudes, R − ( t ) and R + ( t ).In particular, we identified a characteristic switching pat-tern and stability loss delay. Our analytic predictions forthe delay times and the observed robustness with respectto noise are relevant for first experimental investigationsof these effects and provide the basis for analysing similarphenomena in more complex systems. Acknowledgments.
The authors thank M. V. Berry,E.-M. Graefe, M. Koller, D. O. Krimer, A. Mailybaev,A. Neishtadt, and R. Uzdin for useful input. This workwas supported by the European FP7/ITC Project SIQS(600645), project OPSOQI (316607) of the WWTF, theAustrian Science Fund (FWF) through SFB FOQUSF40, SFB NextLite F49-10, project GePartWave No. I1142-N27, and the START grant Y 591-N16.
Appendix A: Non-standard analysis of relativenon-adiabatic transition amplitues
In this appendix we briefly summarise the motiva-tion for the non-standard analysis of quasi-adiabatic non-Hermitian systems [24, 45, 46] and its application to rel-ative non-adiabatic transition amplitudes [34, 49–52]. Inorder to facilitate a simple but rigorous mathematicaltreatment we augment the notation of the paper by in- troducing a dimensionless time s := tT , (A1)where T − is considered ‘small’, and rewrite the govern-ing equation of motion˙ U = − i M ( s ) U (A2)with U (0) = 1 and where U is the evolution operator andthe dot denotes the derivative with respect to time t asusual.
1. Quasi-adiabaticity
We assume M ( s ) to be diagonalisable with eigenvalues λ i ( s ) for all s . Since M ( s ) is non-Hermitian it does notin general have an orthonormal eigenbasis in the sense ofDirac but rather a biorthogonal eigenbasis: a set of righteigenvectors (cid:126)r i ( s ) defined via M ( s ) (cid:126)r i ( s ) = λ i ( s ) (cid:126)r i ( s ) anda set of left eigenvectors (cid:126)l Ti ( s ) defined via (cid:126)l Ti ( s ) M ( s ) = (cid:126)l Ti ( s ) M ( s ) such that (cid:126)l Ti ( s ) (cid:126)r j ( s ) = δ i,j . Ideal adiabaticdynamics may be defined as that for which the dynamicalcoefficients of the instantaneous eigenvectors decouple.In a parallel transported eigenbasis, i.e., (cid:126)l Ti ( s ) (cid:126)r (cid:48) i ( s ) = 0where the prime denotes the derivative with respect to s ,the adiabatic solution, or adiabatic prediction, is U i,j ( t ) = δ i,j (A3)where we have expanded the evolution operator U thus U ( t ) = (cid:88) i,j U i,j e − iT (cid:82) s ds (cid:48) λ i ( s (cid:48) ) (cid:126)r i ( s ) (cid:126)l Tj (0) . (A4)In the adiabatic solution the interaction between thedynamical coefficients of the instantaneous eigenvectorsdue to the finite variation of these eigenvectors is ignored.The full equation of motion for U expanded as above is˙ U p,q = − i (cid:88) i (cid:54) = p T − ˜ f p,i ( s ) e − iT (cid:82) s ds (cid:48) [ λ i ( s (cid:48) ) − λ p ( s (cid:48) )] U i,p (A5)where we have defined˜ f p,i ( s ) := − i(cid:126)l Tp ( s ) (cid:126)r (cid:48) i ( s ) . (A6)The adiabatic solution ignores ˜ f p,q ( s ) for p (cid:54) = q . Assum-ing the system to be initialised to instantaneous eigen-vector q , first order perturbation theory yields that thesolution for the coefficient x p of instantaneous eigenvec-tor p where p (cid:54) = q is x p ( t ) (cid:39) T − ˜ f p,q ( s ) e − iT (cid:82) s ds (cid:48) [ λ q ( s (cid:48) ) − λ p ( s (cid:48) )] λ q ( s ) − λ p ( s ) . (A7)This expression vanishes linearly with | T − / [ λ q ( s ) − λ p ( s )] but diverges exponentially with T (cid:61) (cid:82) s ds (cid:48) [ λ q ( s (cid:48) ) − λ p ( s (cid:48) )] if (cid:82) s ds (cid:48) λ q ( s (cid:48) ) > (cid:82) s ds (cid:48) λ p ( s (cid:48) ). Second order per-turbation theory contains no more information as regards x p but does reveal that x q ( t ) differs from unity with ananalogous scaling. The traditional quantum adiabaticcondition ε p,q ( t ) := (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) T − ˜ f p,q ( s ) λ q ( s ) − λ p ( s ) (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) (cid:28) U for which (cid:61) (cid:82) s ds (cid:48) λ i ( s (cid:48) ) is greatest, i.e., the least dis-sipative instantaneous eigenvectors. It obviously cannotbe the case that every eigenvector is least dissipative,unless all are degenerate, and it is therefore impossiblethat the adiabatic solution Eq. (A3) hold. In the contextof non-Hermitian systems it therefore seems pertinent tocall Eq. (A8) the quasi -adiabatic condition. So long as weinitialise to the least dissipative instantaneous eigenstateand so long as this eigenstate remains the least dissipa-tive, the quasi-adiabatic condition ensures adiabaticity.But if we initialise to an eigenstate that is not the leastdissipative, or the quality of being least dissipative is ex-changed, then perturbation theory breaks down.
2. Relative non-adiabatic transition amplitudes asa slow-fast system
In the main text we argued that for our two-dimensional case the simplest encompasing dynamicaldescription of adiabaticity is afforded by the relative non-adiabatic transition amplitudes—the ratios of the ele-ments of the evolution operator expressed in a paralleltransported eigenbasis. Let us recall the equation of mo-tion for the relative non-adiabatic transition amplitude R ( t ) as defined in the paper:˙ R = − iλ ( s ) R − iT − ˜ f ( s )(1 + R ) . (A9)Treating λ and ˜ f as dynamical variables themselves, inthe limit T − → R = − iλ ( s ) R − iT − ˜ f ( s )(1 + R ) (A10)where s = T − t and t is the initial time. On the otherhand, we may rewrite the equation of motion using s asthe independent variable, T − R (cid:48) = − iλ ( s ) R − iT − ˜ f ( s )(1 + R ) , (A11)whereupon similarly taking T − → − iλ ( s ) R − iT − ˜ f ( s )(1 + R ) . (A12)The difference between Eqs. (A10) and (A12) is that theformer is over a time-scale of order T − and is hence fast , whilst the latter is over a time-scale of order 1 andis hence slow ; we have a slow-fast system. Furthermore,we notice here that the fast time-scale equation is differ-ential whilst the slow algebraic. This is often taken as the definition of a singularly perturbed system and it meansthat any perturbative approach in T − can only be validfor times of order T − . In order to study the long timebehaviour we must turn to non-standard analysis.
3. Slow manifolds
The solutions of the slow time-scale Eq. (A12) are thefixed points of the fast time-scale Eq. (A10) and as suchare known as instantaneous fixed points. We recall theirapproximate expressions from the main text: R ad. ( s ) (cid:39) − T − ˜ f ( s )2 λ ( s ) and R n.ad. ( s ) (cid:39) − λ ( s ) T − ˜ f ( s ) . (A13)Since these are the fixed points of Eq. (A10) we mayuse Eq. (A10) to perform a stability analysis. One findsthat R ad. ( s ) is stable if and only if (cid:61) λ ( s ) <
0, whilst R n.ad. ( s ) is stable if and only if (cid:61) λ ( s ) >
0, and the pos-sible local phase portraits are stable star, stable spiral,centre, unstable spiral, and unstable star [53]. Evidently,the only possible bifurcation is from a stable spiral to anunstable spiral through a centre. The locus of points R ad. ( s ) over s is called a slow manifold and denoted M ad. = { R ad. ( s ) : s } . Similarly for R n.ad. ( s ).
4. Adiabatic manifolds
Due to finite variations in λ ( s ) and ˜ f ( s ) the slow man-ifolds M ad. and M n.ad. are in fact not locally invariant.Nevertheless, a theorem due to N. Fenichel [54] ensuresthe existence of locally invariant manifolds in a T − -neighbourhood of M ad. and M n.ad. . These locally invari-ant manifolds are called adiabatic manifolds and denoted M ad. and M n.ad. respectively. The adiabatic manifoldsdo not obey a simple equation of motion, but we may finda good approximation by considering the particular inte-gral of the N -times linearised equation of motion. We fo-cus on M ad. for clarity. Following the argument of exam-ple 2.1.10 in Ref. [34] one arrives at M ad. = {R ad. ( s ) : s } where R ad. ( s ) (cid:39) N − (cid:88) n =0 T − n (cid:18) − iλ ( s ) dds (cid:19) n R ad. ( s ) (A14)and N is an optimal truncation with a remainder of order e − C/T − for some C >
0. The expression for R n.ad. ( s ) isanalogous. We describe R ad. ( s ) and R n.ad. ( s ) as attrac-tive or unattractive analogously to R ad. ( s ) and R n.ad. ( s )being stable or unstable respectively.1
5. Stability loss delay
At certain critical times t ∗ , or s ∗ , the stability of theinstantaneous fixed points swaps. For example, R ad. ( s )becomes unstable and R n.ad. ( s ) becomes stable at s = s ∗ such that (cid:61) λ ( s ∗ ) = 0 and (cid:61) λ (cid:48) ( s ∗ ) >
0. One mightna¨ıvely suppose an immediate transition between M ad. and M n.ad. at s = s ∗ , but this is not the case. The bi-furcation is dynamical and the type is degenerate Hopf,which in general exhibits the phenomenon known as sta-bility loss delay : the solution R ( t ) continues to follow,say, R ad. ( s ) for a significant time past its loss of stability.Following the argument in Sec. 2 of Ref. [52] one findsthat away from s = s ∗ the solution has the asymptoticexpansion R ( t ) ∼ Ae ˜Ψ( s ) /T − + R ad. ( s ) (A15)where ˜Ψ( s ) = − i (cid:90) ss ∗ ds (cid:48) λ ( s (cid:48) ) , (A16)whilst at s = s ∗ the solution exhibits the discontinuity∆ = − i (cid:90) s ∗ + s ∗− ds ˜ f ( s ) e − ˜Ψ( s ) /T − (A17)where s ∗− < s ∗ and s ∗ + > s ∗ are the intersections of thelevel curve of (cid:60) Ψ that includes the point z ∗ ∈ C suchthat λ ( z ∗ ) = 0. In order to incorporate this discontinu-ity in the asymptotic expansion of the solution for R ( t )we add the term Θ( s − s ∗ )∆ e ˜Ψ( s ) /T − where Θ is theHeaviside step function. This term is proportional to e ( ˜Ψ( s ) − ˜Ψ( z ∗ )) /T − , and therefore diverges as T − → t such that (cid:60) ˜Ψ( s ) − (cid:60) ˜Ψ( z ∗ ) >
0. Thus, thetimes t ∗− < t ∗ and t ∗ + > t ∗ corresponding to s ∗− and s ∗ + are such that:(i) if the solution enter a neighbourhood of M ad. before t ∗− then it must leave at t ∗ + ;(ii) if the solution enter a neighbourhood of M ad. after t ∗− at t − < t ∗ then it must leave at t + > t ∗ suchthat (cid:60) ˜Ψ( s + ) − (cid:60) ˜Ψ( s − ) = 0;(iii) if the solution leave a neighbourhood of M ad. after t ∗ + then it must have entered at t ∗− .Since the third case is sure to be rare, one typically calls t ∗ + the maximal delay time . Note that, this maximal de-lay time is precisely the quasi-adiabatic limit of the delaytime t + calculated via Eq. (26).For an analytic example of a maximal delay time, letus consider the path analysed in Sec. III B. In this casethe level curves of (cid:60) ˜Ψ are e − π (cid:61) z cos( π (cid:60) z ) (cid:39) const . (A18)Evidently, any s − < s − > − / s + = − s − > FIG. 6. Plot of the theoretical maximal delay time (solidline) and numerically observed departure times from a 5%-neighbourhood of R ad. ( s ) (open squares) for case (c) fromSec. III E. delay time is therefore t ∗ + = T /
2, which agrees with themain text.For a numerical example, let us consider case (c) fromSec. III E. Here we calculate the theoretical maximal de-lay time by numerically finding the complex root of λ ( z ),i.e., z ∗ , and then numerically finding where the levelcurve of (cid:60) ˜Ψ on which z ∗ lies intersects the real axis.The results are plotted in Fig. 6. Also plotted in the fig-ure are the results of a numerical solution where we haveinitialised to a neighbourhood of R ad. ( s ) at s = s ∗− andasked for when the numerical solution departs from thisneighbourhood. The agreement between the theoreticalmaximal delay time and the numerically observed depar-ture time is very good. Furthermore, we see qualitativelythe same results as those shown in Fig. 5; the quantita-tive difference is simply due to the finite time required for R ( t ) to leave a small neighbourhood of M ad. and reach | R ( t ) | = 1. Appendix B: Non-adiabatic transitions as amanifestation of the Stokes phenomenon ofasymptotics
The Stokes phenomenon of asymptotics is that sub-dominant exponentials in the asymptotic expansion ofcertain functions disappear and reappear in different sec-tions of the complex plane with different coefficients. InStokes’ words [55],. . . the inferior term enters as it were intothe mist, is hidden for a little while from view,and comes out with its coefficient changed.M. V. Berry and R. Uzdin [27] uncovered the presence ofthe Stokes phenomenon of asymptotics in the solutionsof specific exactly solvable models of quasi-adiabatic non-Hermitian systems and identify non-adiabatic transitionsin such systems as a manifestation. In this appendix wereview one such example.Let us consider again the parametrisation studied inSec. III B: ω ( t ) = r sin φ ( t ), γ ( t ) = 1, and g ( t ) = 1 / r cos φ ( t ) with r (cid:28) φ ( t ) = ˙ φ = 2 π/T = const.An exact solution for this example is presented by M.V. Berry [32] which we now review. Constructing theanalagous quantity to R in the circular basis (cid:126)r (cid:8) = 1 √ (cid:18) i (cid:19) and (cid:126)r (cid:9) = 1 √ (cid:18) − i (cid:19) , (B1)which we denote p , i.e., (cid:126)x ( t ) = c (cid:8) ( t ) (cid:126)r (cid:8) + c (cid:9) ( t ) (cid:126)r (cid:9) and p ( t ) = c (cid:9) ( t ) /c (cid:8) ( t ), one arrives at the equation of motion˙ p (cid:39) re iφ ( t ) + p , (B2)where we have used r (cid:28)
1. Note that, (cid:126)r (cid:8) and (cid:126)r (cid:9) are noteigenvectors. Introducing the new independent variable ζ = ( i/ ε ) − e iφ ( t ) / where ε = π/ (4 √ rT ) and usingthe ansatz p ( t ) = − ( d/dt ) log f ( ζ ) yields ζ d fdζ + ζ dfdζ + ζ f = 0 . (B3)This is the zeroth order Bessel equation and the solutionfrom time t = t is thus p ( t ) = i φζ C ( ζ ) C ( ζ ) (B4)where C n ( ζ ) = c J J n ( ζ ) + c Y Y n ( ζ ) is a linear combinationof order n Bessel functions of the first and second kindwith the ratio c Y c J = − ip ( t ) J ( ζ ) + ˙ φζ J ( ζ )2 ip ( t ) Y ( ζ ) + ˙ φζ Y ( ζ ) (B5) where p ( t ) is the initial condition for p and ζ =( i/ ε ) − e iφ ( t ) / .In order to translate this result into an expression for R ( t ) we have only to transform from the circular basisEq. (B1) to the original basis and then from that to theparallel transported eigenbasis Eq. (5), (cid:126)r − ( t ) = (cid:18) cos ϑ ( t ) / ϑ ( t ) / (cid:19) and (cid:126)r + ( t ) = (cid:18) − sin ϑ ( t ) / ϑ ( t ) / (cid:19) (B6)where tan ϑ ( t ) = − g ( t ) / [ ω ( t ) + iγ ( t ) / p and R as M¨obiustransformations, it is immediately seen that p ( t ) = e iϑ ( t ) iR ( t )1 − iR ( t ) , (B7)where R ( t ) is the initial condition for R , and R ( t ) = i − e − iϑ ( t ) p ( t )1 + e − iϑ ( t ) p ( t ) . (B8)Let us focus on the asymptotic expansion about t = t ∗ where R ad. ( t ) becomes unstable. Recalling that λ ( t ) (cid:39)√ re iφ ( t ) / and (cid:61) λ ( t ∗ ) = 0 we see that with ζ ∗ =( i/ ε ) − e iφ ( t ∗ ) / we have (cid:60) ζ ∗ = 0. Without loss ofgenerality, we suppose (cid:61) ζ ∗ >
0. Using 9.1.35, 9.1.36,9.2.1, and 9.2.2 in Ref. [48] and assuming that the solu-tion is exponentially close to R ad. ( t ) by t = t ∗ yields R ( t ) ∼ i − e ϑ ( t ) p ad. ( t ) + 2 i Θ( −(cid:60) ζ ) e iζ [1 − e ϑ ( t ) p n.ad. ( t )]1 + e ϑ ( t ) p ad. ( t ) + 2 i Θ( −(cid:60) ζ ) e iζ [1 + e ϑ ( t ) p n.ad. ( t )] (B9)for − π/ < arg ζ < π/ p ad. ( t ) and p n.ad. ( t ) cor-respond to R ad. ( t ) and R n.ad. ( t ) respectively and Θ is theHeaviside step function. The discontinuity in this asymp-totic expansion is precisely the Stokes phenomenon ofasymptotics. For t ∈ ( t ∗ − T / , t ∗ + T /
2) the discontin-uous term is subdominant and R ( t ) ∼ R ad. ( t ), whereasfor t > t ∗ + T / R ( t ) ∼ R n.ad. ( t ).The connection to the expansion Eq. (16) from Sec. III Bis more clearly seen by expanding to first order about R ad. ( t ): R ( t ) ≈ R ad. ( t ) − t − t ∗ ) e − exp( iφ ( t ) / / (2 ε ) . (B10) Comparing this expansion to Eq. (16) we find very goodagreement. The small difference that here the disconti-nuity at t = t ∗ is ∆ = − e − / (2 ε ) whereas in Sec. III D wefound ∆ = − πe − / (2 ε ) is principally due to the asymp-totic expansions employed in calculating Eq. (B9), whichdiffer from those used in calculating Eq. (16). [1] M. Born, Z. Phys. A-Hadron Nucl. , 167 (1927).[2] S. G. Johnson, P. Bienstman, M. A. Skorobogatiy, M.Ibanescu, E. Lidorikis, and J. D. Joannopoulos, Phys. Rev. E , 066608 (2002).[3] E. Farhi, J. Goldstone, S. Gutmann, J. Lapan, A. Lund-gren, and D. Preda, Science , 472 (2001). [4] T. Kato, Perturbation Theory for Linear Operators (Springer-Verlag, Berlin, 1995), 2nd ed.[5] A. P. Seyranian and A. A. Mailybaev,
MultiparameterStability Theory with Mechanical Applications (Worl Sci-entific Publishing, Singapore, 2003).[6] M. V. Berry, Czech. J. Phys. , 1039 (2004).[7] A. P. Seyranian, O. N. Kirillov, and A. A. Mailybaev, J.Phys. A: Math. Gen. , 1723 (2005).[8] I. Rotter, J. Phys. A: Math. Theor. , 153001 (2009).[9] N. Moiseyev, Non-Hermitian Quantum Mechanics (Cam-bridge Univerity Press, Cambridge, 2011).[10] W. D. Heiss, J. Phys. A: Math. Theor. , 444016 (2012).[11] D. C. Brody and E.-M. Graefe, Phys. Rev. Lett. ,230405 (2012).[12] N. Bachelard, C. Garay, J. Arlandis, R. Touzani, and P.Sebbah, arXiv:1407.8220.[13] B. Peng, S¸. K. ¨Ozdemir, S. Rotter, H. Yilmaz, M.Liertzer, F. Monifi, C. M. Bender, F. Nori, and L. Yang,Science , 328 (2014).[14] P. Ambichl, K. G. Makris, L. Ge, Y. Chong, A. D. Stone,and S. Rotter, Phys. Rev. X , 041030 (2013).[15] S. Klaiman, U. G¨unther, and N. Moiseyev, Phys. Rev.Lett. , 080402 (2008).[16] C. E. R¨uter, K. G. Makris, R. El-Ganeiny, D. N.Christodoulides, M. Segev, and D. Kip, Nat. Phys. ,192 (2010).[17] M. Brandstetter, M. Liertzer, C. Deutsch, P. Klang, J.Sch¨oberl, H. E. T¨ureci, G. Strasser, K. Unterrainer, andS. Rotter, Nat. Comm. , 4034 (2014).[18] C. Dembowski, H.-D. Gr¨af, H. L. Harney, A. Heine, W.D. Heiss, H. Rehfeld, and A. Richter, Phys. Rev. Lett. , 787 (2001).[19] T. Gao, E. Estrecho, K. Y. Bliokh, T. C. H. Liew, M. D.Fraser, S. Brodbeck, M. Kamp, C. Schneider, S. H¨ofling,Y. Yamamoto, F. Nori, Y. S. Kivshar, A. Truscott, R.Dall, and E. A. Ostrovskaya, arXiv:1504.00978.[20] O. Latinne, N. J. Kylstra, M. D¨orr, J. Purvis, M. Terau-Dunseath, C. J. Joachain, P. G. Burke, and C. J. Noble,Phys. Rev. Lett. , 46 (1995).[21] R. Lefebvre, O. Atabek, M. ˇSindelka, and N. Moiseyev,Phys. Rev. Lett , 123003 (2009).[22] O. Atabek, R. Lefebvre, M. Lepers, A. Jaouadi, O.Dulieu, and V. Kokoouline, Phys. Rev. Lett. , 173002(2011).[23] A. Kvitsinsky and S. Putterman, J. Math. Phys. , 1403(1991).[24] G. Nenciu and G. Rasche, J. Phys. A: Math. Gen. ,5741 (1992).[25] G. Dridi, S. Gu´erin, H. R. Jauslin, D. Viennot, and G.Jolicard, Phys. Rev. A , 022109 (2010).[26] R. Uzdin, A. Mailybaev, and N. Moiseyev, J. Phys. A:Math. Theor. , 435302 (2011). [27] M. V. Berry and R. Uzdin, J. Phys. A: Math. Theor. ,435303 (2011).[28] I. Gilary, A. A. Mailybaev, and N. Moiseyev, Phys. Rev.A , 010102 (2013).[29] A. Leclerc, G. Jolicard, and J. P. Killingbeck, J. Phys.B: At. Mol. Opt. , 145503 (2013).[30] E.-M. Graefe, A. A. Mailybaev, and N. Moiseyev, Phys.Rev. A , 065302 (2014).[32] M. V. Berry, J. Opt. , 115701 (2011).[33] M. V. Berry, Proc. R. Soc. Lond. A , 7 (1989).[34] N. Berglund and B. Gentz, Noise-induced phenomenain slow-fast dynamical systems: a sample-paths approach (Springer-Verlag, London, 2006).[35] M. Diener, Math. Intell. , 38 (1984).[36] A. I. Neishtadt, Discrete Cont. Dyn.-S , 897 (2009).[37] T. Faust, J. Rieger, M. J. Seitner, J. P. Kotthaus, andE. M. Weig, Nat. Phys. , 485 (2013).[38] H. Okamoto, A. Gourgout, C.-Y. Chang, K. Onomitsu,I. Mahboob, E. Y. Chang, and H. Yamaguchi, Nat. Phys. , 480 (2013).[39] X.-W. Xu, Y. Liu, C.-P. Sun, and Y. Li, arXiv:1402.7222.[40] A. B. Shkarin, N. E. Flowers-Jacobs, S. W. Hoch, A. D.Kashkanova, C. Deutsch, J. Reichel, and J. G. E. Harris,Phys. Rev. Lett. , 013602 (2014).[41] H. Jing, S. K. ¨Ozdemir, X.-Y. L¨u, J. Zhang, L. Yang,and F. Nori, Phys. Rev. Lett. , 053604 (2014).[42] R. Uzdin and M. Moiseyev, Phys. Rev. A , 031804(2012).[43] J. Doppler et al. , (to be published).[44] With the choice Eq. (5) for the instantaneous eigenvec-tors there is no geometric phase [45, 46].[45] A. Leclerc, D. Viennot, and G. Jolicard, J. Phys. A:Math. Theor. , 415201 (2012).[46] S. Ib´a˜nez and J. G. Muga, Phys. Rev. A , 033403(2014).[47] I. S. Gradshteyn and I. M. Ryzhik, Table of Integrals,Series, and Products , A. Jeffrey, and D. Zwillinger (Else-vier, London, 2007), 7th ed.[48] M. Abramowitz and I. A. Stegun,
Handbook of Mathe-matical Functions With Formulas, Graphs, and Mathe-matical Tables , tenth printing (Dover, New York, 1972).[49] A. Bohun, MSc Thesis, Imperial College, London (2011).[50] A. I. Neishtadt, Diff. Equat. , 1385 (1987).[51] A. I. Neishtadt, Diff. Equat. , 171 (1988).[52] F. Diener and M. Diener, “Maximal Delay”, DynamicalBifurcations (Springer, Berlin, 1991).[53] P. Glendinning,
Stability, Instability and Chaos: an In-troduction to the Theory of Nonlinear Differential Equa-tions (Cambridge University Press, 1994).[54] N. Fenichel, J. Differ. Equations , 53 (1979).[55] G. G. Stokes, Acta Math.26