# Topology by Dissipation in Atomic Quantum Wires

###### Abstract

Robust edge states and non-Abelian excitations are the trademark of topological states of matter, with promising applications such as “topologically protected” quantum memory and computing. While so far topological phases have been exclusively discussed in a Hamiltonian context, we show that such phases and the associated topological protection and phenomena also emerge in open quantum systems with engineered dissipation. The specific system studied here is a quantum wire of spinless atomic fermions in an optical lattice coupled to a bath. The key feature of the dissipative dynamics described by a Lindblad master equation is the existence of Majorana edge modes, representing a non-local decoherence free subspace. The isolation of the edge states is enforced by a dissipative gap in the p-wave paired bulk of the wire. We describe dissipative non-Abelian braiding operations within the Majorana subspace, and we illustrate the insensitivity to imperfections. Topological protection is granted by a nontrivial winding number of the system density matrix.

Topological properties can protect quantum systems from microscopic details and imperfections. In condensed matter physics this is illustrated by the seminal examples of the quantum Hall effect and the recently discovered topological insulators KaneMele (); Bernevig (); Zhang (); Nayak08 (); FuKane (); HasanKane (). The ground state of the Hamiltonian of such systems is characterized by nonzero values of topological invariants which imply the existence of robust edge states in interfaces to topologically trivial phases. Due to their topological origin, these modes are immune against a wide class of perturbations.

The conceptually simplest example illustrating these phenomena is Kitaev’s quantum wire representing a topological superconducting state supporting Majorana fermions as edge states Kitaev00 (). The pair of Majorana edge modes represents a nonlocal fermion which is a promising building block to encode topological qubits Sau (); Beenakker (); Alicea11 (). Similar to the Majorana excitations near vortices of a superconductor ReadGreen (); Ivanov (), they show nonabelian exchange statistics when braided in 1D wire networks Alicea11 ().

Remarkably, the above described topological features and phenomena not only occur as properties of Hamiltonians, but appear also in driven dissipative quantum systems. Below we will develop such a topological program for a dissipative many-body system parallel to the Hamiltonian case. We will do this for a dissipative version of Kitaev’s quantum wire. This represents the simplest instance exhibiting the key features such as dissipation induced Majorana edge modes, decoupled from the dynamically created p-wave superfluid bulk by a dissipative gap. The dissipation induced topological order is generated in stationary states far away from thermodynamic equilibrium which are not necessarily pure, i.e. described in terms of a wave function only, and is reached exponentially fast from arbitrary initial states. This is in marked contrast to recent ideas of topological order in Hamiltonian systems under non-equilibrium periodic driving conditions Lindner11 (); Kitagawa10 (). Such a system can be realized with cold atoms in optical lattices, where the generation of topological order in the more conventional Hamiltonian settings has been proposed recently in a variety of settings Zhang07 (); Stanescu09 (); Goldman10 (); Bermudez10 (); Jiang10 ().

## I Dissipative edge modes in a fermionic quantum wire

Our goal is to develop a master equation for a dissipative quantum wire which exhibits topological properties including Majorana edge states. To illustrate the analogies and differences to the Hamiltonian case, and in particular to motivate our construction of the master equation with the topological states as dark steady states, we start by briefly summarizing Kitaev’s model of the topological superconductor.

Topological quantum wire – Kitaev considers spinless fermions on a finite chain of sites described by a Hamiltonian

with a hopping term with amplitude , a pairing term with order parameter , and a chemical potential . The topologically non-trivial phase of the model is best illustrated for the choice of parameters and , where the Hamiltonian simplifies to

(1) |

Here we have defined Majorana operators as the quadrature components of the complex fermion operators with properties . The Hamiltonian is readily diagonalized in terms of fermionic Bogoliubov quasiparticle operators , where importantly the pairing of Majoranas is from different physical sites.

The ground state satisfies the condition for all . The bulk of the wire describes a fermionic -wave superfluid with a bulk spectral gap, which here equals the constant dispersion . For a finite wire, the absence of the term indicates the existence of a two-dimensional zero energy non-local fermionic subspace spanned by . While highly delocalized in terms of the original complex fermion, in the real Majorana basis the situation is described in terms of two Majorana edge modes () which are completely localized on the leftmost (rightmost) Majorana site (), describing “half” a fermion each. These edge modes exist in the whole parameter regime , however leaking more and more strongly into the wire when approaching the critical values. Their existence is robust against perturbations such as disorder, which can be traced back to the bulk gap in connection with their topological origin Kitaev00 ().

Dissipative topological quantum wire – Consider the master equation for spinless fermions in a 1D chain with sites,

(2) |

with the system density operator. The two terms on the right hand side are a Hamiltonian and a dissipative term, respectively. Here we concentrate on purely dissipative dynamics which occurs at rate , and set . We choose the Lindblad operators as the Bogoliubov operators defined above,

(3) |

These Lindblad operators are quasi-local superpositions of annihilation and creation operators (see Fig. 1). Due to the fermionic nature of , this choice ensures that the bulk of the system “cools” under the above dynamics to the unique pure state , which by construction agrees with the -wave superfluid ground state of the Hamiltonian (1). Following Diehl08 (); Verstraete09 (); Diehl10b (), this steady state is thus a many-body dark state of the Liouvillian, .

The approach to this steady state is governed by the damping spectrum of the Liouvillian . Diagonality of in the implies a flat damping spectrum in analogy to the excitation spectrum of the Hamiltonian above. While the damping spectrum is always positive semi-definite for fermions, this “dissipative gap” implies exponentially fast approach of all observables to their steady state values.

For a finite wire we find dissipative zero modes related to the absence of the Lindblad operator . More precisely, there exists a subspace spanned by the edge-localized Majorana modes , with the above Fock basis , which is decoupled from dissipation, i.e. with .

Edge modes as nonlocal decoherence free subspace – These dissipative edge modes are readily revealed in solutions of the master equation. Eq. (2) is quadratic in the fermion operators, which implies solutions in terms of Gaussian density operators . Here we have defined a column vector of the Majorana operators, and is a real antisymmetric matrix related to the correlation matrix , which equally is real antisymmetric. Writing the Lindblad operators in the Majorana basis, , such that the Liouvillian parameters are encoded in a hermitian matrix , this covariance matrix is seen to obey the dissipation-fluctuation equation Eisert11 ()

(4) |

with real matrices and . Physically, the matrix describes a drift or damping, while the matrix is related to fluctuations in steady state, determined by . Writing , the approach to steady state is governed by , i.e., the eigenvalues of the positive semi-definite matrix Prosen08 () give the damping spectrum. The “dark” nonlocal subspace of edge modes, decoupled from dissipation, is thus associated with the subspace of zero eigenvalues of the damping matrix . In a spectral decomposition , and identifying by greek subscripts the zero eigenvalues subspace, we can write by partitioning

(5) |

While the bulk ( sector) damps out to the steady state by dissipative evolution, the density matrix in the edge mode subspace ( sector) does not evolve, preserving its initial correlations. The coupling density matrix elements (mixed sectors) damp out according to ; in the presence of a dissipative gap as in the example above, this fadeout of correlations is exponentially fast, leading to a dynamical decoupling of the edge subspace and the bulk. More generally, this structure of the master equation appears whenever there exists a basis in which each Lindblad operator is block diagonal with blocks associated to edge and bulk, and with vanishing entries in the edge block (see appendix).

In summary, we arrive at the physical picture that dissipative evolution cools the bulk into a p-wave superfluid and thereby isolates the edge mode subspace, , providing a highly nonlocal decoherence free subspace Lidar98 (). A physical implementation of the master equation (2) with cold atoms is outlined schematically in Fig. 2. More details on the setup are given in the implementation section, following ideas of Ref. Diehl08 ().

## Ii Stability of edge mode subspace

Here we study the robustness of the edge mode subspace against (i) global parameter changes in the Lindblad operators, while preserving their translation invariance, and (ii) static disorder, which breaks this invariance. We consider two examples of quadratic master equations (2) with Lindblad operators deviating from the ideal case (3),

(6) | |||

(7) |

where the ideal case corresponds to . In the first case, the steady state of the bulk remains pure. In the second case, we find a mixed state while still preserving the properties of the edge subspace. As elaborated on in the appendix, this results from the fact that the first case (c) represents a canonical transformation up to normalization of the ideal Lindblad operators in momentum space, while the second one (n) is not (cf. Eq. (47) in the appendix). In this latter case, the steady state has no counterpart as a ground state of some Hamiltonian. This difference is illustrated in Figs. 2 a), b). There, we plot the purity spectrum in steady state (with the edge mode subspace initialized as pure); a pure state in the bulk is indicated by all eigenvalues of being equal to . Note that static disorder, implemented in terms of small random variations of the Lindblad parameters of range from site to site, makes the first case non-equivalent to a canonical transformation, and, therefore, degrades the purity of the steady state.

While the purity spectrum is qualitatively different for both kinds of parameter deformations, the spectra of damping matrices, , are rather similar. The characteristic features are (i) the existence of a dissipative gap , closing at , and (ii) the existence of two zero modes throughout the parameter space. The associated orthogonal eigenvectors describing the Majorana modes can be constructed explicitly (see appendix), with localization length given by in both cases, in units of the lattice constant . This shows the characteristic exponential edge localization of Majorana modes close to the ideal case, while their extent becomes comparable to the system size close to the gap closing points (see Fig. 3). Adding disorder modifies the bulk spectrum quantitatively of , while the zero mode subspace persists. In fact, the existence of two zero eigenvalues in the presence of disorder can be established by explicit construction (see appendix).

In summary, a two-dimensional decoherence free subspace exists for the entire parameter range notwithstanding disorder or the fact that the bulk steady state may be strongly mixed. A sensible notion of protection of the subspace in terms of dissipative isolation from the bulk, however, only exists sufficiently far away from the point where the damping gap closes.

## Iii Adiabatic parameter changes and dissipative braiding

The above can be generalized to a dissipative quantum wire network of finite chains, as counterpart of the Hamiltonian networks discussed by Alicea et al. Alicea11 (). This results in higher dimensional non-evolving subspaces of dimension , again governed by the structure given by Eq. (5). As we will show below, this leads to the possibility of braiding of the dissipative Majorana edge modes by adiabatic parameter changes in the Liouvillian.

We consider the time evolution of the density matrix in a co-moving basis which follows the decoherence free subspace of edge modes, i.e. preserves the property . Demanding normalization of the instantaneous basis for all times, , this yields

(8) |

with the unitary connection operator and the time evolution in the instantaneous basis. The Heisenberg commutator clearly reflects the emergence of a gauge structure Berry84 (); Simon83 (); WilczekZee83 (); Carollo2003a (); Pachos99 () in the density matrix formalism, which appears independently of what kind of dynamics – unitary or dissipative – generates the physical time evolution, represented by the second contribution to the above equation. The transformation exerted on the zero mode subspace of either Hamiltonian or Liouvillian with an initial condition is then given by , with time-ordered and .

Of central importance for such state transformations to work without losing the protected subspaces is adiabaticity of the parameter changes. Here, this is a requirement on the ratio of the rate of parameter changes versus the bulk dissipative gap . This separation of time scales, naturally provided due to the non-evolving subspace, prevents the protected decoherence-free subspace from ever being left, a phenomenon sometimes referred to as the Quantum Zeno effect Beige00 ().

Equipped with this understanding of the role of the dissipative gap, we now demonstrate how to move a Majorana fermion adiabatically along the wire in our dissipative setup (cf. Fig. 4 a)), which is the key ingredient to perform braiding. As an example, we describe the move of the right unpaired Majorana fermion from site () to site (). To achieve this purpose, we consider an adiabatic change of the last Lindblad operator of the form: , where adiabatically varies with time from where to where . At the end of this evolution, the site is empty (vacuum), and the right Majorana fermion moves one site to the left. This is the analog of locally tuning the chemical potential in the Hamiltonian setting to move Majoranas Alicea11 (). For the evolution of the Majorana mode population induced by a finite ramping velocity for time , we find , describing a weak dephasing of the Majorana mode (see appendix).

Operating this mechanism on a T-junction Alicea11 () in order to exchange the two modes adiabatically while permanently keeping them sufficiently far apart from each other, the unitary braiding matrix describing the process is for two Majorana modes , demonstrating non-abelian statistics since for . While knowledge of the braiding matrix for each pair of Majorana modes is in principle enough to demonstrate non-abelian statistics, braiding of a single pair of Majorana modes is not physically observable since super-selection rules dictate a diagonal density matrix for a single complex fermion, which therefore is insensitive to the acquired phase.

We follow Bravyi (); Freedman () to construct an interferometric experiment that explicitly shows the non-abelian nature of dissipative braiding of Majoranas fermions. The setup, cf. Fig. 4 b), is given by two finite wires which host four unpaired Majorana modes, and for the left and right unpaired Majorana fermions of the first wire, and similarly and for the second one. These four real fermions correspond to two complex fermions with . For the subspace with an even number of complex fermions, we define the following basis and . We now prepare an intial state of the system as (the corresponding density matrix is ) such that . If we now braid the Majoranas and , the initial state transforms into , where . To distinguish the states and in an occupation number measurement, we first make a unitary change of basis via such that

(9) |

and then measure the number of fermions and on the wires. In the first case, we get with equal probabilities either or , while we always get in the second case, cf. Fig. 4 b).

## Iv Topological Order of the Steady State

We now show that the robustness of the edge modes seen above is indeed related to the existence of topological order in the bulk of a stationary state of Liouvillian evolution, and construct a topological invariant, which characterizes topologically different states. This classification does not rely on the existence of a Hamiltonian or on the purity of a state (in contrast to existing constructions involving ground states of Hamiltonians), and can be entirely formulated in terms of a density matrix.

In an infinite system, a stationary state of a Gaussian translationally invariant Liouvillian is described by a density matrix , where is a hermitian unit trace matrix, which describes the momentum mode pair . The topologically relevant information is encoded in the block of the density matrix in the subspace with even occupation of the modes , (see appendix). The matrix is proportional to , where is the vector of Pauli matrices and is a real three-component vector . The pure states correspond to , i.e. for all . We can naturally extend the definition of to negative following the change of resulting from the transformation of the basis vectors in the subspace under : and . Note that the vector is continuous at because and have only -component.

Once the vector is nonzero for all , the normalized vector defines a mapping of a circle (the Brillouin zone with identified end points due to usual periodicity in the reciprocal lattice) into a unit sphere of end points of . This mapping, however, is topologically trivial (the corresponding homotopy group ), since a circle can always be continuously deformed into a point on the sphere. We therefore need an additional constraint on in order to introduce a nontrivial topology. In our setting, motivated by Kitaev’s model Hamiltonian Kitaev00 (), the constraint is provided by the chiral symmetry Altland97 (); Ryu10 (). In terms of the density matrix, the chiral symmetry is equivalent to the existence of a -independent unitary matrix with , which anticommutes with the traceless part of the density matrix ( in our case): . After representing the matrix in the form , where is a constant unit vector, the chiral symmetry condition reads , i.e., the vector is orthogonal to for all . The end point of is now pinned to a great circle on the sphere such that the vector defines a mapping from the Brillouin zone into a circle. The corresponding homotopy group is now nontrivial, , and such mappings are divided into different topological classes distinguished by an integer topological invariant (winding number). The explicit form of this invariant reads (see appendix)

(10) |

Geometrically, counts the number of times the unit vector winds around the origin when goes across the Brillouin zone, cf. Fig. 5. Importantly, the winding number distinguishes topologically different density matrices for translationally invariant Gaussian systems with chiral symmetry without restriction on the purity of the state.

Let us now discuss specific examples of the steady states resulting from quasi-local dissipative dynamics. As shown in the appendix, for momentum space Lindblad operators , with (unnormalized) Bogoliubov functions which do not induce a current in the steady state (), the solution for the vector reads

(11) |

with . For quasi-canonical deformations (cf. Eq. (6)), the additional symmetry property and simplifies the solution to , implying the purity of the steady state, for all . The solution corresponds to the ground state of some Hamiltonian. For with an integer , the system has a damping gap closing point, and the vector has only -component ( for even and for odd ) for all , leading to . For other values of , one has . However, as one cannot define the direction of for , a global definition of for all is not possible, with the consequence that the potential change of the sign of when passes is meaningless. Note that the gap closing points do not correspond to a phase transition here, because one has either a completely empty or completely filled lattice, such that no thermodynamic observables of a phase transition could sensibly be defined.

For non-canonical deformations, Eq. (7), visualized in Fig. 5, the steady state density matrix is mixed, (cf. also right panel in Fig. 5). Importantly, the topological order persists for quite strongly mixed states, and we find again for . The difference to the previous example is that at , not only the direction of but also the topological invariant is not defined: , alined in the -direction for all , has zeroes: , meaning physically that these modes are in a completely mixed state. The “loss” of topology at can be viewed as a non-equilibrium topological phase transition Rudner09 (); Lindner11 (); Kitagawa10 () as a result of changing the Liouville parameters: The system has well defined thermodynamic properties, since it is half filled, for all . The closing of the dissipative gap at leads to critical behavior, which manifests itself via diverging time scales, resulting e.g in an algebraic approach to steady state (as opposed to exponential behavior away from criticality) Diehl08 (); Verstraete09 (); Diehl10a (); Eisert11 (). Importantly, the vector in the steady state has the reflection property for all , where and . Therefore, the symmetry pattern of the steady state is identical on both sides of the transition, ruling out a conventional Landau-Ginzburg type transition and underpinning the topological nature of the transition.

We see that in our system the stationary state of Liouvillian evolution in an infinite system is indeed characterized (for ) by a nontrivial topological order. If the system is finite, the edges separate the topologically non-trivial bulk from a vacuum, which is topologically trivial. Similar to the Hamiltonian case of topological insulators or superconductors HasanKane (), the “jump” in the topological order guarantees the existence of edge states. For a pure stationary state (quasi-canonical case) this follows from the formal analogy with the Hamiltonian case in which it is well-established (see e.g. Ref. gurarie11 ()). More generally, for a mixed state (non-canonical case), the arguments follow the line of Ref. Kitaev00 () using an alternative equivalent form of Eq. (10) for the topological invariant discussed in the appendix, cf. Eq. (B). The topological origin of the edge modes explains their stability demonstrated above.

## V Implementation in Cold Atomic Gases

Here we devise an interacting (quartic in the fermion operators) Liouville operator which at late times reduces to the quadratic Majorana Liouville operator Eq. (2). The late time limit and the steady state here play the role of a low energy limit and the ground state in equilibrium, respectively, where the physics of weakly correlated superconductors is universally described in terms of Bogoliubov quasiparticle excitations for a variety of microscopic models. Our prescription may thus be seen as a microscopic ”parent Liouvillian”, providing one possible microscopic realization of a Liouville operator which generates the desired properties at and close to steady state.

The implementation idea follows closely an earlier proposal for bosons Diehl08 () and is based on a conspiracy of laser driving and engineered dissipation made possible by immersion of the fermionic target system into a superfluid BEC reservoir. The interaction with the bosonic dissipative reservoir is microscopically based on a conventional s-wave fermion-boson density-density interaction and stands in marked contrast to the proximity effect to a BCS superconductor exploited in solid state implementation proposals of the (Hamiltonian) Majorana wire. Our scheme is illustrated in Fig. 2 and explained in more detail in the appendix. It yields the following number conserving Liouville dynamics,

(12) | |||||

These Lindblad operators give rise to dissipative pairing in the absence of any conservative forces, a mechanism based on an interplay of phase locking and Pauli blocking established recently Diehl10b (). The relation to the Majorana operators is apparent in the thermodynamic limit, where it can be shown that the following general relation between fixed number () and fixed phase () Lindblad operators holds (see appendix),

(13) |

Here are creation (annihilation) parts, respectively, which for Eq. (12) read . In consequence, we see that indeed precisely Kitaev’s quasiparticle operators are obtained as Lindblad operators, . The description in terms of fixed phase operators becomes appropriate at late times, i.e., close to the steady state, where an ordering principle is provided by the macroscopic occupation of only a few correlation functions. The explicit calculation shows that Eq. (2) is produced with effective dissipative rate (see appendix). There we also discuss the leading imperfections, showing that they preserve the chiral symmetry necessary to remain in the above described topological class.

## Vi Conclusions

In this work, we established a complete list of topological features familiar from ground state physics of certain Hamiltonians in engineered dissipative dynamics, highlighting the universality of the concept of topological order in quantum mechanical many-body systems. While the present work has focused on the conceptually simplest system of a quantum wire, the idea of dissipatively induced topological order is more general and will be present in other physical systems and in dimensions higher than one.

## Acknowledgments

We thank C. Bardyn, V. Gurarie, A. Imamoglu, C. Kraus and M. Troyer for helpful discussions. We acknowledge support by the Austrian Science Fund (FOQUS), the European Commission (AQUTE, NAMEQUAM), the Institut für Quanteninformation GmbH, and by a grant from the US Army Research Office with funding from the DARPA OLE program.

## Appendix A Appendix

### a.1 Implementation in cold atomic gases

Microscopic Model – Here we specify a physical setup leading to Eq. (12), following an earlier proposal for bosons Diehl08 () and illustrated in Fig. 2. We start from an optical superlattice setting, with lower sites which make up the fermion wire and correspond to the lowest Bloch band of the lattice, and auxiliary sites associated to the second Bloch band located on the links. These two Bloch bands are coupled via driving lasers, whose Rabi frequencies are chosen with opposite sign for each pair of lower sites. This amounts to a commensurability condition of driving and lattice laser, and ensures the relative minus sign in the annihilation part of the Lindblad operators of Eq. (12), . By immersing the whole setting into a BEC reservoir, particle superpositions in the upper band can spontaneously decay back to the lower one by emission of a Bogoliubov phonon into the BEC bath. This process is isotropic and short-ranged for suitable bath parameters, giving rise to the two-site creation part in Eq. (12) with relative plus sign, . Microscopically, this dissipative mechanism requires a standard s-wave fermion-boson interaction between system (fermions in the optical superlattice) and bath (BEC) particles Diehl08 (). An effective single particle microscopic model such as Eq. (12) is obtained by integrating out the upper band under suitable detuning conditions for the driving laser. Using the recently developed tools of single site addressability for optical lattices Greiner10 (); Kuhr10 (), an edge can be constructed by cutting the superlattice in the right place cf. Fig. 2.

Eq. (12) describes an interacting, number conserving Liouville dynamics, which leads to dissipative pairing Diehl10b (), described by a pure BCS-type wavefunction for paired fermions. There are two ways of seeing this, either by explicit construction of the dark state wavefunction (fixed particle number) or by deriving a suitable mean field theory (fixed phase). In the thermodynamic limit both procedures are equivalent. Here we follow the second option, and highlight the more precise connection between the two approaches below.

Mean field theory at late times – At late times, following Ref. Diehl10b () we can make use of the proximity to the steady state to derive a quadratic mean field theory for this dissipative dynamics. In the spirit of BCS theory, the ordering principle behind this approximation is the macroscopic occupation of only a few correlation functions in steady state, which can be evaluated explicitly on the exactly known steady state. To implement the approximation, as usual in BCS type mean field theories, we give up exact particle number conservation and work with a fixed phase.

We thus start from number conserving Lindblad operators , where the creation part and the annihilation part with translation invariant complex position space functions ; in the example Eq. (3) discussed in the text, . For the practical calculation, we switch to momentum space, where we have . The Fourier transforms for the creation and annihilation part are local in momentum space, , with in the example (due to the structure of the Liouville operator, the exponential prefactors are irrelevant and can be omitted). In the following we work with functions with the property , which however do not necessarily obey a normalization constraint.

For a setting with fixed phase, we now make an ansatz for the steady state density matrix of the product form , with density matrix for each mode pair . Inserting this ansatz into Eq. (12) in momentum space, and using the projection prescription on the mode pair , we obtain the equations of motion for the single pair density matrices in the presence of nonzero mean fields. These result from the coupling to other momentum modes, with values determined by the steady state properties in the late time dynamics. The resulting mean field dynamics close to the steady state is given by

(14) | |||||

where are fermionic quasiparticle operators up to a normalization, obeying the anti-commutation relations . In our example, and , i.e. the dynamics has constant damping rate for all modes. In fact, transforming back to position space produces precisely Eq. (3), .

Since , we see that the above ansatz indeed provides the correct steady state solution. In addition, the equivalence of fixed number and fixed phase wavefunctions can now be justified in our nonequilibrium context a posteriori in the thermodynamic limit: the fixed phase BCS state has relative number fluctuations , where is the particle number operator and the number of degrees of freedom.

We emphasize that, remarkably, the late time evolution of the master equation (12) naturally gives rise to quasilocal squeezing of fermions. The mechanism behind this effect is in close analogy to a superconductor: the system here acts as its own reservoir by providing an order parameter, which allows for the appearance of the off-diagonal pair annihilation / creation terms h.c. in Eq. (14).

General relation of fixed number and fixed phase Lindblad operators – Here we show Eq. (13). Suppose we are given a fixed phase Liouvillian defined by a set of Lindblad operators of a form

(15) |

using the above conventions, such that the momentum space Lindblad operators read . The dark state (with property for all or ) is conveniently formulated in momentum space and reads

(16) |

Now we want to show that the number conserving version is given by the following expressions for the Lindblad operators and related dark state wavefunction,

(17) |

To prove it, we proceed in momentum space. For the normal ordered Lindblad operators , the dark state property for all or is equivalent to the vanishing of the following commutator for all ,

This is true if and only if for all , i.e. for the wavefunction as claimed above.

Note that working with the functions for finite requires an infrared momentum cutoff such that ( the density of particles), that is consistent with the thermodynamic limit .

Thus, for a given real space quadratic master equation, we can immediately construct a number conserving version and vice versa by the above arguments, and indicate the respective fixed number of fixed phase exact dark state wavefunctions.

Imperfections – In the above implementation scheme, the annihilation part is well under control since it relies on a locking of lattice and driving laser. The dominant imperfection appears in the creation part: In the case that the phonon wavelength in the BEC bath is not smaller than the lattice spacing (subradiant case), the decay may take place over several lattice sites, such that e.g. , giving rise to fixed phase Lindblad operators with in momentum space, which exhibit fermionic anticommutation relations guaranteeing a pure steady state. The resulting vector still lies in a plane for arbitrary , i.e. the imperfection preserves the chiral symmetry.

### a.2 Properties of the finite size system

Equations of motion in the Majorana basis – For the practical treatment of the finite size system quadratic in the fermion operators and with Gaussian initial conditions, it is convenient to encode the information in the covariance matrix of second moments in the real Majorana basis. For convenience, we repeat and extend here the definitions of the main text. The Majorana operators obey , , , , where labeling the physical sites in a one dimensional lattice. Following Eisert11 (), we start with a quadratic master equation written in the Majorana basis, i.e. with for hermitian matrix and , for component column vectors . The Liouvillian parameters are then encoded in a hermitian matrix . The equation of motion for the covariance matrix with entries reads ()

(19) |

with real matrices , . The spectrum of is positive semidefinite Prosen08 () and a pure state is signaled by the eigenvalues of being all equal to .

Zero modes – The real symmetric matrix is diagonalizable and has real eigenvalues and -vectors. First we show the existence of zero eigenvalues for this matrix for a wide class of finite systems via explicit construction of the corresponding eigenvectors. The ideal Lindblad operators for the dissipative topological 1D quantum wire are characterized by the following vector in the Majorana basis:

(20) |

with and . Let us now consider the more general vector specified with

(21) |

For each , there are two orthogonal vectors given by

(22) |

and by construction, they fulfill . Due to the structure of the matrix as a sum of outer products of the vectors, we find precisely two zero modes for any finite system size given by the component vectors

(23) |

where represent the left / right Majorana modes. The reason behind these two vectors being zero vectors of the matrices and is that they are orthogonal to the linearly independent vectors and in a space of dimension . Hence, , , and form a complete basis in this space.

Let us now focus on the deformations studied in the main text. Here, the Majorana modes are characterized by vectors

(i) Canonical deformation:

(24) |

(ii) Non-canonical deformation:

(25) |

with . The associated localization length is then given by , with the lattice constant.

In the case of a deformation being non-homogenous but with value changing randomly from site to site, the zero modes are still present, but the expression for and is given by

(26) |

Dissipative isolation of the subspace – As argued in the main text, in addition to the existence of a zero mode subspace, the dissipative evolution must ensure its isolation from the bulk. The conditions for such a situation are readily formulated generally, without reference to the quadratic setting. For a set of Lindblad operators making up the total Liouvillian , we may introduce projectors on the edge (zero mode) and bulk subspaces, and , respectively. A decoupled edge subspace appears if the Lindblad operators are block diagonal, , with the edge block identical to zero, . We then obtain a dissipative evolution for the density matrix in this projection,

The bulk dissipative evolution has Lindblad form. The density matrix in the subspace is a constant of motion. The coupling density matrix elements h.c. damp out according to , i.e. exponentially fast in the presence of a dissipative gap.

This situation is indeed present in our quadratic setting. We switch to the spectral representation , where are the nonvanishing eigenvalues and the corresponding eigenvectors. In addition, we have zero modes