Quench dynamics and parity blocking in Majorana wires

We theoretically explore quench dynamics in a finite-sized topological fermionic p-wave superconducting wire with the goal of demonstrating that topological order can have marked effects on such non-equilibrium dynamics. In the case studied here, topological order is reflected in the presence of two (nearly) isolated Majorana fermionic end bound modes together forming an electronic state that can be occupied or not, leading to two (nearly) degenerate ground states characterized by fermion parity. Our study begins with a characterization of the static properties of the finite-sized wire, including the behavior of the Majorana end modes and the form of the tunnel coupling between them; a transfer matrix approach to analytically determine the locations of the zero energy contours where this coupling vanishes; and a Pfaffian approach to map the ground state parity in the associated phase diagram. We next study the quench dynamics resulting from initializing the system in a topological ground state and then dynamically tuning one of the parameters of the Hamiltonian. For this, we develop a dynamic quantum many-body technique that invokes a Wick's theorem for Majorana fermions, vastly reducing the numerical effort given the exponentially large Hilbert space. We investigate the salient and detailed features of two dynamic quantities - the overlap between the time-evolved state and the instantaneous ground state (adiabatic fidelity) and the residual energy. When the parity of the instantaneous ground state flips successively with time, we find that the time-evolved state can dramatically switch back and forth between this state and an excited state even when the quenching is very slow, a phenomenon that we term"parity blocking". This parity blocking becomes prominently manifest as non-analytic jumps as a function of time in both dynamic quantities.


Introduction
Of late, two different concepts in quantum many-body theory have elicited a surge of active research, partly stemming from experimental advances in condensed matter and cold atomic systemsthe concepts of quench dynamics  and topological order [52][53][54][55][56][57]. Quenching, or ramping, concerns initializing a system in its equilibrium configuration at some point in parameter space followed by inducing non-equilibrium behavior via dynamic tuning of one of the parameters. When the tuning occurs through a critical point separating two phases of matter, no matter how slow the tuning rate 1/τ, the diverging time scale associated with the critical point and critical exponent z always results in out-of-equilibrium dynamics in its vicinity. The quantum version of such Kibble-Zurek physics, initially studied as a thermal quench during the formation of the early Universe [4][5][6][7], offers a mine of valuable information about the critical point in question.
In the realm of topological systems, while quantum Hall systems have been hailed for their topological properties for over three decades [52], the recent attention on other systems has also been spectacular [53][54][55][56][57]. On the theoretical front, among others, two paradigm low-dimensional models have been avidly studied for their topological propertiesthe Majorana wire proposed by Kitaev [58,59], which is effectively a lattice version of a spinless p-wave superconducting wire, and the two-dimensional Kitaev honeycomb model [60].
Here we explore the synergy of these two concepts, namely quench dynamics 4 and topological order, a rich study that is still in its infancy. Such a marriage is exciting from at least two perspectivescan quench dynamics act as a probe for topological order? Can the existence of topological order lead to a different realm in nonequilibrium dynamics? Studies of such synergy [30,31,38,39,61] have just begun to explore diverse and exciting phenomena with regards to dynamic evolution of topological features. In previous work by two of the authors of this article and co-workers, reference [62], the term 'topological blocking' was coined with regards to the role played by a highlighting feature of topological order-ground state degeneracies-in quench dynamics. It is known that a system can have several topological sectors which are associated with these degeneracies and are distinguished from each other by an invariant based on a discrete global symmetry [63,64]. In quenching between a topological and non-topological phase, if the ground states in the two phases belong to different topological sectors and the Hamiltonian commutes with the global symmetry at all times, the system never reaches the ground state of the final phase. As a result, the usual expectation that if the quench is sufficiently slow (i.e., almost adiabatic), the system always remains in the instantaneous ground state is violated. This topological blocking effect was demonstrated in [62] for the Majorana wire and the Kitaev honeycomb model constrained to periodic boundary conditions.
The goal here is to explore quench dynamics in a finite-sized Majorana wire having open boundary conditions, a system that has come into the limelight for topological features that we expect to affect dynamics in a profound way. These features concern the presence of isolated zero energy Majorana fermion bound states at the wire ends within the topological phase; their possible experimental detection in the context of spin-orbit coupled wires [65][66][67][68][69][70] has garnered much attention in terms of fundamental physics as well as implications for topological quantum computation [71]. In the thermodynamic limit, these Majorana end modes together form a Dirac fermion state that can either be occupied or empty. Thus, degenerate topological sectors are identified by fermion parity. In the finite-sized system, tunnel coupling between these end modes splits the degeneracy in a manner that can be tuned by changing the parameters of the system. Here we explore the quench dynamics of tuning through a succession of parity flips of the ground state. The investigation involving open boundaries requires the formulation of new dynamic quantum many-body techniques, which we develop here. By investigating measures commonly studied in quench dynamics, we demonstrate that topological order drastically affects non-equilibrium behavior, the most dramatic signature stemming from quench-dependent switching of topological sectors.
Our presentation is as follows. In section 2, we present a brief description of the salient features of topological blocking and the highlights of this work, including our results regarding parity switching and blocking in a finite-sized Majorana wire. In section 3, we begin our detailed exposition by reviewing the p-wave superconducting wire given by Kitaev's lattice Hamiltonian. We outline a derivation of its bulk spectrum and description of its topological phase diagram based on the presence or absence of Majorana end modes. In section 4, based on a transfer matrix formalism, we analyze the fate of the Majorana end modes for a finite-sized lattice. We obtain an approximate form for the degeneracy splitting and exact solutions for contours in the phase diagram where the splitting vanishes. We then invoke Kitaev's argument based on Pfaffian methods to determine the ground state parity of the system and confirm that parity switches occur at these degeneracy points. In section 5, we begin our discussion of the quench dynamics associated with varying a parameter of the underlying Hamiltonian linearly in time. We summarize the known results in the case of periodic boundary conditions, in particular, the Kibble-Zurek scaling of the post-quench excitations and the topological blocking phenomenon. In section 6, we describe the real space time-dependent formalism that we use to study this problem numerically. We focus on two measures, namely the overlap between the time evolved state and the instantaneous ground state, which we refer to as adiabatic fidelity as in previous work [2], and the residual energy, which is the difference between the expectation value of the Hamiltonian in the time evolved state and the instantaneous ground state energy. In section 7, we extensively discuss our numerical results in detail, pinpointing the effect of topological blocking associated with multiple parity switches. In section 8, we present an overview of our study and connect it with related phenomena, such as the fractional Josephson effect in junctions of Majorana wires, as well as to experiments.
The stage is set by the concept of topological blocking, which, as mentioned in the previous section, was studied in [62]. The study involved quench dynamics in topological systems elicited by changing a parameter of the Hamiltonian to tune from one quantum phase to another. In going between a topological phase and a trivial phase, as in the Majorana wire, or between two topological phases, as in the Kitaev honeycomb model, the focus was on mismatch of degeneracies. It was shown that if the system was initialized in the ground state in a phase with higher degeneracy and tuned to one with a lower degeneracy, two to one in the former case, and four to three in the latter, the phenomenon of topological blocking would occur. Certain topological sectors characterized by topological quantum numbers, for instance, fermion parity, would inhabit the ground state in the initial phase but would have no partner in the ground state of the final phase. As a result, in tuning through a quantum critical point separating the two phases, these sectors would evolve so as to have null overlap with the final ground state no matter how slow the tuning rate, in stark contrast with Kibble-Zurek physics, where only a rate-dependent fraction of the time-evolved state overlaps with the excitation spectrum above the (gapped) final ground state. Moreover, it was shown that even if one took overlap with the instantaneous sectoral ground state, different topological sectors would show quantitatively different dynamic behavior, particularly in wave function overlap.
Here, we explore this notion of topological blocking with regards to a different but related aspect-the switching of topological sectors due to quench dynamics within a topological phase. In a Majorana wire with open boundaries, the topological degeneracy is associated with the presence of Majorana zero modes at the edges. The degeneracy is split due finite-size coupling of these edge modes, which induces fermionic parity sectors within the topological region of the phase diagram. Here we build on the notion of topological/parity blocking arising from tuning through these parity sectors.
In figure 2, we present some of our key results. Initializing the system in the ground state corresponding to a specific on-site chemical potential, and thus some fixed parity, we sweep the chemical potential to undergo several parity switches. As a measure of how closely the time evolved state tracks the instantaneous ground state, we evaluate the wave function overlap (adiabatic fidelity) associated with these two states. As seen in figure 2, the adiabatic fidelity plummets down to zero in certain chemical potential intervals that exactly correspond to the parity switched regions. The initial ground state, while being able to track some of the dynamic evolution, is thus forced to remain within its parity sector, an attribute of the topological phase. This multiple parity blocked dynamics is a dramatic, topologically induced deviation from the continuous evolution expected in quench dynamics.
In what follows, we detail several aspects leading up to this quench behavior, including the formulation of the Majorana wire model, parity switching due to coupling of Majorana wavefunctions and a real space formalism to compute the many-body dynamics.

p-wave superconducting wire
The system that forms our subject of study is a lattice version of the spinless fermionic one-dimensional p-wave superconducting wire with spinless fermions, also referred to as the Kitaev chain [58]. This system can be mapped exactly to a spin-1/2 XY model in a transverse field (i.e., a magnetic field applied along the z direction) via the Jordan-Wigner transformation [72]. The parameters of the system are the nearest-neighbor hopping amplitude w, the superconducting pairing amplitude between nearest neighbors Δ, and a chemical potential μ. This model is a paradigm system for demonstrating numerous interesting topological properties, including the existence of Majorana modes at the ends of an open chain in the topological phase. Figure 1. A finite size superconducting wire in the topological phase characterized by Majorana fermionic end modes. While their wave functions may or may not oscillate (red curves), they all decay into the bulk (envelope) over a characteristic length scale that depends on system parameters. The overlap between the decaying oscillatory wave functions of these two end modes gives rise to a tunnel splitting of the otherwise doubly degenerate zero energy states.
The Hamiltonian of such a system with N sites and open boundary conditions is given by where the λ j are non-negative real numbers. This can be done through a transformation of the form where = ⋯ a a a ā ( , , , ) are not paired with any other operators in the system and therefore do not appear in the Hamiltonian. These isolated modes correspond to the zero energy eigenvectors localized at the ends. The existence of these modes is robust even away from this extreme limit and they only disappear with the closing of the bulk gap.
The ground state of the system in the topological phase is thus doubly degenerate and has two zero energy eigenvalues corresponding to the Majorana modes. These Majorana modes can be combined to form a complex Dirac fermion state, which can be either empty or occupied. Hence, each of the degenerate ground states has a specific fermion parity and the system can be characterized by a related Z 2 -valued topological invariant. This ground state parity will play an important role in the subsequent sections. Numerical results for the (a) adiabatic fidelity  t ( ) and (b) parity of the instantaneous ground state for an even number of sites. The times at which the parity switches its sign are exactly the points where parity blocking occurs, resulting in the adiabatic fidelity plummeting down to zero. Depending on the parameters chosen, the parity after crossing the quantum critical point changes from the initial ground state parity thereby leading to parity blocking for the entire topologically trivial region.
There are three phase boundaries, indicated by dark lines in figure 3, where the bulk gap vanishes. These are the quantum critical lines across which there is a topological phase transition. In the thermodynamic limit (infinite wire) or for a closed chain, one can transform the Hamiltonian into Fourier space, and the single particle energy spectrum takes the form This spectrum has a finite superconducting gap in all the phases; the gap vanishes as one crosses one of the critical lines and reopens upon entering another phase. In the spin language, the topological phases correspond to the ferromagnetic phases of the transverse field XY model (where either the x or the y component of the spins has long range order), and the trivial phases are in the paramagnetic phase. These characteristic features of the system, namely, the topological invariant, the spectrum of the bulk and end modes, and the wave functions of the Majorana end modes have been discussed extensively in previous work (see, for example, [34,56,58]). So far, all the above mentioned characteristics of the model assume the size of the system to be much larger than the decay length of the Majorana end modes. In the next section we consider the case when the Majorana modes at the two ends have a finite overlap, giving rise to consequences such as parity blocking in quench dynamics.

Tunneling between Majorana end modes
In a Majorana wire of finite length, the two Majorana end modes are no longer completely decoupled since there is some overlap between their wave functions (figure 1); the overlap shifts their energies slightly away from zero. In the extreme limit of the topological phase with Δ = w, the end modes are exactly localized at the ends with no overlap between them; hence the effective Hamiltonian governing these modes is with J = 0, But in general, the effective Hamiltonian has an expression in terms of the (almost) zero energy eigenvectors localized at the ends of a finite length wire [58] Here J is in general a function of Δ w, μ w and N. Due to this coupling the two zero modes split in energy into a particle-hole symmetric pair of eigenenergies = ± E J , and the ground state is no longer degenerate. The Majorana end modes can be combined into non-local Dirac fermions as . The Hamiltonian can now be expressed as f † The occupation number c c † can either be zero or 1. Thus we see that the energies ±E come with corresponding eigenstates with opposite fermion parities. The sign of J decides which of these states is the ground state. The parity of the states in the bulk being fixed, the overall ground state parity is then decided by the lower one of the two split energy levels (which lie inside the bulk gap). The coupling J is a function of the chemical potential μ and oscillates, switching its sign at specific values of μ. Therefore the split energy levels cross zero at certain points in time if μ is varied linearly in time. This leads to oscillations in the overall parity of the ground state. Even though these split energy levels are exponentially smaller (for large system size N) than the energies in the bulk, the parity oscillations play a key role in the time evolution of the ground state in the quenching dynamics.
The energy splitting due to the coupling of the two Majorana end modes can be derived by evaluating the overlap between their associated wave functions. For a very long wire, the Majorana modes have zero energy; we calculate their exact wave functions by using the Heisenberg equations of motion = H a [ , ] 0 n and H b [ , ] n , where a n , b n are the two Majorana operators at site n which were denoted by − a n 2 1 and a 2n above. We obtain the difference equations for these operators as: for ⩽ ⩽ − n N 2 1 . These difference equations can be solved exactly by using z-transform methods and taking into account the form of the difference equations at the ends. The wave functions α n and β n for the a n and b n modes on site n, respectively, are of the form Using these wave functions, we can now approximately calculate the energy splitting by assuming that the wire has a finite length N and computing the overlap of the wave functions α L at the left end and α R at the right end. This gives an expression of the form Figure 4(a) shows the oscillations in the splitting as calculated above. This is an approximate calculation because we have assumed the energies to be zero in the Heisenberg equations of motion for the Majorana operators and then calculated their overlap (which shifts the energies slightly away from zero). But one can see that it qualitatively agrees with the exact numerical calculation in figure 4(b). Figure 4(c) shows the oscillations of the ground state parity as a function of μ which is linearly varied with time. One sees an excellent correspondence with the zero crossings of the energy splitting in figure 4(b). Further, although these oscillations appear due to the degeneracy splitting in the topological phase, they do not exist in the entire topological region in the phase diagram. To understand this, we stress the fact that the key ingredient in getting these oscillations is the oscillatory component in the wave functions of the Majorana modes given in equations (9). One can see that the oscillatory functions ωn sin( ) and ωn cos( ) become hyperbolic if ω given by equation (10) becomes imaginary. The boundary at which this happens is given by the circle μ Δ = − w 4( ) 2 2 2 . Beyond this circle the Majorana wave functions have only a decaying (but not oscillatory) component. This implies that there would not be any parity oscillations in this region. This can be clearly seen in figure 5(a).
We remark here that expressions for the energy splitting of the Majorana end modes are known for a continuum model [73,74]. These results are consistent with the splitting being both oscillatory and decaying exponentially with increasing length.
Calculation of ground state fermion parity:to obtain a rigorous knowledge of the parity and its switching as a function of the chemical potential, as has been used in other Majorana wire contexts [75], we employ the measure introduced in Kitaev's well-known work [58]. Given a Hamiltonian of the form in equation (1), the transformation B, which reduces the Hamiltonian to the canonical form, can be represented as a conjugation by a parity preserving unitary operator if B has the form = B e D i.e., if = B det( ) 1. Otherwise B changes the parity. Therefore, the parity of the system is given by In appendix B we illustrate this result with a simple problem of a two-site effective Hamiltonian. This illustration is particularly useful in our case since the effective Hamiltonian for only the coupled Majorana modes is in fact a two-site problem given by Within the topological phase, the dynamics of only these end modes and their associated Dirac fermions determine the overall parity as we saw above.  In terms of the Majorana operators, the parity of a N-site system is given by We note that P is both Hermitian and unitary, and it commutes with the Hamiltonian in equation (2); since = P I 2 , the eigenvalues of P must be ± 1. In the extreme case of Δ = w for the Majorana wire in equation (2), the terms of the form + a a j j 1 are equal to 1. Thus only the term a a N 1 2 remains, which is in fact the term in the effective Hamiltonian. Now, within the topological phase, small deformations of the parameters should not change this fact. Since the parity operator also commutes with the Hamiltonian, we can see that in the topological phase the end modes alone determine the parity of the ground state. The parity equivalence of ground state sectors has been studied in [76].

Exact expression for zero energy contours
While we obtained an approximate result for the tunneling amplitude that splits the degeneracy between Majorana end modes in the previous subsection, a formalism involving transfer matrices [34,77] enables to track the exact points in the parameter space at which this energy changes sign, restoring the zero energy degeneracy and resulting in a parity switch. Previous work has presented similar derivations and results using the equivalent method of chiral decomposition [78].
For a mode with energy exactly equal to zero, equations (8) and (9) are applicable. We see that the a and b modes are decoupled; for definiteness, let us consider a zero energy mode involving the a n ʼs. Equation (8) shows that for ⩽ ⩽ − n N 2 1 , the a n ʼs are related to each other by a transfer matrix M a , n n a n n a If λ 1 and λ 2 are the two eigenvalues of M a , the general solution of equation (8) is regardless of the relative values of w and Δ) and N is odd.
We get the same conditions if we look for a zero energy mode involving the b n ʼs.
In terms of ω defined in equations (10) and (16) is equivalent to saying that ω We see that this differs somewhat from the approximate condition that α α L R given in equation (11) should be equal to zero. The solutions of equation (16) are given by where p is an integer equal to ⋯ N 1, 2, , 2 if N is even and 2 if N is odd. In terms of the variables μ w and Δ w, we observe that equation (17) defines a number of ellipses, which are labeled by the integer p; these are shown in figure 6(a) for N even and figure 6(b) for N odd. Note that all the ellipses pass through the two points given by μ = 0 and Δ = ±w. Figure 6(b) for N odd also contains a zero energy line lying at μ = 0 for all values of Δ w.
Equation (17) can be understood in a simple way for the special case Δ = 0. Equation (17) then reduces to We can understand this as follows. For Δ = 0, equation (1) describes a non-superconducting tight-binding model whose single particle spectrum is given, for an open chain with N sites, by , where = ⋯ q N 1, 2, , . One of these energies vanishes whenever μ satisfies the condition given in equation (18), in particular, when p is equal to the smaller of the two integers q and . This is where the ranges of p mentioned above, namely, = ⋯ p N 1, 2, , 2 for N even and ⋯ − N 1, 2, , ( 1) 2 for N odd come from. In addition, if N is odd (but not if N even), we have a zero energy state at μ = 0 corresponding to = + q N ( 1) 2. Having found all the zero energy lines in the plane defined by μ w and Δ w, we observe that the parity of the fermion number of the ground state flips sign whenever we cross one of these lines. As a check, this is again easy to see for the case Δ = 0. The number of energy levels which are occupied in the ground state changes by 1 and hence the fermion parity changes sign whenever one of the single particle energies E q given above crosses zero. For Δ = 0 and very large negative values of μ, we can see that the ground state of equation (1) contains no fermions; hence, according to equation (13), the fermion parity is +1 for any value of N. For very large positive values of μ, the ground state of equation (1) is completely filled with N fermions; hence the fermion parity according to equation (13) is − ( 1) N . We remark that the oscillations in the parity, which are related to the Kitaev's Pfaffian, map to spin-spin correlations in the transverse spin chain, and as with much of the literature on Majorana wires, these oscillations have been discussed in depth in the spin context [79].

Quenching dynamics in the Majorana wire
Previous work involving the dynamics of quenching in the Majorana wire described above has focused on tuning through quantum critical points separating topological and trivial phases [31,34,61,62]. There have been recent works on the effect of quenching on Majorana modes, signatures of Majorana modes in quenching dynamics and Kibble-Zurek scaling [80][81][82][83][84]. (As mentioned earlier, while the terms quench and ramp are sometimes used to distinguish between instantaneous change of parameter versus a time-dependent change at some given rate, respectively, here we use quench in a more general sense to encompass all such dynamic tuning. Our actual studies are restricted to the ramp case.) While our analysis also explores non-equilibrium dynamics within a particular topological phase, we use similar protocols for changing parameters of the Hamiltonian to tune from one phase to another.
Specifically, we consider a linear variation with time of the chemical potential of the system so as to go across the critical line μ = w 2 , i i Here μ i is the initial chemical potential at t = 0 and τ 1 is the quench rate. Due to the finite rate of variation of μ, the system cannot remain exactly in its ground state and will exhibit non-equilibrium behavior. Namely, excitations (or defects) will be produced in the ground state; this lead to an excess energy of the system and may also lead to a ground state which is in a different topological sector than the initial ground state. These effects can be characterized by the following quantities.
• Defect density: the number of defects produced in the ground state configuration, which is given by the sum over all the excitations.
• Adiabatic Fidelity  t ( ): this is the inner product of the instantaneous ground state ψ | 〉 t ( ) ins of the timedependent Hamiltonian with the time evolved initial ground state Ψ ins • Residual energy E res : this is the energy in excess of the instantaneous ground state of the system. We will define this as the dimensionless quantity where E G is the energy of the ground state at time t.
Previous work: most of the analytical results for the above quantities obtained in earlier work are in the limit of very large system size or with periodic boundary conditions, where one can Fourier transform the Hamiltonian to momentum space. This in fact reduces the calculation to a well known problem of a transition between two states for each value of the momentum k in the Brillouin zone. This is the famous Landau-Zener-Majorana-Stueckelberg problem [85][86][87][88] which can, under a few assumptions, be solved exactly to obtain the probability of excitation from the ground state to the excited state. Using this probability, we can obtain expressions for the defect density, adiabatic fidelity and residual energy. In the limit of long time t, all these quantities have a universal power law scaling as a function of the quench rate which is related to the post-quench excitations. This is the well studied Kibble-Zurek scaling . Given the quench rate τ 1 , the defect density and residual energy scale as τ 1 and the adiabatic fidelity scales as τ c exp( ) in this one-dimensional Majorana wire.
For the Majorana wire with open boundary conditions, there have been investigations of the behavior of single particle states under a quench. While it has been found that the single particle bulk states still obey the Kibble-Zurek scaling for the defect density, the quench for an initial state with a Majorana end mode has been found to be non-universal and dependent on the topological features of the system [31]. The end states are not robust with respect to the quench and they delocalize to merge with the bulk states. This leads to a scaling of the defect density as τ 0 (i.e., independent of the quenching rate), which is very different from the Kibble-Zurek scaling.
Another investigation for the open chain looks into a quantity called the Loschmidt echo, which is the survival probability of the Majorana end modes under a quench [61,80]. Upon quenching across the critical point, the probability decays to extremely small values as the end modes merge with the bulk states when the system is near the critical point where the gap between the end and bulk states vanishes. But interestingly, when the system is quenched to exactly the critical point, although the probability of survival goes to zero initially, it revives itself completely at regular intervals of time. This is attributed to the nearly equal spacing of the low-lying energy levels in the bulk near the critical point; this spacing is of the order of N 1 while the scaling of the gap at the critical point is also of order N 1 . This leads to oscillations in the survival probability with a time period which is proportional to the system size N. We will see below the remnant of this effect in the residual energy for an extremely slow quench.
In this work we will study the quenching dynamics in the above mentioned quantities for an open Majorana wire. We will mainly focus on the many-body states rather than the single particle states and the novel parity switching mechanism discussed in the previous section, which comes about in the topological phase due to the coupling between the Majorana end modes.

Real space formalism for studying quenching dynamics for open boundary conditions
In comparison to the translationally invariant systems with periodic boundary conditions studied in previous work, a challenge encountered in these finite-sized systems with open boundary conditions is that one cannot exploit the momentum basis, which in previous works highly simplified the quench dynamical problem. Here,in principle, we are faced with the full-fledged 2 N -dimensional Hilbert space associated with fermions on a N-site lattice associated with the Fock space formed by fermion occupancy on each site. Here, we develop and present a dynamic many-body technique to reduce the problem to a numerically tractable form. Our technique hinges on two principles in calculating expectation values or overlaps between states in this time-dependent setting. The first is to use the Heisenberg picture so that the crux of the information on the time evolution is given by the relation between fermionic creation/annihilation operators at different times. The second is to invoke an analog of Wicks theorem for Majorana operators. The computation then reduces to dealing with time-dependent × N N 2 2 matrices, allowing us to embark on an exhaustive analysis of adiabatic fidelity and residual energies and to pinpoint attributes of parity blocked dynamics.
Let us start with a general time-dependent Hamiltonian which is quadratic in Majorana operators a j = ⋯ j N ( 1, 2, , 2 ), Here M(t) is a real antisymmetric matrix; its elements will be functions of w, Δ and μ for the Majorana wire Hamiltonian in equation (2). This can be converted to the canonical form . The N (2 )-dimensional matrix B(t) comprises of the eigenvectors of H(t) and it belongs to the group U N (2 ) with = ± det B ( ) 1. The eigenvalues of the Hamiltonian are λ ± j . Adiabatic fidelity calculation: as defined in a previous section the adiabatic fidelity is given by . The corresponding annihilation operators of these many-body states satisfy the relations Here  is the evolution operator, with  denoting time ordering. The two sets of annihilation operators are related by The key idea underlying the calculation in real space is to express the quantities of interest to us in terms of objects which can be calculated numerically in a simple way. Given the form of the initial Hamiltonian H(0) and the time-dependent H(t), we can see that the quantities B(0), B(t), S t ( , 0) and G(t) can be easily computed. Given these and the annihilation operators for the ground states, the calculation of the adiabatic fidelity reduces to a computation of the determinant of an antisymmetric matrix A given by We now directly state an important result, deferring the detailed derivation to appendix. The adiabatic fidelity defined in equation (20) is given by 1 4 Given this relation and the Hamiltonian H(t), we can numerically calculate the adiabatic fidelity as a function of time for a system with open boundary conditions. This can naturally be extended to periodic/antiperiodic boundary conditions as well.
Residual energy calculation: another quantity of interest, the residual energy, defined in equation (21), measures the excess energy contained in the time-evolved quench dependent state compared to the instantaneous ground state energy. This quantity can also be calculated with the real space formalism developed in this section. Following the same strategy as for the adiabatic fidelity, the final expression can simply be expressed in terms of the matrix G(t), The details of the derivation are given in appendix.

Results
We now present the results that we obtain for an open Majorana wire with N sites with the Hamiltonian given in equation (2), where μ varies linearly in time as shown in equation (19). Given this specific form of H(t), we numerically calculate all the quantities B(0), B(t) in equation (24), S t ( , 0), G(t) in equation (25) and finally the adiabatic fidelity  t ( ) in equation (26) and the residual energy E Res in equation (27). In what follows, we provide an in-depth discussion of the novel phenomenon of topological and associated parity blocking as elucidated in section 2.
For periodic boundary conditions, the effect of the number of sites on the fermion parity of the ground state and the consequent phenomenon called topological blocking on the quenching dynamics of the ground state has been discussed in detail in [62]. For an open chain, we saw in section 4 that the Majorana end modes play an important role in determining the parity. For a fixed Δ and N, the parity changes sign as we sweep across the topological phase by varying μ. Here we will explicitly see this parity blocking effect in the evolution of the ground state within the topological phase. The initial parity of the system is set by the value of μ i and the number of sites N. As we will see below, the choice of μ i can have drastic consequences, especially for an odd number of sites. Figures 2 and 7 show the numerical results for the adiabatic fidelity  t ( ) along with the parity of the instantaneous ground state for an open chain with an even and odd number of sites, respectively. The case of the initial value μ = 0 i for an odd number of sites will be discussed later. We can see from figures 2 and 7 that for both even and odd number of sites, the system starts in a particular fermion parity sector, and as it moves within the topological phase the instantaneous ground state switches parity regularly. On crossing the critical point it can either have opposite parity from the initial state or the same parity, depending on the initial parity sector. On the other hand, as we are dealing with parity conserved systems, the state which is time evolved from the initial ground state continues to have the same fermion parity. Thus the overlap of the time evolved state with the instantaneous ground state plunges to zero at times when the instantaneous parity becomes opposite to the initial parity. We call these parity oscillations, which occur for an open Majorana wire, as the parity blocking effect. The initial ground state is blocked from having any non-zero overlap with the instantaneous ground state for certain values of μ. Finally, on crossing the quantum critical point it becomes zero at all times if the parity is flipped from the initial ground state; this is also a manifestation of parity blocking. Hence, in figure 2, the case of an even number of sites, the parities of the initial and final ground states are the opposite and the system shows parity blocking for the entire topologically trivial phase, while in figure 7, the parities are matched and there is some residual overlap in the trivial phase.

Adiabatic fidelity and parity blocking
Domain with no parity blocking: as we showed in section 4, the oscillations in the parity do not occur throughout the parameter space corresponding to the topological phase. Namely, there are no oscillations if 2 . This implies that there ought to be no parity blocking in the adiabatic fidelity. We indeed see this in the numerical results shown in figure 8.
Fermion parity degeneracy for an odd number of sites: for an open chain with an odd number of sites, and for states which belong to the odd fermion sector, the μ = 0 i point is special in that it has two degenerate ground states. To see this, let us define an operator where the last term on the right-hand side is given by 2 , we find that C generates the particle-hole transformation In this case the system has the same parity as the initial ground state on crossing the quantum critical point (figure (b)) and therefore has a non-vanishing overlap with it.
. (29) n n N n n n N n This is a symmetry of the Hamiltonian in equation (1) if μ = 0. We now note that the parity and chargeconjugation operators P and C satisfy = − PC C P ( 1) N . Thus P and C anticommute if the number of sites N is odd. Since P and C both commute with H if μ = 0, every energy of a system with an odd number of sites will have a two-fold degeneracy with the two eigenstates having opposite fermion parities. (This can be shown as follows. If ψ | 〉 + is an eigenstate of H and P with eigenvalues E and +1 respectively, the relations is an eigenstate of H and P with eigenvalues E and −1 respectively.) Therefore, starting from μ = 0 i means starting with parity states whose degeneracy is not split. The time evolution of an arbitrarily chosen ground state would therefore be that of a linear combination of both parity states. Even though the parity of the instantaneous ground state would keep switching as one sweeps through the topological phase, the overlap would be finite as the time evolved state will be in a superposition of both parities. But the value of the adiabatic fidelity will be smaller than in the parity blocked case as the amplitude is split between the two superposed states. Figure 9 shows the results for this unique case of quenching. We can clearly see in figure 9 that the splitting is zero at time t = 0. This change in the initial condition drastically affects the evolution of the ground state and we do not see complete parity blocking in this case.
Comparison with a closed chain: the case of parity blocking in closed chains with periodic and antiperiodic boundary conditions has been studied in detail in a previous work [62]. While the previous analysis involved momentum modes, simplifying the problem to a set of two-level Landau-Zener systems, the real space formalism is easily extended here to compute all the quantities for a closed chain. The numerical results are shown in figure 10. Even though there is parity blocking in the case of an open chain, we may still expect the envelope of the adiabatic fidelity to be comparable with that of the periodic case. In figure 10, we see a good match between the envelopes in the two cases. For a closed chain, the final parity can flip from the initial state depending on the boundary conditions. As can be seen, the parity does not change for the periodic closed chain and therefore one still has a finite overlap. The match between the envelopes in the open and closed chain cases is very close. This suggests that overall the open chain would also respect the Kibble-Zurek behavior for the defect production and excitation density that was found in the closed chain case. The crucial difference between open and closed chains is the parity blocking and switching due to the coupling of the two Majorana end modes.

Residual energy
The numerical results for the variation of the residual energy with time using the full many-body formulation in equation (27) is shown in figure 11 for both open and closed chains. The two cases have the same average behavior. Both show a rapid increase in the energy of the system above the instantaneous many-body ground state as they approach the critical point. This rapid rise is due to the system falling out of equilibrium upon approaching the critical point and thus losing adiabaticity. Far beyond the critical point, we find that the energy asymptotes to a fixed average value. The scaling analysis of this quantity for the transverse field Ising chain has been studied numerically in [89].
Effect of Majorana end modes: the crucial difference between the open and closed chains in figure 11 is the presence of the abrupt jumps at small times in the case of the open chain. Reflecting the behavior of adiabatic for the odd sector. This is the special case where the initial state is in a superposition of the odd and even parity states.Thus the time evolved states will not be completely 'parity blocked' but the amplitude of adiabatic fidelity will be reduced. As we go to smaller N the splitting is exponentially enhanced and one can clearly see the effect of it in 'skewing' the superposition towards the state which contributes to the ground state.
fidelity, these jumps correspond to switching back and forth between the ground state and excited states due to parity blocking. Upon comparing the behavior of the residual energy with the energy splitting and parity switching of figure 4(b), we find that there is a complete match between the points at which the jumps in the residual energy take place and the points where the parity switching occurs.
Signatures of Loschmidt echoes: in addition to parity blocking, in our many-body system, we find evidence for the Majorana-mode related physics found in previous work on single particle dynamics in quenching [61]. The Loschmidt echo studied in this work calculates the probability for an initial Majorana mode to become a single-particle bulk excitation as a function of time as one sweeps through the critical point. If one varies μ t ( ) extremely slowly so as to be close to the adiabatic limit, the gap between the mid-gap end states and bulk states scales as N 1 on approaching the critical point since the dynamical critical exponent is equal to 1. The level spacing of the low lying bulk states also scales as N 1 . Hence the Loschmidt echo turns out to be a periodic function with period N as one quenches to or across the critical point.
In our full-fledged many-body treatment, these echoes appears as 'chirps' of excitations whose occurrence has a period of N. These excitations contribute to the overall energy of the state. Note that this is true only if the quench rate is extremely slow. As shown in figure 12, one can see these oscillations in the numerical results at a Figure 10. Numerical calculation of  t ( ) for a closed chain with 34 sites. The periodic and antiperiodic closed chains represent even and odd fermion sectors, respectively. One can see that there is blocking in the second case, whereas a small amount of overlap persists in the first case after crossing the critical point. Also the envelope of the adiabatic fidelities for closed chains is compared with that of the open chain for the same number of sites. Even though there is no 'parity blocking' within the topological phase in the case of closed chain due to absence of the edge modes, the overall the behavior remains qualitatively the same. Figure 11. Residual energy plots with the critical point occurring at t = 2. One can notice in the case of open chain the oscillations before crossing the critical point, which arise due to the oscillation of mid-gap states. One can see that the steps arising due to the splitting scales inversely with the system size. low quench rate like τ = 1 0.1. The frequency of occurrence of these chirps in the excitations is indeed halved when the system size is doubled. These chirps also appear in the adiabatic fidelity.
To summarize, open and closed chains broadly show similar average behaviors as expected for quench dynamics, However, both in the adiabatic fidelity and residual energy, distinct non-analytic features arise in the form of jumps only for an open chain, and these can be attributed directly to the presence of end modes and their associated fermion parity.

Discussion
In conclusion, our study of non-equilibrium behavior in finite-sized Majorana wires demonstrates that the presence of topological order can dramatically alter quench dynamics. Previous work involving Majorana wires having periodic boundary conditions brought to light the notion of topological blocking in tuning between different quantum phases [62]. In contrast, we have seen here that the coupling between Majorana end modes and the associated ground state parity flips as a function of a tuning parameter gives rise to a more drastic manifestation of topological blocking due to a succession of switches between topological sectors within a single topological phase. As a result, some common measures studied in the quench dynamics literature, such as wave function overlaps between time-evolved states and instantaneous ground states (the adiabatic fidelity), and residual energies, show a series of non-analytic structures in the form of characteristic jumps which are not observed in standard Kibble-Zurek physics.
Our work has shown that there is a much richer texture in the phase diagram of the Kitaev chain or the Majorana wire than has been presented earlier. The circle which separates the regions of oscillating and purely decaying wave functions of the Majorana end modes exists in the thermodynamic limit. In addition, we have shown that a coupling between the end modes, due to a finite length of the chain, leads to further divisions within the circle in the form of ellipses, each division corresponding to a particular fermion parity of the ground state. Although the energy splitting between the Majorana modes goes to zero as we increase the chain length, the number of fermion parity switches increases linearly with the length. This has a dramatic consequence for the adiabatic fidelity under a quench, namely, parity blocking occurs more frequently as we increase the system size. (This is very different from conventional finite size effects which typically vanish in the thermodynamic limit.) We therefore see that parity blocking is not merely a finite size effect, but is a relevant manifestation of the physics of Majorana modes and their topological nature in any real system. This study shows that quench dynamics serves well as a probe of topological order. While blocking features need not be unique to topological systems in that quantum invariants in other systems can possibly have similar effects, they are necessary conditions under appropriate circumstances (for instance, open boundary condition in the case studied here). Moreover, unlike in most other systems, such as ferromagnets having local order, we expect this blocking phenomenon to be robust against local perturbations. In the case of a finite-sized Majorana wire, the succession of parity switches associated with topological sectors is a crucial aspect of topological order; while studies of the static behavior have been extensive (see, for instance, [90][91][92][93]), here we have performed the nearly unexplored study of their effect on quench dynamics. In fact, the issue of parity forms the basis of several discussions and proposals for Majorana wires, particularly in light of the potential experimental discovery of isolated Majorana end modes and their implications for topological braiding and quantum computing. Several schemes involve changing the on-site chemical potential at specific locations as a means of manipulating and dynamically moving the isolated end modes. A popular study regarding the end modes is the fractional Josephson effect (see, for instance, [58,[94][95][96][97][98]), which involves parity switches between two finite-sized Majorana wires connected to each other at their ends and their effect on Josephson physics (in principle, other zero energy end bound states could mediate such an effect). Our study here is highly relevant to these lines of investigation and it provides a dynamic quantum many-body formulation that goes beyond quasi-static approximations.
While our study primarily aims to understand the effects of parity switching on issues typically studied in the literature on quench dynamics, an experimental setup probing the predictions would be remarkable. While the arena of cold atomic gases is more ideally suited for measurements of residual energy and adiabatic fidelity, realizing topological order in these systems is still in its initial phases [99]. In the setting involving spin-orbit coupled wires, where the isolated end modes have potentially been observed, in principle, ground state parity switches can be observed by coupling the wire to another system. For instance, a possible read-out could involve tunnel-coupling to a quantum dot or STM tip (see, for example, [100,101]). Further studies would involve pinpointing ways of measuring the behavior of adiabatic fidelity predicted here in such a setup.
Finally, the study of topological/parity blocking in quench dynamics presented here and the associated quantum many-body formulation offer wide scope for further exploration. Several aspects of this initial study require more detailed investigation, for instance, more involved studies of system size, and further connections with Kibble-Zurek scaling and single-particle physics, including anomalous scaling due to boundary effects and appearance of Loschmidt echoes. Oscillations have been found in the derivative of the Renyi entropy with respect to the chemical potential in [102] and it may be interesting to see if there are such effects related to the oscillations in the ground state parity. A host of open issues related to topological blocking in Majorana wires include constraints on thermalization imposed by topological order, effects of external potentials, such as quasiperiodic potentials and disorder, and higher dimensional analogs, such as the Kitaev honeycomb model. The key idea underlying the calculation in real space is to express all the quantities of interest in terms of those which can be calculated numerically. We can see that given the form of the initial Hamiltonian H(0) and the time-dependent H(t), the quantities B(0), B(t) and S t ( , 0) can be easily computed. Given these and the annihilation operators of the ground states, we will now derive the final expression for the adiabatic fidelity.
Consider the Fock space of 2 N states ϕ | 〉 a , = ⋯ a 1, 2, , 2 N . Their fermionic occupation numbers are given We will make use of this fact to reduce the calculation in the 2 N dimensional Fock space to a matrix computation in N 2 dimensions as follows. Consider the quantity