Stability of dynamical quantum phase transitions in quenched topological insulators: From multiband to disordered systems
SStability of dynamical quantum phase transitions in quenched topological insulators:From multiband to disordered systems
Christian B. Mendl
1, 2, ∗ and Jan Carl Budich † Technische Universit¨at Dresden, Institute of Scientific Computing,Zellescher Weg 12-14, 01069 Dresden, Germany Technische Universit¨at M¨unchen, Department of Informatics and Institutefor Advanced Study, Boltzmannstraße 3, 85748 Garching, Germany Institute of Theoretical Physics, Technische Universit¨at Dresden, 01062 Dresden, Germany (Dated: December 26, 2019)Dynamical quantum phase transitions (DQPTs) represent a counterpart in non-equilibrium quan-tum time evolution of thermal phase transitions at equilibrium, where real time becomes analogousto a control parameter such as temperature. In quenched quantum systems, recently the occurrenceof DQPTs has been demonstrated, both with theory and experiment, to be intimately connected tochanges of topological properties. Here, we contribute to broadening the systematic understandingof this relation between topology and DQPTs to multi-orbital and disordered systems. Specifically,we provide a detailed ergodicity analysis to derive criteria for DQPTs in all spatial dimensions, andconstruct basic counter-examples to the occurrence of DQPTs in multi-band topological insulatormodels. As a numerical case study illustrating our results, we report on microscopic simulations ofthe quench dynamics in the Harper-Hofstadter model. Furthermore, going gradually from multi-band to disordered systems, we approach random disorder by increasing the (super) unit cell withinwhich random perturbations are switched on adiabatically. This leads to an intriguing order oflimits problem which we address by extensive numerical calculations on quenched one-dimensionaltopological insulators and superconductors with disorder.
I. INTRODUCTION
Motivated by experimental progress on realizing quan-tum matter far from equilibrium in various physical sys-tems including ultracold atomic gases [1, 2], trapped ions[3–5], nitrogen-vacancy centers in diamond [6] and light-driven condensed matter systems [7, 8], investigating the(coherent) quench dynamics of quantum many-body sys-tems has become a broad frontier of current research [9].A prominent example allowing for a systematic studyof intriguing non-equilibrium phenomena is provided bydynamical quantum phase transitions (DQPTs) [10–24],a counterpart of thermal phase transitions in coherentquantum time evolution, where the role of a control pa-rameter is replaced by real time.The formal analog of a (boundary) partition functionis in the context of DQPTs played by the Loschmidt am-plitude G ( t ) = (cid:104) ψ | e − iHt | ψ (cid:105) = r ( t )e iφ ( t ) , (1)with | ψ (cid:105) denoting the initial state and H denoting theHamiltonian governing the non-equilibrium time evolu-tion, i.e., | ψ (cid:105) is far from being an eigenstate of H . Therole of a free energy density is assumed by the so-calledrate function g ( t ) = − log( |G ( t ) | ) /N , where N is the sizeof the system, i.e., in our present context the number oflattice sites. Further following this formal analogy tothermal systems, DQPTs are then simply hallmarked by ∗ [email protected] † [email protected] non-analytical behavior of g as a function of real time,manifesting in characteristic cusps in g ( t ) or one of itstime-derivatives. These cusps are accompanied by zerosof G ( t ), known in statistical physics as Fisher zeros of thepartition function [25].Taking a closer look at the analytical origin of DQPTs, π -phase slips of the Pancharatnam geometrical phase φ G ( t ) [26, 27] (see Fig. 1a for an illustration) havebeen identified as a generic phenomenon behind the non-analytical behavior of g ( t ) [16]. The phase φ G ( t ) is ob-tained from the total phase φ ( t ) (see Eq. (1)) of the com-plex Loschmidt amplitude by subtracting the dynamicalphase φ G ( t ) = φ ( t ) − φ dyn ( t ) (2)with the dynamical phase φ dyn ( t ) = − (cid:82) t d s (cid:104) ψ ( s ) | H | ψ ( s ) (cid:105) . Now, when G ( t ) goes through aFisher zero, its total phase φ ( t ) generically jumps by π ,as for any zero crossing of a complex-valued function.Since the dynamical phase φ dyn ( t ) is always continuousin time, this jump must occur in the geometricalphase φ G ( t ). For the simple case of a time-dependenttwo-level system – which is immediately relevant for theexperimentally realized two-band models – φ G ( t ) maybe readily visualized using a Bloch sphere representation(see Fig. 1a). In this picture, φ G ( t ) is simply given byhalf of the area bounded by the time evolution trajectorybetween times τ = 0 and τ = t , which is augmented toa closed path by a geodesic connecting its end points.At a Fisher zero, | ψ (cid:105) and | ψ ( t ) (cid:105) then correspond toantipodal points of the Bloch sphere which renders theirgeodesic connection (and with that φ G ( t )) ill-defined.This provides a simple picture of how jumps in φ G ( t ) a r X i v : . [ c ond - m a t . qu a n t - g a s ] D ec (a) geometrical phase φ G ( t )(b) φ Gk ( t ) for the three-band Hofstadter model ( t ) - ( t ) (c) corresponding rate function g ( t ) and closeup of g (cid:48) ( t ) FIG. 1. (a) Interpretation of the Pancharatnam geometricalphase on the Bloch sphere as half of the surface area en-closed by the trajectory up to time t , and the geodesic curveleading back to the initial wavefunction. (b) Time evolutionsnapshots of the geometrical phase for the three-band Hofs-tadter model after a quench. Phase vortices are circled in red,and the marked areas show the admissible region accordingto the criterion (5) (cross-hatched) and the complement ofthe exclusion (6) (dotted). (c) Corresponding rate functionand its derivative. Cusps of g (cid:48) ( t ) hallmark DQPTs, i.e., the(dis-)appearance of Fisher zeros and phase vortex pairs. occur at Fisher zeros hallmarking DQPTs.Among many other intriguing applications (seeRef. [20] for a review), DQPTs have become an impor-tant diagnostic tool for identifying topological insula-tor phases [28, 29] in systems far from equilibrium, ashas been demonstrated in recent experiments on var-ious physical platforms, ranging from ultracold atomicgases [19], over superconducting qubit systems [30], andquantum walks in photonic systems [31, 32], to nanome-chanical settings [33]. The underlying conceptual in-sight is that changes in the topological properties overa quench generically imply the occurrence of DQPTs[15, 16]. Moreover, a one-to-one correspondence distin-guishing such topology-induced DQPTs from accidentalones has been derived by identifying a dynamical topolog-ical order parameter for DQPTs [16]. Shortly thereafter,generalizing the relation between DQPTs and topologi-cal properties, the occurrence of DQPTs in the quenchdynamics in multiband topological insulators has beeninvestigated [34].Our present work is aimed at further generalizing theunderstanding of the interplay between topology andDQPTs. More concretely, the purpose of our analy-sis is twofold: First, we revisit the quench dynamics inmultiband systems, going beyond Ref. [34] by providinga comprehensive ergodicity analysis resulting in criteriafor DQPTs that depend on the spatial dimension of thesystem, and by constructing basic counter-examples tothe occurrence of DQPTs in multi-band topological in-sulator models, where not all individual bands are topo-logically nontrivial (see Sec. II). Furthermore, our resultson multi-band models are supported by numerical simu-lations of the quench dynamics in the Hofstadter model(see Sec. III). Second, we connect the theory of DQPTsin multi-band and disordered systems, by approachingdisorder from an angle of increasing the (super) unit cellwithin which random perturbations are switched on adi-abatically. This leads to an intriguing order of limitsproblem, and to settle the question of whether topology-induced DQPTs generically survive up to a finite disorderstrength, we present extensive numerical simulations onquenches in a disordered one-dimensional (1D) topologi-cal insulator model (see Sec. IV). II. ERGODICITY ANALYSIS FOR MULTIBANDSYSTEMS
We consider free fermions on a (hypercubic) d -dimensional lattice with unit lattice constant and n de-grees of freedom per site. For the quantum quench, thesystem is assumed to be prepared in an insulating stateof a filled lowest Bloch band, forming the ground state ofsome initial Hamiltonian H i , before the system Hamil-tonian is quenched at time t = 0 to a final Hamiltonian H . A. Loschmidt amplitude in multiband latticemodels
Assuming lattice translation invariance, the conser-vation of lattice momentum allows us to factorize theLoschmidt amplitude as G ( t ) = (cid:81) k G k ( t ) with G k ( t ) = (cid:104) ψ k | e − iH ( k ) t | ψ k (cid:105) = r k ( t ) e iφ k ( t ) , (3)where H ( k ) denotes the n × n post quench Bloch Hamil-tonian in reciprocal space and | ψ k (cid:105) is the occupied Blochstate of the initial Hamiltonian.Denoting the eigenvalues and eigenvectors of the post-quench Hamiltonian H f ( k ) by E k,α and | u k,α (cid:105) , respec-tively, Eq. (3) can be written as G k ( t ) = n (cid:88) α =1 |(cid:104) u k,α | ψ k (cid:105)| e − iE k,α t . (4)As mentioned, note that this formula holds for the specialcase of a single filled band. B. General criteria for Fisher zeros
Because of the generalized triangle inequality in thecomplex plane, the occurrence of a Fisher zero at mo-mentum k , i.e., G k ( t ) = 0 for some time t , then requires[34] |(cid:104) u k,α | ψ k (cid:105)| ≤
12 for all α = 1 , . . . , n. (5)This condition affords a simple geometric interpretationwhen thinking of the sum in Eq. (4) as a polygonalchain in the complex plane, the edges of which havelength |(cid:104) u k,α | ψ k (cid:105)| that rotate with independent frequen-cies E k,α : A violation of Eq. (5) then simply means thatone edge dominates in length over all others such thatconcatenating all edges can never lead to a closed poly-gon, independent of their direction.Another relevant criterion for the (non-)occurrence ofFisher zeros at a fixed time t is whether the points { e − iE k,α t } α =1 ,...,n , all lie within a minor arc of the unitcircle; equivalently, whether the convex polygon with ver-tices { e − iE k,α t } α =1 ,...,n (as points in the complex plane)contains the origin. In other words, if there exists a ω ( t ) ∈ R such thatcos( E k,α t − ω ( t )) > α = 1 , . . . , n, (6)then the sum in Eq. (4) cannot be zero.Note that the condition (5) only depends on the ini-tial state and the eigenvectors of H f ( k ), whereas the dy-namical criterion (6) solely depends on the eigenvalues of H f ( k ) and time t . C. Abundance of Fisher zeros
In Ref. [34], it has been shown that quenches from atrivial initial state into a post-quench Hamiltonian, allindividual bands of which have non-zero Chern number,there must be a momentum for which Eq. (5) is satis-fied. Basic ergodicity arguments then imply that G k ( t )must come arbitrarily close to zero at some finite time t . However, these important insights do not yet pro-vide a sufficient condition for the actual occurrence of aFisher-zero, i.e., an exact zero crossing of G k ( t ) at anyfinite time. In the following, we fill this gap by perform-ing an additional dimensional analysis, revealing also thegeneric dependence of the abundance of Fisher zeros onthe spatial dimension d . We note that zeros of the parti-tion function in the complex plane have been studied inthe context of phase transitions for more than 50 years,including the analysis of the dimensional dependence ofcritical exponents [25, 35–40].We start by observing that Eq. (5) for n > admissible region of spatialdimension d , i.e., in a whole neighborhood in momentumspace. Therefore, as a subset of the ( d + 1)-dimensionalmomentum-time space (where momentum space is con-strained to the admissible region), the dimension of themanifold of Fisher zeros G k ( t ) = 0 is generically givenby ( d + 1) − d −
1, since both the real and imag-inary parts of G k ( t ) have to be tuned to zero. Thisdimensional counting is independent of n for n > − iE k, t , . . . , e − iE k,n t ) is ergodic on the n -dimensional torus as long as the energies are rationallyindependent. As a consequence, in a one-dimensionalsystem ( d = 1), the Fisher zeros are expected to occur atisolated points in time-momentum space, while for d = 2,the set of Fisher zeros are curves in the three-dimensionalmomentum-time space, in agreement with microscopicsimulations on the quench-dynamics of two-band modelsin d = 2 [41].We now elaborate on the somewhat exceptional but ex-perimentally highly relevant case n = 2. There, Eq. (5)implies |(cid:104) u k, | ψ k (cid:105)| = |(cid:104) u k, | ψ k (cid:105)| = 1 / d − d -dimensional neighbor-hood found for n >
2. However, this reduction in di-mension of the set of admissible momenta for n = 2is exactly compensated by the fact that then G k ( t ) =e − it ( E k, + E k, ) / cos( t ( E k, − E k, ) /
2) in the admissibleregion which requires only tuning of a single real con-dition (the argument of the cos-function) in order toachieve zeros. Hence, Fisher zeros are now guaranteedto occur at all admissible momenta, namely at the times t k,l = (2 l + 1) π/ ( E k, − E k, ) , l = 1 , , . . . such that theyafter all still form a ( d − n > D. Avoided DQPTs in quenched Chern insulators
Quenches from trivial states to Chern insulator Hamil-tonians imply DQPTs, at least when assuming that allindividual bands of the post quench Hamiltonian havenon-zero Chern number [34]. To demonstrate that thisquite strong assumption is indeed necessary, we constructa basic counter-example, where the post quench Hamilto-nian is in a Chern insulator regime, but where no Fisherzeros or DQPTs occur as not all individual bands havenon-vanishing Chern number. To this end, consider asystem with three bands, where we quench from an ini-tial Hamiltonian with only topologically trivial bands toa Chern insulator which has Chern numbers (1 , , − > / √ III. QUENCHED HOFSTADTER MODEL
In this section, we practically verify our general ergod-icity analysis by time-dependent simulations of DQPTsin multi-band systems. For concreteness, we considerthe q -band (magnetic flux 2 π/q per unit cell) Hofstadtermodel [42, 43] defined in the Landau gauge by the mo-mentum space (Bloch) Hamiltonian H ( k ) = k x ) 1 e − ik y k x − πq ) 1. . .e ik y k x − π ( q − q ) . (7)Because of conservation of lattice momentum, theLoschmidt amplitude factorizes (see Eq. (3)). We con-sider two scenarios for the topologically trivial initialstate | ψ k (cid:105) : (i) the initial state occupies the first orbital(in the basis of Eq. (7)) for each k , i.e., | ψ k (cid:105) = | e (cid:105) , and(ii) | ψ k (cid:105) is equal to a fixed complex-random state (inde-pendent of k ). The second scenario will exemplify theabsence of symmetries beyond lattice momentum conser-vation.Fig. 1b shows snapshots of φ Gk ( t ) at several points intime, for scenario (i) and q = 3. Fisher zeros, at which φ Gk ( t ) is ill-defined, appear as phase vortices at isolated k -points (circled in red) which contain the whole range ofphases, [ − π, π ], in any (arbitrarily small) neighborhood.This is in line with the dimension analysis in Sec. II C:the Fisher zeros should describe a d − t . Indeed the phase vortices stay inside bothregions, as required. Note that the stripe-shaped patternof the static admissible region (for the present model pa-rameters) implies that phase vortex–antivortex pairs areconstrained to remain within a single stripe. At t = 2 . g ( t )and a closeup of its derivative. Because of the factoriza-tion G ( t ) = (cid:81) k G k ( t ), the rate function equals g ( t ) = − N log (cid:0) |G ( t ) | (cid:1) = − | BZ | (cid:90) BZ d k log (cid:0) |G k ( t ) | (cid:1) (8)with BZ = [ − π, π ] being the Brillouin zone of thepresent model. The (weak) log-singularity of the in-tegrand at Fisher zeros leads to a cusp in the deriva-tive g (cid:48) ( t ) at their (dis-)appearance, as visible in Fig. 1c.Specifically, Fisher zeros occur for the first time around t = 1 . t = 2 . φ Gk ( t ) with respect to lattice momentum,first note that H ( k x , − k y ) = H ( k ) T according to Eq. (7),such that for real-valued | ψ k (cid:105) , G ( k x , − k y ) ( t ) = G k ( t ). Inparticular, this mirror symmetry holds in the first sce-nario. Moreover, central inversion ( k → − k ) can be ex-pressed as unitary transformation: Let P q be the q × q permutation matrix which sends the j -th entry of a vec-tor (counting from zero) to − j mod q ( j = 0 , . . . , q − U ( k y ) = P q · e − ik y . (9)Then U ( k y ) † H ( k ) U ( k y ) = H ( − k ) . (10)It follows that G − k ( t ) = (cid:104) U ( k y ) ψ − k | e − iH ( k ) t | U ( k y ) ψ − k (cid:105) . (11)Since U ( k y ) | e (cid:105) = e − ik y | e (cid:105) and since the phase factore − ik y cancels in G − k ( t ), this explains the inversion sym-metry apparent in Fig. 1b.In contrast, for the second scenario (ii) of a com-plex random initial state, our analysis does not predictany momentum symmetry. Fig. 2 shows the geometrical (a) φ Gk ( t ) for the Hofstadter model with random initial state ( t ) - ( t ) (b) corresponding rate function g ( t ) and closeup of g (cid:48) ( t ) FIG. 2. Pancharatnam geometrical phase and rate functionfor the three-band Hofstadter model as in Fig. 1, but for a( k -independent) initial state with complex random entries.The lack of momentum symmetry of the geometrical phase isexpected (see main text). Note that the dynamical criterionfor Fisher zeros (dotted areas) only depends on the spectrumof H ( k ) and thus agrees with Fig. 1b. phase and rate function for the second scenario, and in-deed momentum symmetry is now absent. Nevertheless,the dynamical exclusion criterion in Eq. (6) only dependson the spectrum of H ( k ) and time, and thus agrees forboth scenarios. In particular, it disallows any Fisher ze-ros at t = 2 .
5, as in the first scenario.
IV. DISORDERED SYSTEMSA. General framework
We now gradually extend the framework of Pancharat-nam geometric phase vortices leading to DQPTs frommulti-orbital to disordered systems. To this end we con-sider systems that are still periodic, but with respectto a super cell containing (cid:96) (cid:29) µ c † j c j to µ j c † j c j with µ j = ¯ µ + ∆ µ j and ∆ µ j the perturbation. For sufficientlylarge (cid:96) , the system resembles a disordered system (with-out any periodicity), as the relevant physical propertiesare expected to be negligibly changed when matchingdistant coefficients, i.e., ∆ µ j + (cid:96) = ∆ µ j . The momen-tum representation of the Hamiltonian is now based on asupercell of size (cid:96) . For example, an unperturbed Hamil-tonian in Bogoliubov-de Gennes form H = 12 π (cid:90) T d k (cid:16) ˆ c † k ˆ c − k (cid:17) (cid:16) (cid:126)d ( k ) · (cid:126)σ (cid:17) (cid:18) ˆ c k ˆ c †− k (cid:19) (12)(with (cid:126)σ the vector of Pauli matrices) is changed to H (cid:96) = 12 π (cid:90) T d k ( ˆ χ (cid:96)k ) † h (cid:96) ( k ) ˆ χ (cid:96)k (13)with ˆ χ (cid:96)k = (cid:16) ˆ c k, ˆ c †− k, · · · ˆ c k,(cid:96) − ˆ c †− k,(cid:96) − (cid:17) T (14)and h (cid:96) ( k ) being a 2 (cid:96) × (cid:96) matrix depending on the dis-order realization. The index α in ˆ c k,α appearing in (14)may be interpreted as orbital index.The Loschmidt amplitude defined in (3) becomes inthe supercell representation G (cid:96)k ( t ) = det (cid:16)(cid:10) ψ k,j | e − ih (cid:96) ( k ) t | ψ k,j (cid:48) (cid:11)(cid:17) (cid:96)j,j (cid:48) =1 (15)with the orthonormal ψ k,j , j = 1 , . . . , (cid:96) defining the ini-tial state as Slater determinant | ψ k, · · · ψ k,(cid:96) (cid:105) of occupiedmodes. We denote the complex phase of G (cid:96)k ( t ) by φ (cid:96)k ( t ).Note that one recovers the special case of zero noise as G (cid:96)k ( t ) = (cid:96) (cid:89) j =1 (cid:10) ψ k,j | e − ih (cid:96) ( k ) t | ψ k,j (cid:11) (16)since the matrix in (15) can then be canonically diago-nalized due to translation invariance. In particular, thecorresponding phase φ (cid:96)k ( t ) is then given by the following (cid:96) -fold superposition of phases: φ (cid:96)k ( t ) = (cid:96) (cid:88) j =1 φ k,j ( t ) mod 2 π. (17)The dynamical phase reads in the supercell representa-tion φ dyn ,(cid:96)k ( t ) = − t (cid:96) (cid:88) j =1 (cid:104) ψ k,j | h (cid:96) ( k ) | ψ k,j (cid:105) mod 2 π, (18)and analogously φ G,(cid:96)k ( t ) = φ (cid:96)k ( t ) − φ dyn ,(cid:96)k ( t ). FIG. 3. Pancharatnam geometrical phase for the disordered Kitaev chain with period (cid:96) and increasing disorder strength ∆ µ max .Each row corresponds to a fixed disorder strength (starting from zero disorder in the top row), and each column to a fixedsupercell size (cid:96) . The dashed vertical lines mark the critical momentum k c of the ordered system (∆ µ max = 0). Since G (cid:96)k ( t ) is a real-analytic function of the noisecoefficients (such as ∆ µ j in the example), the non-analytic points of the Pancharatnam geometrical phase(i.e., Fisher zeros of the Loschmidt amplitude) cannotinstantaneously disappear when continuously increasingthe noise strength; instead, the non-analytic points willcontinuously move in the k - t -plane, potentially annihi-lating or being created in pairs. B. Disordered Kitaev chain
As a specific example, we investigate the Kitaev chain[44, 45] described by the Hamiltonian H = (cid:88) j ∈ Z (cid:104) − t (cid:16) c † j c j +1 + h.c. (cid:17) + µ (cid:16) c † j c j − (cid:17) + (cid:16) ∆ c j c j +1 + h.c. (cid:17) (cid:105) (19)where t is the hopping amplitude, µ is the chemical po-tential and ∆ is the superconducting gap.Switching to the Bogoliubov-de Gennes momentumrepresentation of the Hamiltonian, H ( k ) = (cid:126)d ( k ) · (cid:126)τ (20)with (cid:126)d ( k ) = (cid:0) , ∆ sin( k ) , µ − t cos( k ) (cid:1) and (cid:126)τ the Nambupseudospin, one obtains the Pancharatnam geometricalphase defined in (2), which allows to identify singularpoints of the Loschmidt amplitude [16]. We now employthe supercell representation to investigate the effects ofdisorder (see also Appendix A for technical details): Forsimplicity, we solely let the chemical potential in (19) besite-dependent, i.e., µ j = ¯ µ + ∆ µ j with independent andidentically distributed random variables ∆ µ j chosen fromsome interval [ − ∆ µ max , ∆ µ max ] (uniformly distributed);we retain periodicity with period (cid:96) ∈ N , i.e., µ j + (cid:96) = µ j for all j ∈ Z . The specific parameters for the followingare t = 1, ¯ µ = 6 and ∆ = 1. We checked, however,that different disorder scenarios, such as adding noiseto the hopping amplitudes or superconducting gap pa-rameters instead of the potentials, lead to qualitativelysimilar findings regarding the physics of DQPTs. In par-ticular, the cusps of the rate function (see below) remainintact. This holds even though these scenarios differ re-garding their effectiveness in localizing the eigenstates ofthe Hamiltonian.Fig. 3 shows the Pancharatnam geometrical phase fora fixed noise realization but increasing noise strength,and various supercell sizes (cid:96) . According to Eq. (17), thesupercell representation effectively folds back the phasealong the momentum direction. Accordingly, the geomet-rical phase assumes a stripe-like pattern with increasing (cid:96) , i.e., it varies less as a function of momentum.Using the supercell representation, the rate functionfor the present model reads g (cid:96) ( t ) = − π (cid:90) π d k log (cid:2)(cid:12)(cid:12) G (cid:96)k ( t ) (cid:12)(cid:12)(cid:3) /(cid:96). (21)Thus the zeros of G (cid:96)k ( t ) result in (weak) log-singularitiesof the integrand and corresponding cusps of g (cid:96) ( t ).Fig. 4 visualizes g (cid:96) ( t ) for the disordered Kitaev chainand random disorder realizations, illustrating (a) the ef-fect of increasing disorder at fixed supercell size and (b)increasing supercell size at fixed disorder strength. Thecritical time t c of the first Fisher zero for the case withoutdisorder has been obtained semi-analytically [16]. Oneobserves in Fig. 4a that the rate function is continuouslydeformed with increasing noise strength, and while thetime points of the cusps (i.e., Fisher zeros) shift, thecusps do not instantaneously disappear (see also the timederivative around t c on the right). This is expected dueto the real-analytic dependence of the Loschmidt ampli-tude on the noise coefficients, as detailed above. Visually,the perseverance of the cusps can be understood basedon the geometrical phase in Fig. 3. Namely, the cuspscorrespond precisely to the phase vortices, and thus the Δμ = Δμ = / Δμ = Δμ = / t c ( t ) / t c - - g ′ ( t ) (a) g (cid:96) ( t ) and dd t g (cid:96) ( t ) around t c for (cid:96) = 4 ℓ = ℓ = ℓ = ℓ = / t c ( t ) / t c - - g ′ ( t ) (b) g (cid:96) ( t ) and dd t g (cid:96) ( t ) around t c for ∆ µ max = FIG. 4. Rate function and its derivative around the first criti-cal time point t c , for the disordered Kitaev chain with randomdisorder realizations in the supercell representation. (dis-)appearance of cusps and vortex–antivortex pairs atmomenta k and − k with increasing disorder strength isequivalent. This does not happen instantaneously whenturning on disorder at finite (cid:96) , since the vortex positionsdepend continuously on the disorder strength and have afinite distance in momentum at zero disorder.However, in the limit (cid:96) → ∞ the size of the effec-tive Brillouin zone associated with the supercell shrinksto zero, leading to a non-trivial order of limits problemfor the stability of DQPTs against disorder. To set-tle this issue, we performed extensive numerical simu-lations on systems with finite disorder strength and large (cid:96) . Our results, summarized in Fig. 4, give strong numer-ical evidence that the non-analyticities in the rate func-tion hallmarking DQPTs persist up to significant disor-der strength even in the large (cid:96) limit, i.e., when approach-ing the disordered case without residual translational in-variance.Having investigated instances of disorder realizationsso far, we will now analyze averaging effects as the su-percell size increases. Fig. 5 shows the histogram andcorresponding variance of the rate function g (cid:96) ( t ) eval-uated at time point t = 1, for various supercell sizes (cid:96) . The observed ∼ /(cid:96) scaling of the variance is likelydue to the disorder contributions from individual latticessites being (almost) independent, analogous to the sumof independent random variables in the central limit the-orem. This becomes plausible when assuming that (16)holds approximately for weak disorder, and by inserting ( t = ) ℓ = ℓ = ℓ = ℓ = (a) histogram ● ● ● ● ● ● ● ● ~ ℓ - ℓ - - - Var ( g ( t = )) (b) variance FIG. 5. Histogram and corresponding variance for randomdisorder realizations (∆ µ max = ) of the rate function g (cid:96) ( t )evaluated at t = 1. The variance exhibits a ∼ /(cid:96) scaling,with (cid:96) the supercell size. (16) into (21) we get: g (cid:96) ( t ) ≈ − (cid:96) (cid:96) (cid:88) j =1 π (cid:90) π d k log (cid:104)(cid:12)(cid:12)(cid:10) ψ k,j | e − ih (cid:96) ( k ) t | ψ k,j (cid:11)(cid:12)(cid:12)(cid:105) . (22)Now a lattice to momentum transformation applied to h (cid:96) ( k ) may be understood as an orthogonal transforma-tion of the random disorder coefficients, and if these aremultivariate normal distributed, the transformed coeffi-cients will remain independent. V. CONCLUDING DISCUSSION
We investigated the stability and topological proper-ties of dynamical quantum phase transitions going be-yond the minimal setting of lattice translation invari-ant two-band models in two somewhat related directions.First, building up on recent results [34] on the occurrenceof DQPTs in multi-orbital systems, we demonstrated howthe phenomenology of DQPTs depends on the spatialdimension of the system by means of a more in depthergodicity analysis of the Loschmidt amplitude. We em-phasize that our analysis (and Ref. [34]) was based onthe assumption of a single filled band. Hence, the deriva-tion of strict criteria for the occurrence of DQPTs inmulti-band systems with more than one occupied bandsremains an interesting subject of future research.Second, we considered random potential fluctuationswithin a (super) unit cell of increasing size as a routetowards understanding the stability of DQPTs in disor-dered systems. This approach yielded clear analyticalinsights supporting for the considered settings the stabil- ity of DQPTs for finite unit cells with random potential.However, a non-trivial order of limits problem renders ananalytical proof for the truly disordered case of an infi-nite spatial period of the random potential elusive. To fillthis gap, at least for the considered model systems, wepresented numerical simulations for systems with largeunit cells, thus corroborating the existence of DQPTs ashallmarked by non-analyticities of the rate function upto significant disorder strength.The numerical simulations presented in this work en-courage accompanying theoretical investigations: Specif-ically, a promising direction could be a perturbationanalysis (with respect to disorder strength) applied toEq. (15), which should result in (16) as a lowest orderterm. Also, the question of whether and to what extentdisorder contributions to the rate function can indeed betreated as independent (as conjectured in Sec. IV B) maybe settled in future work.We close by briefly discussing the relation of ourpresent analysis to recent other studies on the combi-nation of disorder and DQPTs. In Ref. [46], the inter-play between quasi-periodic potentials and DQPTs hasbeen investigated, demonstrating the existence of Fisherzeros in certain limits of quasi-disorder, and identifyingthe value of the Loschmidt echo as a marker for local-ization. Shortly after, in Ref. [47], DQPTs have beenexemplified to serve as a tool for diagnosing Andersonlocalization transitions in certain disordered 1D and 3Dmodels. Very recently, the effect of disorder on DQPTsin extended toric code models has been analyzed [48].Approaching the fate of DQPTs in disordered systemsby following vortices in the geometric phase in systemswith a growing disordered super-cell, however, is uniqueto our present work.
ACKNOWLEDGMENTS
We acknowledge helpful discussions with Markus Heyl.J.C.B. acknowledges financial support from the GermanResearch Foundation (DFG) through the CollaborativeResearch Centre SFB 1143 (Project No. 247310070) andthe W¨urzburg-Dresden Cluster of Excellence on Com-plexity and Topology in Quantum Matter – ct.qmat(EXC 2147, Project No. 39085490).
Appendix A: Generalized Kitaev chain with periodicsupercell structure
We consider the Kitaev chain as in (19), generalized tosite-dependent coefficients, i.e., H = (cid:88) j ∈ Z (cid:104) − t j (cid:16) c † j c j +1 + h.c. (cid:17) + µ j (cid:16) c † j c j − (cid:17) + (cid:16) ∆ j c j c j +1 + h.c. (cid:17) (cid:105) . (A1)The Hamiltonian may formally be represented inBogoliubov-de Gennes form as H = (cid:16) · · · c † j c j c † j +1 c j +1 · · · (cid:17) × . . . . . . . . . B † j − A j B j B † j A j +1 B j +1 . . . . . . . . . ... c j c † j c j +1 c † j +1 ... with 2 × A j = 12 (cid:18) µ j − µ j (cid:19) and B j = 12 (cid:18) − t j − ∆ ∗ j ∆ j t j (cid:19) . (A2)In the following, we assume periodicity with period (cid:96) ∈ N , i.e., t j + (cid:96) = t j , µ j + (cid:96) = µ j and ∆ j + (cid:96) = ∆ j for all j ∈ Z . Thus we may subsume the creation and annihilationoperators in a spinor χ (cid:96)n = (cid:16) c (cid:96)n c † (cid:96)n · · · c (cid:96)n + (cid:96) − c † (cid:96)n + (cid:96) − (cid:17) T (A3)and represent the Hamiltonian as H = (cid:88) n ∈ Z (cid:2) ( χ (cid:96)n ) † h (cid:96) local χ (cid:96)n + (cid:0) ( χ (cid:96)n ) † h (cid:96) hop χ (cid:96)n +1 + h.c. (cid:1)(cid:3) (A4)with h (cid:96) local = A B B † A B . . . . . . . . . B † (cid:96) − A (cid:96) − (A5)and h (cid:96) hop = B (cid:96) − . (A6) To arrive at a momentum representation of the Hamil-tonian, we use Fourier transformation χ (cid:96)n = 12 π (cid:90) T d k e ikn ˆ χ (cid:96)k (A7)with ˆ χ (cid:96)k = (cid:16) ˆ c k, ˆ c †− k, · · · ˆ c k,(cid:96) − ˆ c †− k,(cid:96) − (cid:17) T . (A8)Here the index α in ˆ c k,α may be interpreted as orbitalindex. The first Brillouin zone is equal to the interval T = [ − π, π ] with periodic boundary conditions. Inserting(A7) into (A4) yields H = 12 π (cid:90) T d k ( ˆ χ (cid:96)k ) † (cid:2) h (cid:96) local + (cid:0) e ik h (cid:96) hop + h.c. (cid:1)(cid:3) ˆ χ (cid:96)k . (A9)Note that the conventional momentum representationof the Kitaev chain is recovered for (cid:96) = 1: in this case, h = A and h = B , such that (for real-valued∆ ) H = 12 π (cid:90) T d k (cid:126)d ( k ) · (cid:126)τ (A10)with (cid:126)d ( k ) = (cid:0) , ∆ sin( k ) , µ − t cos( k ) (cid:1) and (cid:126)τ being theNambu pseudospin.From a slightly different perspective, for the specialcase A = · · · = A (cid:96) − and B = · · · = B (cid:96) − we mayagain use Fourier transformation applied to the orbitals:ˆ χ (cid:96)k,α = (cid:32) ˆ c k,α ˆ c †− k,α (cid:33) = 1 (cid:96) π ( (cid:96) − (cid:88) q =0 e iα ( k + q ) /(cid:96) (cid:32) ˆ c ( k + q ) /(cid:96) ˆ c †− ( k + q ) /(cid:96) (cid:33) . (A11)If this is inserted into (A9) it results in H = 12 π (cid:90) T d k (cid:96) π ( (cid:96) − (cid:88) q =0 (cid:16) ˆ c † ( k + q ) /(cid:96) ˆ c − ( k + q ) /(cid:96) (cid:17) × (cid:104) A + (cid:0) e i ( k + q ) /(cid:96) B + h.c. (cid:1)(cid:105) (cid:32) ˆ c ( k + q ) /(cid:96) ˆ c †− ( k + q ) /(cid:96) (cid:33) , (A12)and with the substitution p = ( k + q ) /(cid:96) ∈ T , one arrivesat H = 12 π (cid:90) T d p (cid:16) ˆ c † p ˆ c − p (cid:17) (cid:2) A + (cid:0) e ip B + h.c. (cid:1)(cid:3) (cid:18) ˆ c p ˆ c †− p (cid:19) . (A13)Thus we have again recovered (A10) (for real-valued ∆),as expected.0 REFERENCES [1] I. Bloch, J. Dalibard, and W. Zwerger, Many-bodyphysics with ultracold gases,
Rev. Mod. Phys.
80, 885–964 (2008).[2] N. Goldman, J. C. Budich, and P. Zoller, Topologicalquantum matter with ultracold gases in optical lattices,
Nat. Phys.
12, 639–645 (2016).[3] R. Blatt and C. F. Roos, Quantum simulations withtrapped ions,
Nat. Phys.
8, 277–284 (2012).[4] P. Jurcevic, B. P. Lanyon, P. Hauke, C. Hempel, P. Zoller,R. Blatt, and C. F. Roos, Quasiparticle engineering andentanglement propagation in a quantum many-body sys-tem,
Nature
Nature
Phys. Rev. B
Rev. Mod. Phys.
82, 1489–1537 (2010).[8] T. Byrnes, N. Y. Kim, and Y. Yamamoto, Exciton-polariton condensates,
Nat. Phys.
10, 803–813 (2014).[9] A. Polkovnikov, K. Sengupta, A. Silva, and M. Vengalat-tore, Colloquium: Nonequilibrium dynamics of closed in-teracting quantum systems,
Rev. Mod. Phys.
83, 863–883(2011).[10] M. Heyl, A. Polkovnikov, and S. Kehrein, Dynamicalquantum phase transitions in the transverse-field Isingmodel,
Phys. Rev. Lett.
Phys. Rev.B
87, 195104 (2013).[12] E. Canovi, P. Werner, and M. Eckstein, First-order dy-namical phase transitions,
Phys. Rev. Lett.
Phys. Rev. B
92, 104306 (2015).[14] M. Heyl, Scaling and universality at dynamical quantumphase transitions,
Phys. Rev. Lett.
Phys. Rev. B
91, 155127(2015).[16] J. C. Budich and M. Heyl, Dynamical topological orderparameters far from equilibrium,
Phys. Rev. B
93, 085416(2016).[17] A. A. Zvyagin, Dynamical quantum phase transitions(Review Article),
Low Temperature Physics
42, 971–994(2016).[18] P. Jurcevic, H. Shen, P. Hauke, C. Maier, T. Brydges,C. Hempel, B. P. Lanyon, M. Heyl, R. Blatt, and C. F.Roos, Direct observation of dynamical quantum phasetransitions in an interacting many-body system,
Phys.Rev. Lett.
NaturePhysics
14, 265–268 (2018).[20] M. Heyl, Dynamical quantum phase transitions: a re-view,
Rep. Prog. Phys.
81, 054001 (2018).[21] B. ˇZunkoviˇc, M. Heyl, M. Knap, and A. Silva, Dynamicalquantum phase transitions in spin chains with long-rangeinteractions: Merging different concepts of nonequilib-rium criticality,
Phys. Rev. Lett.
Phys. Rev. B
97, 060304(2018).[23] X. Qiu, T.-S. Deng, G.-C. Guo, and W. Yi, Dynamicaltopological invariants and reduced rate functions for dy-namical quantum phase transitions in two dimensions,
Phys. Rev. A
98, 021601 (2018).[24] N. Sedlmayr, P. Jaeger, M. Maiti, and J. Sirker, Bulk-boundary correspondence for dynamical phase transi-tions in one-dimensional topological insulators and su-perconductors,
Phys. Rev. B
97, 064304 (2018).[25] M. E. Fisher, The theory of equilibrium critical phenom-ena,
Rep. Prog. Phys.
30, 615–730 (1967).[26] S. Pancharatnam, Generalized theory of interference, andits applications,
Proceedings of the Indian Academy ofSciences - Section A
44, 247–262 (1956).[27] J. Samuel and R. Bhandari, General setting for Berry’sphase,
Phys. Rev. Lett.
60, 2339–2342 (1988).[28] M. Z. Hasan and C. L. Kane, Colloquium: Topologicalinsulators,
Rev. Mod. Phys.
82, 3045–3067 (2010).[29] X.-L. Qi and S.-C. Zhang, Topological insulators and su-perconductors,
Rev. Mod. Phys.
83, 1057–1110 (2011).[30] X.-Y. Guo, C. Yang, Y. Zeng, Y. Peng, H.-K. Li, H. Deng,Y.-R. Jin, S. Chen, D. Zheng, and H. Fan, Observationof a dynamical quantum phase transition by a supercon-ducting qubit simulation,
Phys. Rev. Applied
11, 044080(2019).[31] K. Wang, X. Qiu, L. Xiao, X. Zhan, Z. Bian, W. Yi, andP. Xue, Simulating dynamic quantum phase transitionsin photonic quantum walks,
Phys. Rev. Lett. arXiv e-prints arXiv:1808.03930 (2018).[33] T. Tian, Y. Ke, L. Zhang, S. Lin, Z. Shi, P. Huang,C. Lee, and J. Du, Observation of dynamical phase tran-sitions in a topological nanomechanical system,
Phys.Rev. B
Phys. Rev. Lett.
Prog. Theor. Phys.
38, 72 (1967).[36] R. Abe, Singularity of specific heat in the second orderphase transition,
Prog. Theor. Phys.
38, 322 (1967).[37] M. Suzuki, Note on the singularity of specific heat in thesecond order phase transition,
Prog. Theor. Phys.
Z. Phys.
Z. Phys.
449 (1969).[40] W. van Saarloos and D. A Kurtze, Location of zerosin the complex temperature plane: absence of Lee-Yangtheorem,
J. Phys. A: Math. Gen.
17, 1301 (1984).[41] M. Heyl and J. C. Budich, Dynamical topological quan-tum phase transitions for mixed states,
Phys. Rev. B
Phys. Rev. B
14, 2239–2249 (1976).[43] M. Aidelsburger, M. Atala, M. Lohse, J. T. Barreiro,B. Paredes, and I. Bloch, Realization of the HofstadterHamiltonian with ultracold atoms in optical lattices,
Phys. Rev. Lett.
Phys. Usp.
44, 131 (2001).[45] C.W.J. Beenakker, Search for Majorana fermions in su-perconductors,
Annu. Rev. Condens. Matter Phys.
Phys. Rev.B
95, 184201 (2017).[47] H. Yin, S. Chen, X. Gao, and P. Wang, Zeros ofLoschmidt echo in the presence of Anderson localization,
Phys. Rev. A
97, 033624 (2018).[48] V. Srivastav, U. Bhattacharya, and A. Dutta, Dynamicalquantum phase transitions in extended toric-code mod-els,