Disorder-protected topological entropy after a quantum quench

TL;DR

This paper investigates how disorder can protect topological entropy in quantum systems after a sudden quench, demonstrating that Anderson localization enhances the robustness of topological phases.

Contribution

It introduces an analytical method to study the time evolution of topological entropy post-quench and shows disorder-induced localization preserves topological order.

Findings

Disorder in the Hamiltonian couplings enhances topological entropy resilience.

Anderson localization prevents the decay of topological order after a quantum quench.

Topological phases can be protected by disorder even out of equilibrium.

Abstract

Topological phases of matter are considered the bedrock of novel quantum materials as well as ideal candidates for quantum computers that possess robustness at the physical level. The robustness of the topological phase at finite temperature or away from equilibrium is therefore a very desirable feature. Disorder can improve the lifetime of the encoded topological qubits. Here we tackle the problem of the survival of the topological phase as detected by topological entropy, after a sudden quantum quench. We introduce a method to study analytically the time evolution of the system after a quantum quench and show that disorder in the couplings of the Hamiltonian of the toric code and the resulting Anderson localization can make the topological entropy resilient.

Figures

Equations

$$H_{TC}=-\sum _s{A}_s-\sum _p{B}_p$$

$$H(J,\lambda )=-\sum _s{J}_s{A}_s-\sum _p{J}_p{B}_p-\sum _{\begin{array}{c}i\in \text{odd}\\ \text{rows}\end{array}}{\lambda }_i^z{\sigma }_i^z-\sum _{\begin{array}{c}j\in \text{even}\\ \text{rows}\end{array}}{\lambda }_j^x{\sigma }_j^x.$$

$$S_2^T(t,J,\lambda )=lo{g}_2\left(\frac{P^{ABC}(t,J,\lambda )P^B(t,J,\lambda )}{P^{AB}(t,J,\lambda )P^{BC}(t,J,\lambda )}\right)$$

$$\tilde{H}(J,\lambda )$$

$$P(t)=$$

$$\Delta S_2^T\equiv S_2^T(0)-\underset{t\to \infty }{\lim }{S}_2^T(t).$$

$$\Delta S_2^T(r_m,W_J)=g(r_m)f(W_J{r}_m^{-b}),$$

$$\tilde{\mathcal{H}}=-\sum _{l=1}^N\left(J_l{\tau }_l^z+{\lambda }_l{\tau }_l^x{\tau }_{l+1}^x\right).$$

$$\tilde{\mathcal{H}}=$$

$$\displaystyle M=\left(\begin{array}[]{cc}A&B\\ -B&-A\\ \end{array}\right)$$

$$\displaystyle\left(\begin{array}[]{cc}A&B\\ -B&-A\\ \end{array}\right)=\left(\begin{array}[]{cc}g^{T}&h^{T}\\ h^{T}&g^{T}\\ \end{array}\right)\left(\begin{array}[]{cc}\omega&0\\ 0&-\omega\\ \end{array}\right)\left(\begin{array}[]{cc}g&h\\ h&g\\ \end{array}\right).$$

$$G_{ij}(t)$$

$$G_{ij}^A(t)$$

$$G_{ij}^B(t)$$

$$a_i^H(t)=\sum _j{\varphi }_{ij}^{\ast }(t)c_j^{†}+{\varphi }_{ij}(t)c_j,$$

$$b_i^H(t)=\sum _j{\psi }_{ij}^{\ast }(t)c_j^{†}-{\psi }_{ij}(t)c_j.$$

$$G_{ij}(t)$$

$$G_{ij}^A(t)$$

$$G_{ij}^B(t)$$

$$\displaystyle\!\!\!\!\!i\frac{d}{dt}\left(\!\begin{array}[]{cc}g(t)&h^{\ast}(t)\\ h(t)&g^{\ast}(t)\\ \end{array}\!\!\right)=\left(\!\!\begin{array}[]{cc}A&B\\ -B&-A\\ \end{array}\!\!\right)\left(\!\begin{array}[]{cc}g(t)&h^{\ast}(t)\\ h(t)&g^{\ast}(t)\\ \end{array}\!\!\right).$$

$$\displaystyle\!\left(\!\!\begin{array}[]{cc}g(t)&\!h^{\ast}(t)\\ h(t)&\!g^{\ast}(t)\\ \end{array}\!\!\!\right)\!=\!\left(\begin{array}[]{cc}g^{T}&h^{T}\\ h^{T}&g^{T}\\ \end{array}\right)\left(\begin{array}[]{cc}e^{-i\omega t}&0\\ 0&e^{i\omega t}\\ \end{array}\right)\left(\begin{array}[]{cc}g&h\\ h&g\\ \end{array}\right).$$

$$\varphi (t)$$

$$\psi (t)$$

Full text

Introduction.— Novel quantum phases in many-body systems that feature topological order are of extreme importance in both condensed matter physics wenbook and in quantum information nayak:2008 . They possess gapped energy spectrum and robust ground-state degeneracy, which is supposed to be a promising candidate of the self-correcting quantum memory dennis . These novel quantum phases can not be described by the Landau paradigm of symmetry breaking and are not characterized by local order parameters. Instead, they are characterized by a long-range pattern of entanglement dubbed topological entropy (TE) hiz1 ; hamma:2005b ; kitaevpreskill ; levin:2006 that serves as nonlocal order parameter hammahaas .

In order to exploit topological order for realistic applications such as the robust quantum memory, the system needs to be robust not only in the ground state degeneracy but must also feature robustness at both the dynamical level and at finite temperature chesitherm ; iblisdir ; chamon3d . Topologically ordered systems in two and three dimensions based on local Hamiltonians with commuting operators are not stable both at finite temperature rmp2016 ; Hastings2011 ; finiteT ; chamon3d or when cast away from equilibrium kay2009 ; Yu2016 . On the other hand, both the topological phase and its self correcting quantum memory are robust in four or greater spatial dimensions, which, unfortunately, is not realistic for implementation dennis ; alicki ; rmp2016 ; hammamazac . The depletion of both topological entropy and topological quantum memory is due to the diffusion of defects that ultimately destroy both features, as they are intimately connected, although not exactly the same thingchamon3d ; hammamazac ; nussinov . Several schemes have been proposed to overcome these shortcomings, from the introduction of long-range interactions between the excitations toricboson ; chesi-mem ; Bardyn2016 , to models that feature membrane condensation membrane together with the absence of string-like excitations haah , and the introduction of localization through disordered couplings disorderTCM1 ; disorderTCM2 ; bravyi_majorana , the latter showing that disorder can increase the lifetime of quantum memory, also see rmp2016 for extended references.

In this paper, we study how disorder can improve the resilience of topological order after a quantum quench. To this end, we study the time evolution of TE in the toric code with randomized couplings. While for the clean system the TE will self-thermalize after the quantum quench Yu2016 , we show that we can obtain a stable TE with a error that can be made arbitrarily small by a disorder increasing nearly with the square root of the system size.

Quantum quench and time evolution of TE.— We consider the two-dimensional toric code model (TCM) introduced by Kitaev kitaev:2003 defined on a periodic $M\times N$ rectangular lattice, with spins $1/2$ on the bonds. The TCM Hamiltonian is given by

[TABLE]

where the star operators $A_{s}\equiv\prod_{i\in s}\sigma_{i}^{x}$ and the plaquette operators $B_{p}\equiv\prod_{i\in p}\sigma_{i}^{z}$ are stabilizer operators which belong to stars($s$) and plaquettes ($p$) on the lattice containing four spins each (see Fig. 1). The Hamiltonian is exactly solvable since any two stabilizer operators commute. So the ground space of the Hamiltonian is made of simultaneous eigenstates of all stabilizer operators with eigenvalue $+1$. Considering the global constrains $\prod_{s}A_{s}=\prod_{p}B_{p}=\mathbbm{1}$, one can see that the ground-state manifold $\mathscr{L}$ is fourfold degenerate. The logical operators encoding the topological qubits are given by ($W^{x}_{1},W^{z}_{1}$) and ($W^{x}_{2},W^{z}_{2}$) where $W^{\alpha}_{a}$ is defined as $W^{\alpha}_{a}=\prod_{j\in\gamma^{\alpha}_{a}}\sigma^{\alpha}_{a}$ ($\alpha=x,z$ and $a=1,2$) with each $\gamma^{\alpha}_{a}$ a non-contractible string winding around the torus (see Fig. 1). By defining the reference state $|0\rangle=1/\sqrt{2^{NM-1}}\prod_{s}(1+A_{s})|\Uparrow\rangle$, with $|\Uparrow\rangle$ being the all spin-up state in the $\sigma^{z}$ basis, a generic state in the ground state manifold can be expressed as $|\xi\rangle=\sum_{i,j=0}^{1}a_{ij}(W_{1}^{x})^{i}(W_{2}^{x})^{j}|0\rangle$ with $\sum_{i,j=0}^{1}a_{ij}^{2}=1$.

As the TCM system is gapped, the bipartite entanglement in the ground state satisfies the area law for the entanglement arealaw ; amico . However, there is a $O(1)$ correction of topological origin hiz1 that can be extracted by a clever linear combination of different entropies corresponding to different subsystems kitaevpreskill ; levin:2006 . This topological correction is a stable feature of the topological quantum phase hammahaas ; santra2014 . Remarkably, the same linear combination obtained by the Rényi entropies also serves as order parameter for the topological phase flammia:2009 ; halasz and it is easier to compute and experimentally accessible in principle Abanin2012 ; Zoller2012 , especially for the case of 2-Rényi entropy nature15750 : $S^{R}_{2}=-\log_{2}\mbox{Tr}[\rho^{2}_{R}]$, where $P=\mbox{Tr}[\rho^{2}_{R}]$ is the purity of the state reduced to the subsystem $R$, namely $\rho_{R}$.

One can directly calculate the von Neumann entropyhamma:2005b and the 2-Rényi entropy flammia:2009 ; santra2014 for a simply connected region $R$ of an arbitrary state in ground-state manifold of TCM to obtain $S^{R}_{von}=S^{R}_{2}=L_{R}-1$. Here $L_{R}$ is the boundary length of the region $R$ while the sub-leading constant is topological entropy characterizing the topological phase. Following Ref.levin:2006 ; flammia:2009 , the topological Rényi entropy is defined as $S^{T}_{\alpha}=S^{AB}_{\alpha}+S^{BC}_{\alpha}-S^{B}_{\alpha}-S^{ABC}_{\alpha}$, where the subsystem A, B and C are illustrated in Fig.2. In the ground state of TCM, we have $S^{T}_{2}=2$.

In order to explore the non-equilibrium time evolution of $S^{T}_{2}$ at zero temperature, we study the scenario of the sudden quantum quench. We prepare the initial quantum memory as one ground state of TCM at t=0: $|\Psi(0)\rangle=1/\sqrt{2}(1+W^{x}_{1})|0\rangle$ note1 , then we suddenly change the Hamiltonian to

[TABLE]

$J_{s}$ and $J_{p}$ are the stabilizer strengths depending on their space positions, which is random in general. The external fields are arranged in a special fashion (see Fig. 1). Superimposing together the original lattice and the dual lattice, we place the field in the $z$ direction with magnitude of $\lambda^{z}_{i}$ on the odd rows (black dots) and the field in the $x$ direction with $\lambda^{x}_{j}$ on the even rows (white dots). The wave function will undergo a unitary time evolution $|\Psi(t)\rangle=e^{-itH(J,\lambda)}|\Psi(0)\rangle$. Note that $|\Psi(t)\rangle$ is always in the sector of $W^{x}_{1}=1$ and $W^{z}_{2}=1$ note1 . The time-dependent topological 2-Rényi entropy is given by

[TABLE]

and therefore the behavior of $S^{T}_{2}(t)$ is determined by the purity for each subsystem, namely, $P^{R}(t)=\mbox{Tr}_{R}\left[\mbox{Tr}_{\bar{R}}\left(e^{-iHt}|\Psi(0)\rangle\langle\Psi(0)|e^{iHt}\right)\right]^{2}$.

The Hamiltonian Eq.(2) can be divided into two mutual commutative parts ${H=H_{1}+H_{2}}$, where $H_{1}=-\sum_{s}J_{s}A_{s}-\sum_{i\in\text{odd}~{}\text{rows}}\lambda^{z}_{i}\sigma^{z}_{i}$ and $H_{2}=-\sum_{p}J_{p}B_{p}-\sum_{j\in\text{even}~{}\text{rows}}\lambda^{x}_{j}\sigma^{x}_{j}$. We can map the stabilizer operators to the effective spins living on the lattice and dual-lattice sites, which means $A_{s}\mapsto\tau^{z}_{s}$ and $B_{p}\mapsto\tau^{z}_{p}$. In this spin ‘$\tau$-picture’, the external fields $\sigma^{z}_{i}$ and $\sigma^{x}_{j}$ flip their two neighbour effective spins, thus $\sigma^{z}_{i}\mapsto\tau^{x}_{s}\tau^{x}_{s^{\prime}}$ and $\sigma^{x}_{j}\mapsto\tau^{x}_{p}\tau^{x}_{p^{\prime}}$, where i labels the bond between two neighboring sites $\langle s,s^{\prime}\rangle$ on the lattice while j labels the bond between $\langle p,p^{\prime}\rangle$ on the dual-lattice. The corresponding Hamiltonian of Eq.(2) in the ‘$\tau$-picture’ is the sum of total $2M$ independent Ising chains with periodic boundary conditions in either odd or even rows:

[TABLE]

In presence of random couplings $J_{s},J_{p}$, this model features Anderson localization for the $\tau$ degrees of freedom, which correspond to the anyonic excitations of the toric code Anderson1958 ; kitaev:2003 . Mapping into the ‘$\tau$-picture’, the initial state $|\Psi(0)\rangle\mapsto|\tilde{\Psi}(0)\rangle$ turns out to be the all spin-up state and the time-evolution state $|\Psi(t)\rangle\mapsto|\tilde{\Psi}(t)\rangle$ results in the tensor-product state of each rows.

Although in the $\tau$-picture the state can be expressed as the tensor product of the $2M$ rows, that does not imply that the system is not entangled in two dimensions in the original spin degrees of freedom, namely the $\sigma$-picture. We employ the $\tau$-picture as a tool to compute correlation functions that enter in the computation of $S^{T}_{2}$, as it was shown in Yu2016 . In the $\tau$-picture, the formula for the purity of any subsystem $R$ reads as

[TABLE]

In the above formula, $\Phi_{1}(t)$ and $\Phi_{2}(t)$ describe the time evolution of the system in odd rows and even rows respectively. The operators $\tilde{g},\tilde{z},\tilde{h},\tilde{x},\partial\tilde{g}\partial\tilde{x}(\partial\tilde{g})$ and $\partial\bar{g}\partial\bar{x}(\partial\bar{g})$ represent string operators in ‘$\tau$-picture’ operating with the Pauli algebra in either subsystem $A$ or its complement. The definition of these operators is detailed in the supplemental material. The phase factor takes into account whether the two operators commute or anti-commute. The constant $C_{P}$ is irrelevant when we compute the $S^{T}_{2}$ by Substituting Eq.(5) into Eq.(3). As one can see, the calculation requires all the knowledge of the many-spin correlation functions in the summation. In the supplemental material, one can find the details for the calculation of all these correlation functions, even for the disordered system. Through mapping to free fermions, the calculation of correlation functions is mapped onto the evaluation of a Pfaffian, which can be reduced to a determinant whose maximal dimension is $2(R+1)$ barouch:1971 . Only the evaluation of the Pfaffian has to be performed by a computer. The complexity of the problem only resides on the fact that Eq.(5) contains a number of correlation functions that is exponential in $(R+r)$, but each of these correlation functions can be evaluated very efficiently.

Main result.— After this effort, one can obtain the value of the purity Eq.(5) and the topological 2-Rényi entropy Eq.(3). For the clean case, that is $J_{s}=J_{p}=1$ and $\lambda_{s}=\lambda_{p}\neq 0$, the system is integrable, and $S^{T}_{2}(t)$ vanishes in long time evolution after a quantum quench therefore reaching the value of thermal equilibrium Yu2016 . At long times, the encoded logical qubit will be lost. Moreover, the Loschmidt echo $|\langle\Psi(t)|\Psi(0)\rangle|^{2}$ decays Yu2016 together with the effective magnetization in TCM kay2011 .

We have mentioned that previous results rmp2016 indicated an improvement of the lifetime for the topological encoded qubit. Now we set out to show our main result, that the dynamical localization of the anyons due to the disorder in the couplings makes the TE resilient. As we shall see, it is very important to see how disorder must scale with the system size and with the amount of protection desired.

The coupling strengths $J_{s}=1+\delta J_{s}$ and $J_{p}=1+\delta J_{p}$ are randomized as $\delta J_{s}$ and $\delta J_{p}$ are uniformly distributed in $[-W_{J},W_{J}]$. The external-field strengths are set to be fixed $\lambda_{s}=\lambda_{p}=0.6$. We set the system size as $N=200$, and the subsystem with extension $R=8$ and thickness $r=2$. On the other hand, $r_{m}$ is made to scale and $M$ can be any number proportional to $r_{m}$. Given a disorder strength $W_{J}$, the resulting TE as a function of time $S^{T}_{2}(t)$ were averaged over 100 realizations of the disorder. A simplification (without loss of generality) is used here, that each row has the same arrangement of stabilizer strengths in every single realization, so to reduce the total number of correlation functions to compute.

The results are displayed in Fig. 3, where we show the time evolution of $S^{T}_{2}(2)$ for $r_{m}=9$, and different disorder strengths $W_{J}$ from [math] to $50$ in interval of $5$. We can see that after a very short time the average $S^{T}_{2}(t)$ tends to equilibrate with small fluctuations. As the disorder strength $W_{J}$ increases, the equilibrium value $S^{T}_{2}(t)$ is closer to the initial value. The main message of the result is that Anderson localization induced by the disordered stabilizer strength makes the topological order resilient after a quantum quench.

One here needs to be careful with the scaling of these quantities with the system size. It turns out that the dependence is on the size of the ‘hole’ $r_{m}$ in the way the subsystem is partitioned. We define the quantity

[TABLE]

indicating the drop of TE after infinite time. As $r_{m}$ increases, this drop also increases, as it is shown in Fig. 4. However, we can see that as the disorder is increased, one can make this drop arbitrarily small. To determine the behavior of $\Delta S^{T}_{2}$ with large $r_{m}$ and disorder strength $W_{J}$, we perform a scaling collapse with a scaling function taking the form

[TABLE]

where $b$ is the scaling parameters, and $f$ is an undetermined function. We choose $g(r_{m})=\frac{1}{2}\Delta S^{T}_{2}(r_{m},W_{J}=0)$. Numerical fitting shows that $g(r_{m})$ satisfies a exponential function: $g(r_{m})=\frac{1}{2}(q_{0}e^{-p_{0}r_{m}}+a_{0})$, where $q_{0}=-2.73$, $p_{0}=0.91$, and $a_{0}=2$. The inset to Fig. 3 shows the collapsed data and the obtained scaling parameter, where $x_{r_{m}}=W_{J}r^{-b}_{m}$ and $y_{r_{m}}=\Delta S^{T}_{2}/g(r_{m})$. The collapsed data shows that in order to keep TE unchanged when $r_{m}$ increases, one needs to increase disorder strength grows as $W_{J}\sim r^{b}_{m}$, where $b=0.48$.

For $x_{r_{m}}\gg 1$, the fitted function $f$ takes the form $f(x_{r_{m}})=16.9/(x_{r_{m}}^{1.23}+7.23)$, which shows that one can obtain an arbitrarily good protection of the TE by setting $\Delta S^{T}_{2}(r_{m},W_{J})<\epsilon$ and then scaling the disorder strength with $r_{m}^{b}h(\epsilon)$, where $h(\epsilon)=(16.9/\epsilon-7.23)^{0.81}$ from the fitted funtion.

Conclusions and Outlook.— In this paper, we investigated the resilience of the topological phase after a quantum quench by means of dynamical localization induced by disorder in the couplings of the Hamiltonian. We use the topological entropy as a order parameter for the phase. In the toric code with random couplings and a special arrangement of the external fields, the system can be cast in the form of free fermions, and disorder will induce Anderson localization for the anyonic excitations of the system. We have shown that the phase can be protected arbitrarily well at arbitrarily long times by scaling the disorder strength with the square root of the system size.

If we let the external fields in Eq.(2) take a more general form, an interacting terms for the fermions appears. Thus strong disorder should give rise to many-body localization (MBL) huse1 ; gogolin ; Bloch15 ; Choi16 ; Monroe16 . MBL features a very slow dynamic and absence of thermalization at reasonable times abanin2 ; abaninlog ; junyang and preserves initial information to some extent. With low disorder, the system is supposed to thermalize. However, if the interaction strength corresponding to the external field is small enough, the localization still trumps the propagation Basko2005 . So it is possible that the protection of topological phase will hold for small enough general external field. The competition of the disorder strength in the stabilizer couplings and interaction strength from external field leads to the MBL phase transition. The crossover study on topological order and MBL phase transition will be the scope of future research.

I Supplemental material

I.1 Definition of the operators in the purity formula

In this section we provide the definition of the groups of string operators in the purity formula Eq.(5) in the main text, namely $\partial G^{\prime}_{R}$, $G^{\prime}_{R}$, $Z_{R}$, $H^{\prime}_{R}$, $X^{\prime}_{R}$ and $\partial G^{\prime}_{\bar{R}}$. Considering the mapping from the physical spins, namely ‘$\sigma$-picture’, to the ‘$\tau$-picture’ is isotropic in the Hilbert space we concern, we will not distinguish the notations in the two pictures. We first define the operators in the ‘$\sigma$-picture’, then map them to the ‘$\tau$-picture’ for further calculation Yu2016 .

Given a subsystem $R$, $X^{\prime}_{R}$ is the group generated by the $\sigma^{x}$ operators on the lattice bonds belonging to $R$ (on the even rows). $G^{\prime}_{R}$ is formed by the products of star operators, such that every element of the group acts only on the subsystem $R$. For a compact lattice, One can also define $X^{\prime}_{\bar{R}}$ and $G^{\prime}_{\bar{R}}$ for the complement subsystem of $R$ in the same fashion. $H^{\prime}_{R}$ is formed by all the products of plaquette operators which commute with any element in group $X^{\prime}_{\bar{R}}$. $Z_{R}$ is generated by the $\sigma^{z}$ operators on the odd rows which commute with any element in $G^{\prime}_{\bar{R}}$. As an example, we show the generators of $X^{\prime}_{R}$ and $Z_{R}$ in Fig. 5 for the subsystem being $AB$, where each white (black) dot represents a $\sigma^{x}$ ($\sigma^{z}$) operator. So far, as the element in each group, $\tilde{z}\in Z_{R}$, $\tilde{x}\in X^{\prime}_{R}$, $\tilde{g}\in G^{\prime}_{R}$ and $\tilde{h}\in H^{\prime}_{R}$ are well defined.

The definitions of $\partial\tilde{g}\partial\tilde{x}(\partial\tilde{g})$ and $\partial\bar{g}\partial\bar{x}(\partial\bar{g})$ are more complicated, so let us do it step by step. First, we denote by $\partial\tilde{g}$ a product of star operators that does not belong to either $G^{\prime}_{R}$ or $G^{\prime}_{\bar{R}}$. This means that $\partial\tilde{g}$ acts on both $R$ and $\bar{R}$. Second, $\partial\tilde{x}$ is defined as an element of $X^{\prime}_{\bar{R}}$. Third, $\partial\tilde{x}$ is a function of $\partial\tilde{g}$, such that $\partial\tilde{g}\partial\tilde{x}(\partial\tilde{g})$ acts only on the subsystem $R$. And finally, $\partial\tilde{g}\partial\tilde{x}(\partial\tilde{g})$ can not belong to the group $G^{\prime}_{R}\times X^{\prime}_{R}$. The set of the $\partial\tilde{g}\partial\tilde{x}(\partial\tilde{g})$’s satisfying all the above four conditions forms a group. Because $\partial\tilde{x}$ is a function of $\partial\tilde{g}$, This group is mapped to $\partial G^{\prime}_{R}$ by $\partial\tilde{g}\partial\tilde{x}(\partial\tilde{g})\mapsto\partial\tilde{g}$. We can get $\partial\bar{g}\partial\bar{x}(\partial\bar{g})$ and $\partial G^{\prime}_{\bar{R}}$ by changing $R$ to $\bar{R}$.

I.2 Calculation of the correlation functions

In this section, we show the details of the calculation of the correlation functions necessary for the evaluation of the purity Eq.(5). The time evolved state after the quantum quench is expressed as a tensor-product over all the rows in the $\tau$-picture. Now we show how to find the state for every row, by mapping each Ising chain into free fermions through Jordan-Wigner transformation. We consider an arbitrary row and label the ‘$\tau$’ spins with $l$ from $1$ to $N$. The corresponding Ising Hamiltonian with random couplings randomising is

[TABLE]

The Jordan-Wigner transformation maps the spin operators into the fermion operators by $\tau^{z}_{l}=1-2c^{\dagger}_{l}c_{l}$ and $\tau_{l}^{+}=(\tau_{l}^{-})^{\dagger}=\prod_{j=1}^{l-1}(1-2c^{\dagger}_{j}c_{j})c_{l}$. We can rewrite the above Hamiltonian as

[TABLE]

Note that in the fermion representation the all spin-up state is mapped to the vacuum state denoted as $|\tilde{\psi}(0)\rangle$, which correspond to the initial state. The state evolves as $|\tilde{\psi}(t)\rangle=e^{-it\tilde{\mathscr{H}}}|\tilde{\psi}(0)\rangle$, and the parity of the number of fermions is conserved in the sector we consider, so $\text{exp}(i\pi\sum_{l=1}^{N}c_{l}^{\dagger}c_{l})=1$. This quadratic fermion Hamiltonian can be diagonalized by the canonical transformation. We write Eq.(9) in matrix form as $\tilde{\mathscr{H}}=\frac{1}{2}C^{\dagger}MC$, where $C^{\dagger}=(c^{\dagger}_{1},c^{\dagger}_{2},\cdots,c^{\dagger}_{N},c_{1},c_{2},\cdots,c_{N})$, and

[TABLE]

where $A$ and $B$ are $N\times N$ matrix with elements given by $A_{i,i}=2J_{i}$, $A_{i,i+1}=-\lambda_{i}$, $A_{i+1,i}=-\lambda_{i}$, $B_{i,i+1}=\lambda_{i}$, $B_{i+1,i}=-\lambda_{i}$. It is worth noting that the sign of boundary terms needs to be changed caused by the even parity of particle numbers. The matrix $M$ can be diagonalized numerically after the orthogonal transformation

[TABLE]

The spectrum of diagonal matrix $\omega$ consists of the elementary excitations of $\tilde{\mathscr{H}}$.

To calculate the many-spin correlation function in Eq.(5), we can define ‘majorana like’ operators $a_{j}=c_{j}^{\dagger}+c_{j}$ and $b_{j}=c_{j}^{\dagger}-c_{j}$ with $a^{2}_{j}=1$ and $b^{2}_{j}=-1$. Then we represent the spin operators as $\tau^{z}_{j}=a_{j}b_{j}$ and $\tau^{x}_{j}\tau^{x}_{j+1}=b_{j}a_{j+1}$. It is revealed that the modula square of the correlation function has the form of $|\langle\tilde{\psi}(t)|\cdots a_{l}\cdots a_{s}\cdots b_{u}\cdots b_{v}\cdots|\tilde{\psi}(t)\rangle|^{2}$, or $|\langle\tilde{\psi}(0)|\cdots a^{H}_{l}(t)\cdots a^{H}_{s}(t)\cdots b^{H}_{u}(t)\cdots b^{H}_{v}(t)\cdots|\tilde{\psi}(0)\rangle|^{2}$ written in Heisenberg picture. Owing to Wick’s theorem, this many-fermion correlation function can be expressed as a Pfaffian barouch:1971 . So All we need to know, in the end, are three types of two-point correlation function:

[TABLE]

As mentioned before, $|\tilde{\psi}(0)\rangle$ is the vacuum state, namely $c_{j}|\tilde{\psi}(0)\rangle=0$ for any $j$. So we can expand $a^{H}_{i}(t)$, $b^{H}_{i}(t)$ by

[TABLE]

Substitute the above equations in to Eq.(I.2), we get

[TABLE]

The matrixes $\phi(t)$ and $\psi(t)$ can be evaluated by solving the Heisenberg equation $i\frac{d}{dt}c^{H}_{i}(t)=\left[c^{H}_{i}(t),\tilde{\mathscr{H}}\right]=\sum_{j}A_{ij}c^{H}_{j}(t)+B_{ij}c^{H}_{j}(t)^{\dagger}$. Together with the expansion $c^{H}_{i}(t)=\sum_{j}g_{ij}(t)c_{j}+h^{\ast}_{ij}(t)c^{\dagger}_{j}$, we get the matrix equation

[TABLE]

We can solve the above equation by substituting Eq.(21) into it and combining with the initial condition $g(t)=I$ and $h(t)=0$. The solution is

[TABLE]

Applying $\phi(t)=g(t)+h(t)$ and $\psi(t)=g(t)-h(t)$, we have

[TABLE]

where $\phi=g+h$ and $\psi=g-h$. Combining with Eq.(I.2) and Eq.(I.2) we can compute the pfaffian of each may-spin correlation function by reducing it to a determinant.