Braiding errors in interacting Majorana quantum wires

TL;DR

This paper investigates braiding errors in Majorana bound states within Kitaev wires, highlighting how non-adiabaticity, hybridization, and interactions influence the fidelity and phase accuracy of braiding operations.

Contribution

It provides a detailed numerical analysis of braiding errors considering the full quasi-particle background and interaction effects, advancing understanding of topological quantum computation.

Findings

Fidelity loss mainly due to incomplete adiabaticity and MBS hybridization.

Explicit braiding phase extraction reveals breakdown of independent-particle picture.

Nearest-neighbor interactions can improve or impair braiding depending on their sign and strength.

Abstract

Avenues of Majorana bound states (MBSs) have become one of the primary directions towards a possible realization of topological quantum computation. For a Y-junction of Kitaev quantum wires, we numerically investigate the braiding of MBSs while considering the full quasi-particle background. The two central sources of braiding errors are found to be the fidelity loss due to the incomplete adiabaticity of the braiding operation as well as the hybridization of the MBS. The explicit extraction of the braiding phase in the low-energy Majorana sector from the full many-particle Hilbert space allows us to analyze the breakdown of the independent-particle picture of Majorana braiding. Furthermore, we find nearest-neighbor interactions to significantly affect the braiding performance to the better or worse, depending on the sign and magnitude of the coupling.

Figures

Equations

$$H$$

$$H_{\mathrm{Y}}^{}={\mu }_0^{}{c}_0^{†}{c}_0^{}$$

$${\mu }_j^{(n)}(\tau )=2{\mu }_{c}m\left(\frac{\tau }{T}[1+\alpha (\ell -1)]-\alpha (\ell -j)\right),$$

$$m(q)$$

$$r(q)$$

$$|\psi (t+T)⟩$$

$$w_{\mathrm{loss}}(t_{\mathrm{f}})=1-F(t_{\mathrm{f}}{)}^2=1-|⟨\psi (t_{\mathrm{f}})|\psi (t_0)⟩{|}^2,$$

$${\varphi }_{\text{g}}[C_0]=\arg ⟨\psi (t_0)|\psi (t_{\text{f}})⟩-\Im \underset{t_0}{\overset{t_{\text{f}}}{\int }}⟨\psi (t)|\dot{\psi (t)}⟩\mathrm{d}t,$$

$${\varphi }_{\mathrm{g}}[{\tilde{C}}_0]$$

$$\Phi (t)={\varphi }_{\mathrm{g}}^o(t)-{\varphi }_{\mathrm{g}}^e(t).$$

Full text

Introduction.

Topological quantum states of matter have become one of the most vibrant fields of contemporary condensed matter physics. One subbranch thereof, the search for realizations of topological quantum computation, has decisively fueled the interest both from theory and experiment to address manifold scenarios of this kind. Among them, Majorana bound states (MBS) have been attracting pivotal attention in recent years Alicea (2012); Leijnse and Flensberg (2012); Beenakker (2013); Sarma et al. (2015). While Majorana-type quasiparticles have been previously addressed in the context of the fractional quantum Hall effect Moore and Read (1991); Read and Green (2000), unconventional superconductivity Read and Green (2000); Kitaev (2001); Ivanov (2001), and spin liquids Kitaev (2006); Greiter and Thomale (2009), the proposal by Fu and Kane Fu and Kane (2008) to employ conventional superconducting proximity effect to stabilize MBS in vortex cores at the surface of a topological insulator unleashed a remarkable body of research that has brought the detection and manipulation of MBS much closer to reality. After crucial progress towards simpler realizations of proximity-induced topological superconducting phases Lutchyn et al. (2010); Oreg et al. (2010); Nadj-Perge et al. (2013), several experimental groups have reported increasingly compelling evidence for the observation of Majorana zero modes in semiconductor nanowires and iron atomic chains that are coupled to a bulk superconductor Mourik et al. (2012); Rokhinson et al. (2012); Das et al. (2012); Deng et al. (2012); Churchill et al. (2013); Nadj-Perge et al. (2014); Albrecht et al. (2016); Deng et al. (2016); Pawlak et al. (2016); Feldman et al. (2017); Zhang et al (2016). The unambiguous detection of MBS, however, has so far remained elusive, and parity-changing quasiparticle poisoning through single-particle tunneling into the wire or chain likewise constitutes a major challenge. One crucial experimental finding in favor of the existence of Majorana zero modes would be a braiding experiment Alicea et al. (2011); van Heck et al. (2012); Li et al. (2016); Aasen et al. (2016), revealing their non-trivial braiding statistics, and showing the path to more complex multi-MBS protocols.

A plethora of approaches is currently underway towards the realization of topological quantum computation in networks of semiconductor nanowires proximity-coupled to a superconductor Hyart et al. (2013); Aasen et al. (2016); Landau et al. (2016); Plugge et al. (2016, 2017); Vijay and Fu (2016a); Karzig et al. (2016). Alternative routes include the aforementioned proximitized topological-insulator surfaces Fu and Kane (2008); Vijay et al. (2015); Vijay and Fu (2016b), and atomic chains Nadj-Perge et al. (2013); Li et al. (2016) or cold-atom setups Kraus et al. (2013); Laflamme et al. (2014) where the originally envisioned lattice realizations of topological superconductivity Kitaev (2001) are more directly accessible. Given the amount of different obstacles still to overcome, it is hard to predict which direction will succeed. At the current stage of the field, it thus appears worthwhile following up on several of these directions at the same time: New results on any given approach will offer guidance for inevitable refinements of the others.

In this paper, we numerically analyze the braiding process of MBSs in a minimal model of an interacting wire network. We define a Y-junction of three Kitaev chains Kitaev (2001) in which we investigate the controlled time-evolution of a single pair of MBSs as the most elementary braiding operation. We employ time-dependent on-site potentials to locally drive the system in and out of the topologically non-trivial phase, which hence allows us to spatially move the MBSs Alicea et al. (2011). The potential-manipulation protocol is analyzed to yield an optimal performance as to minimize non-adiabatic effects. The small system sizes accessible through exact diagonalization allow us to systematically identify the different sources of braiding errors. Since we keep the full Hilbert space during the braiding operation, we can correlate the taken protocol with the fidelity of the MBS sector in which we operate. In addition, this allows us to extract explicitly the geometric exchange “braiding” phase acquired during the execution of the braid, directly encoding the non-Abelian statistics of MBSs. By varying the superconducting wire bulk gap, we then find how the localization length of the individual Majorana zero modes as well as the resulting MBSs hybridization affect the braiding operation. Finally, interactions in the wires, modeled by a nearest-neighbor density-density coupling Gangadharaiah et al. (2011); Stoudenmire et al. (2011); Sela et al. (2011), are found to predominantly affect the braiding operation through its impact on the wire bulk gap. Depending on the non-interacting gap versus hopping strength, for weak attractive (repulsive) interaction, the braiding operation is more stable due to an increase (decrease) of the effective superconducting gap and reduction of the MBS localization length. For strong attractive or repulsive nearest-neighbor coupling, the interactions have a negative impact on the braiding performance, where ultimately the topological phase becomes inaccessible altogether Gangadharaiah et al. (2011); Stoudenmire et al. (2011); Sela et al. (2011).

Kitaev chain.

The (interacting) Kitaev model is the most elementary system exhibiting MBSs Kitaev (2001); Gangadharaiah et al. (2011). It is given by the Hamiltonian Hassler and Schuricht (2012); Thomale et al. (2013)

[TABLE]

where $c^{\dagger}_{i}$ ($c^{\phantom{\dagger}}_{i}$) is a creation (annihilation) operator of the spinless fermion on site $i$ and $n_{i}=c^{\dagger}_{i}c^{\phantom{\dagger}}_{i}$ is the density operator. $J$ denotes the nearest neighbor hopping amplitude, $\Delta=|\Delta|{\mathrm{e}}^{i\varphi}_{\phantom{0}}$ the $p$-wave superconducting pairing amplitude with phase $\varphi$, $V$ the strength of the nearest neighbor interaction, and $\mu_{i}$ the site-dependent potential. The superconducting pairing term breaks $U(1)$ fermion number symmetry down to $\mathbb{Z}_{2}$ (parity conservation), where the two parity sectors will be denoted even $(e)$ and odd $(o)$ below. Unless stated otherwise, we initially constrain ourselves to the non-interacting limit $V=0$. The Kitaev chain features a topological trivial and non-trivial phase for $|\mu|>2J$ and $|\mu|\leqslant 2J$, respectively, so long as $|\Delta|>0$ Kitaev (2001); Alicea (2012). For the ideal parameter set $\mu=0$, $|\Delta|=J$ and hard wall boundary conditions, hence residing in the non-trivial regime, Eq. (1) takes a particularly simple form in terms of a Majorana fermion representation $H=-J\sum_{j=1}^{L-1}i\gamma^{\phantom{\dagger}}_{2j}\gamma^{\phantom{\dagger}}_{2j+1}$, where $\gamma^{\phantom{\dagger}}_{2j-1}={\mathrm{e}}^{i\varphi/2}_{\phantom{0}}c^{\phantom{\dagger}}_{j}+{\mathrm{e}}^{-i\varphi/2}_{\phantom{0}}c^{\dagger}_{j}$ and $\gamma^{\phantom{\dagger}}_{2j}=({\mathrm{e}}^{i\varphi/2}_{\phantom{0}}c^{\phantom{\dagger}}_{j}-{\mathrm{e}}^{-i\varphi/2}_{\phantom{0}}c^{\dagger}_{j})/i$ are real Majorana operators $\gamma^{\dagger}_{i}=\gamma^{\phantom{\dagger}}_{i}$ obeying a Clifford algebra, $\{\gamma^{\phantom{\dagger}}_{i},\gamma^{\phantom{\dagger}}_{j}\}=2\delta^{\phantom{\dagger}}_{ij}$. We observe that the first and last Majorana modes, residing on the first and last site, respectively, are decoupled from the rest of the fermionic chain. The corresponding Majorana operators, $\gamma^{\phantom{\dagger}}_{1}\equiv\gamma^{\phantom{\dagger}}_{\mathrm{L}}$ and $\gamma^{\phantom{\dagger}}_{2L}\equiv\gamma^{\phantom{\dagger}}_{\mathrm{R}}$, do not appear in the Hamiltonian, $[H,\gamma^{\phantom{\dagger}}_{\mathrm{L}}]=[H,\gamma^{\phantom{\dagger}}_{\mathrm{R}}]=0$, and form a fermionic zero-energy mode comprised by two individual Majorana (zero) modes localized on single sites at each end of the chain. After introducing a non-local fermionic operator $f=\frac{1}{2}(\gamma_{\mathrm{L}}+i\gamma_{\mathrm{R}})$, these zero-energy states can be identified with MBS eigenstates $f|{0}\rangle=0$ and $f_{\phantom{0}}^{\dagger}|{0}\rangle=|{1}\rangle$, which reside in the middle of the effective bulk gap $4|\Delta|/J$. The states $|{0}\rangle$ and $|{1}\rangle$ also coincide with the parity eigenstates of the full fermionic chain, and can be associated to the ground state wavefunctions of the Kitaev Hamiltonian (1) in the even ($e$) and odd ($o$) parity sector, respectively Greiter et al. (2014).

Away from this particular parameter set but still within the non-trivial regime, i.e. for $|\Delta|<J$ or $0<\mu<2J$, the Majorana zero modes develop an exponential tail extending into the bulk, yielding a small energy split of the two eigenstates previously located at zero energy (in different parity sectors). This splitting energy, usually referred to as MBS hybridization energy, decays exponentially with the distance between MBSs on a length-scale governed by $\Delta$. The protection offered by the bulk energy gap and spatial isolation of MBSs are the main stimuli for $|{0}\rangle$ and $|{1}\rangle$ to serve as a basis for the sought-after topological qubit. As we show below, targeting the entire Hilbert space – and not only a subspace of sub-gap Majorana states in terms of an effective low-energy theory as conventionally assumed Alicea et al. (2011); Aasen et al. (2016); Cheng et al. (2011); Clarke et al. (2011); Halperin et al. (2012); Hyart et al. (2013); Karzig et al. (2013); Scheurer and Shnirman (2013); Karzig et al. (2015a, b); Kraus et al. (2013); Laflamme et al. (2014); Landau et al. (2016); Li et al. (2016); Plugge et al. (2016, 2017); Karzig et al. (2016); Vijay and Fu (2016a); van Heck et al. (2012); Vijay et al. (2015); Vijay and Fu (2016b) – is essential when manipulations such as braiding are performed on finite time scales (non-adiabatically) and in finite-size systems. Preceding works undertook steps in this direction Amorim et al. (2015); Pedrocchi and DiVincenzo (2015); Knapp et al. (2016), but attempted to ground the errors in an effective low-energy (Majorana) picture.

In the following, via a sequence of predefined exact manipulations of potenitals $\mu_{i}$ of Hamiltonian (1), we perform a braiding operation imposed on an initially prepared set of MBSs Alicea et al. (2011). Such manipulations employ unitary but in general non-adiabatic dynamics. This implementation requires a way to (i) dislocate MBSs without destroying them or creating new ones, (ii) realize a system geometry that allows to exchange Majorana modes while keeping their overlap exponentially small during the entire braiding process, and (iii) provide an exchange process that gives fidelity, taken as the overlap between initial and final states, close to unity.

At the elementary level of our description, the goal is not to define the optimal exchange protocol applied to a realistic setup for topological quantum computation (though we made sure to select a favorable implementation among several “simple” test routines). Rather, for the simplest braiding operation possible, we intend to identify the time scale and conditions under which it is feasible to perform a braiding with the desired fidelity and braid-phase accuracy. For an analysis of tailor-made braid protocols and possible drawbacks or limitations thereof, based on effective models of the low-energy (Majorana) sector of the system, cf. Refs. Karzig et al. (2015a, b); Knapp et al. (2016).

Y-junction.

The minimal geometry required for the exchange “braiding” of MBSs is a junction with at least three legs (tri-junction) Alicea et al. (2011); Alicea (2012). It can be formed by connecting three Kitaev chains (later referred to as junction legs) to an additional connecting junction site at $i=0$ 111The site indices in each junction leg are taken as increasing from the central tri-junction site to the free edge., which otherwise witnesses parameters sets identical to the respective adjacent leg:

[TABLE]

Here $H^{(n)}$ for each leg corresponds to the Kitaev Hamiltonian (1) with parameters $J^{(n)}$, $\Delta^{(n)}$, $V^{(n)}$, and $\mu^{(n)}_{i}$. In general, one may choose different $p$-wave pairing terms and hopping amplitudes in each leg of the Y-junction. For simplicity, we here assume that $J^{(n)}=1$ (thus setting the reference scale in plots and equations below), $\Delta^{(1)}=-\Delta$, $\Delta^{(2)}=\Delta\,e^{i\varphi_{2}}$, and $\Delta^{(3)}=\Delta\,{\mathrm{e}}^{i\varphi_{3}}$, with a real pairing amplitude $\Delta>0$. The superconducting pairing phase $\varphi_{n}$ is chosen to coincide with the corresponding geometric angle of the $n$-th leg of the Y-junction, cf. Fig. 1. This type of triplet pairing could, for example, be realized in a Y-junction placed in proximity to a chiral $p$-wave superconductor Halperin et al. (2012). The connection of two Kitaev chain legs, both of which are in the topological (trivial) phase, produces a single segment residing in the topological (trivial) phase; unless the mutual phase difference between the pairing terms of these chains equals exactly $\varphi_{i}-\varphi_{j}=\pi$, MBSs on the connecting site fuse, and the entire segment ends up with the two left-over MBS residing on the free ends of the new, larger topological segment Alicea et al. (2011). If the mutual phase, however, is exactly $\pi$, MBSs on the connecting site do not fuse and the joint segment maintains all four MBSs (one from each end of the initial chains) – a case we intentionally avoid in the following.

By ramping up the on-site potentials over the critical value of $\mu_{c}=2J$ starting from one end of the topological segment (i.e., at the location of a Majorana zero mode), locally one drives the system into the trivial phase. As the ramp of potentials is successively carried out along the chain, one thus continuously displaces the MBS alongside with the trivial-topological domain wall. This action can be undone by again lowering the local potentials below the critical value, and switching the chosen sites back from the trivial to topological phase. Finally, concerting such ramping procedures in sequences across all three junction legs, one can create the desired Majorana zero mode exchange or braiding operation Alicea et al. (2011). This process is schematically depicted in Fig. 1. We will only study Y-junctions with legs of equal length $L^{(n)}=\ell$, and take $\varphi_{2}=-\varphi_{3}=\varphi$. The initially homogeneous potentials in the junction legs are $\mu^{(1)}=\mu^{(2)}=0$ and $\mu^{(3)}=2\mu_{c}$, meaning that legs $1$ and $2$ reside in the topological, and leg $3$ in the trivial phase. In passing we note that since the time evolution of an initial state prepared like this fixes a time direction, there is a subtle time-reversal symmetry breaking appearing in the braiding setup as seen by reversing the sign of $\varphi$ in the time evolution. We thus only consider $\varphi>0$ in the following.

Ramping protocol.

We investigated several ramping protocols for the successive time variation of the on-site potentials $\mu_{i}^{(n)}(t)$. The most stable results were found by employing a sine-squared ramp defined as

[TABLE]

where $\tau\in[0,T]$ denotes time during the respective ramp-up protocol step in Fig. 1, and we used the scalar function

[TABLE]

Here, $r(q)$ is a linear ramp of unit height and duration

[TABLE]

We also show results obtained by the simpler, linear “guillotine” ramp, $m(q)=r(q)$, where the smoothening by the sine-squared function in Eq. (4) is switched off. Both procedures raise the on-site potentials in the $n$th leg, containing $\ell$ sites, from $\mu^{(n)}_{\text{init}}=0$ to $\mu^{(n)}_{\text{f}}=2\mu_{c}=4J$ within the time period $T$. The time required to lift the on-site potential on each individual site is $T/(1+\alpha(\ell-1))$, and the ramping is delayed by $\alpha T/(1+\alpha(\ell-1))$ between consecutive sites. For a guillotine ramp, the modulations of local potentials are reminiscent of a guillotine knife passing a rectangular window (hence the name), and $\alpha$ is the inclination of the knife.

The delay coefficient (inclination) can be varied between $\alpha=0$ and $\alpha=1$, corresponding to a simultaneous ramp of all sites or a consecutive site-by-site ramping, respectively. In the following, we set $\alpha=0.025$ and vary the time step $T$, where $\alpha^{-1},T\gg\ell$ at all times, i.e., a small ramp delay (inclination) and slow ramp protocol. Finally, the ramp-down is implemented as an exact time-reverse of the ramp-up protocol, with $\tau\to T-\tau$ in Eq. (3). The ramping (and hence braiding) protocol is accomplished numerically by considering the piece-wise constant, time-independent Hamiltonian

[TABLE]

where $\mathcal{T}$ is the time ordering operator, and $H^{\phantom{\dagger}}_{\mathrm{Y}}(t)$ corresponds to (2) with time-dependent local potential $\mu^{(n)}_{j}\rightarrow\mu^{(n)}_{j}(t)$ given in (3). The time-step discretization $\Delta t/T\ll 1$ governs the accuracy of approximation of the continuous time-dependent Hamiltonian by a piece-wise constant one, and is chosen sufficiently small against all energy scales of the problem, yielding a total of $N=T/\Delta t$ steps per braiding sweep in Fig. 1. As such, the outlined braiding procedure demands a total (final) braid time $t_{\text{f}}=6T$, and performs a complete cyclic evolution in the parameter space of the Hamiltonian.

As braiding errors are detected, there are several parameters that can be used to control and alter the braiding process. The ramping speed can be controlled by the time period $T$ required for each step of the protocol, and by the ramp delay/inclination $\alpha$. Since the latter also determines the steepness of the trivial-topological domain wall propagating through the chain, generally, for steeper inclinations, a slower braid velocity $\sim\ell/T$ with larger time period $T$ is required Karzig et al. (2013); Scheurer and Shnirman (2013). Furthermore, the rate $\Delta t/T$ above determines the smoothness of the ramping potential. In the following, we take $\Delta t/T\ll 1$ small against all relevant energy scales and simulate an effectively smooth ramp.

Fidelity and geometric braiding phase.

In order to achieve an ideal braiding process, it is necessary (but not sufficient) that the fidelity between the initial and obtained final state equals unity. Since numerical simulations of the full Hilbert space are costy in general, and sampling the fidelity for many input states – in particular if one is interested in longer braid sequences – is not feasible, we propose to analyze the topological character and protection of braid operations for the ground states of the two parity sectors. For this purpose, the braiding statistics are directly encoded in the relative geometric phase factor (the braid phase) between final and initial states in the even and odd parity sector. Braiding errors then manifest both in deviations of the braid fidelity from unity and of the braiding phase from $\pi/2$. Both quantities are accessible by considering input states $|{0}\rangle$ and $|{1}\rangle$ of the topological qubit, chosen as the ground states of $H_{\mathrm{Y}}(t_{0}=0)$ (2) in the respective total-parity sector. (Note that transitions between even $(e)$ and odd $(o)$ total-parity states are forbidden under unitary time evolution (9) with a parity-conserving Hamiltonian (2).)

In order to quantify the fidelity loss during the braiding, we introduce the loss function

[TABLE]

where $F(t_{\mathrm{f}})=|\langle{\psi(t_{\mathrm{f}})}|{\psi(t_{0})}\rangle|$ is the fidelity for pure states Nielsen and Chuang (2000), taken between the initial and final, time-evolved state. Note that zero loss (unit fidelity) does not necessarily imply that the overall process is adiabatic, however, a cyclic and completely adiabatic process should always yield $w_{\mathrm{loss}}=0$. Similar measures, but with an adiabatically evolved state as reference, have been used in Refs. Karzig et al. (2015a); Amorim et al. (2015). In the simple two-Majorana braiding setup of Fig. 1, however, superpositions of states of different total parity are unphysical, and fidelities in a fixed parity sector lack information about the braid phase. Further, in the full Hilbert space of the interacting Kitaev chain (1), obtaining general analytic results, even under adiabatic evolution, appears impossible.

Our braiding protocol describes a cyclic evolution in the space of Hamiltonian parameters, where the Hamiltonian takes its original form at the end of the process, i.e., at time $t_{\text{f}}=6T$. Assuming adiabaticity, one can measure a Berry phase Berry (1984) acquired during this cyclic evolution of the system, where the time-evolving state should always correspond to an eigenstate of the instantaneous Hamiltonian. The fidelity after execution of the braid protocol, however, is not exactly unity (we find $w_{\mathrm{loss}}>0$) because of non-adiabatic spectral losses due to the finite ramping time. This in fact will manifest in a non-constant dynamical evolution of the acquired phase even after the system is evolved back to the initial Hamiltonian, cf. Fig. 2. As a consequence, the Berry phase is not applicable to realistic, finite-time braiding protocols. The unitary time-evolution of states under a continuous time-dependent Hamiltonian, however, still describes a smooth curve in the Hilbert space. Since quantum states that differ merely by a phase factor give rise to the same physics, one may then employ the gauge- and parametrization-invariant functional Samuel and Bhandari (1988); Mukunda and Simon (1993)

[TABLE]

which measures the geometrical phase for the smooth open curve of normalized (unit) vectors in the Hilbert space $\mathcal{C}_{0}=\{|{\psi(t)}\rangle\in\mathcal{N}_{0}\ |\ t\in[t_{0},t_{\text{f}}]\subset\mathbb{R}\}\subset\mathcal{N}_{0}$. Here $\mathcal{N}_{0}$ denotes the space of normalized vectors and $C_{0}$ is the projection of $\mathcal{C}_{0}$ to the so called ray space (or projective Hilbert space), where members that differ only by a phase factor are regarded as equivalent 222 An equivalence class of states in ray space is a projection operator $|\psi\rangle\langle\psi|$ to the equivalence class represented by $|\psi\rangle$. This is a natural projection that maps each vector to the ray on which it lies.. We note that if the fidelity between the initial and final states is exactly one (zero loss), $C_{0}$ becomes a closed curve in ray space, and $\phi_{\text{g}}[C_{0}]$ corresponds to the Aharonov-Anandan geometrical phase Aharonov and Anandan (1987). If the process is also adiabatic, one recovers the conventional Berry phase Berry (1984) (for discussion, cf. Ref. Mukunda and Simon (1993)). In our time-dependent many-particle problem, the ramping (and hence braiding) protocol is accomplished numerically by considering the piece-wise constant, time-independent Hamiltonian in Eq. (9). Therefore, $C_{0}$ is a continuous but generally not a smooth curve. To determine the geometric phase numerically, we hence employ a discretized variant of the above functional Mukunda and Simon (1993)

[TABLE]

which measures the geometric phase for the closed ${N+1}$-sided polygon of unit vectors in the Hilbert space $\tilde{\mathcal{C}}_{0}=\left\{|{\psi(t_{j})}\rangle\in\mathcal{N}_{0}\ |\ t_{j}=t_{0}+j\Delta t,~{}j=0,\dots,N\right\}$, with$N=t_{\mathrm{f}}/\Delta t$. The closed polygon $\tilde{C}_{0}$ parametrizes the braid protocol and corresponds to a projection of time-ordered set $\tilde{\mathcal{C}}_{0}$ onto ray space. This expression is also referred to as the Bargmann vertex formula for the geometric phase. In the continuous limit $\Delta t\rightarrow 0$ in Eq. (9), the expression in Eq. (12) reduces to the continuous-time version in Eq. (11). Note also that the final expression above is cyclically symmetric. In terms of a physical interpretation, the first and second term in Eqs.(11)/(12) correspond to the complete and the local, dynamic parts of the acquired phase, respectively. We thus explicitly exclude purely dynamic contributions, for example those generated by a finite Majorana hybridization Alicea et al. (2011); Clarke et al. (2011); Cheng et al. (2011), from our considerations.

Exact diagonalization.

The initial states $|{e}\rangle$ and $|{o}\rangle$ are chosen as the ground states of $H_{\mathrm{Y}}(t_{0})$ in even and odd total-parity sectors, obtained through Davidson exact diagonalization Sadkane and Sidje (1999). The braiding is then carried out iteratively by the Arnoldi approximation of the product of matrix exponentials (9) with a time-evolving state vector Moler and Loan (2003). Since the time-dependent Hamiltonian conserves fermion parity, the simulation for each total-parity sector can be carried out separately. Further, while strictly speaking the geometric phase is only defined after the braid process is carried out completely and a closed loop in parameter space is accomplished, it is convenient to define the exchange phase at arbitrary time $t$ of the braiding process,

[TABLE]

Here, $\phi^{e}_{\mathrm{g}}(t)$ and $\phi^{o}_{\mathrm{g}}(t)$ are the individual phases acquired during the time evolution of states $|{e}\rangle$ and $|{o}\rangle$. We note ${|{e}\rangle~{}\rightarrow~{}{\mathrm{e}}^{i\phi^{e}_{\mathrm{g}}(t)}|{e}\rangle}$ and ${|{o}\rangle~{}\rightarrow~{}{\mathrm{e}}^{i\phi^{o}_{\mathrm{g}}(t)}|{o}\rangle}$, and in the following drop the sub-script “$\mathrm{g}$” for brevity.

We now are equipped to study the symmetric Y-junction with equal-size legs $\ell=5$, where qualitatively similar results were found for leg sizes $\ell\leqslant 4$. At the ideal Kitaev point $\Delta=J=1$, the lowest-energy (initial) states in the even- and odd-parity sector correspond to $|{0}\rangle$ and $|{1}\rangle$ states of the topological qubit, i.e., they differ in terms of the Majorana sector only. After the exchange “braiding” of Majorana modes, in general, both states will acquire different geometrical phases $\phi^{\mathrm{e}}$ and $\phi^{\mathrm{o}}$.In an ideal case the accumulated many-body geometrical phase for $|{0}\rangle$ and $|{1}\rangle$ (and hence for $|{e}\rangle$ and $|{o}\rangle$) will only differ by a non-trivial phase stemming from the protected braid exchange of Majoranas $\gamma_{\mathrm{L}}$ and $\gamma_{\mathrm{R}}$. In Fig. 2 we show the time evolution of phases $\phi^{\mathrm{e}}(t,\varphi)$ and $\phi^{\mathrm{o}}(t,\varphi)$, as well as the phase difference $\Phi(t,\varphi)$ in Eq. (13), for sine-squared and guillotine ramp protocols and several superconducting pairing phases. Here, after finishing the braiding protocol at time $t_{\mathrm{f}}=6T$, we continue the time evolution for yet another period $T$. In the fully adiabatic case, the geometrical phases in each parity sector – and hence the difference between them – should remain unchanged for times $t>t_{\mathrm{f}}$, i.e., when the braiding exchange has been executed completely. In our case, however, the phases keep evolving approximately linearly in time (insets in Fig. 2). We wish to emphasize that this phase evolution is not due to purely dynamical phases as, e.g., caused by finite Majorana hybridization Clarke et al. (2011); Cheng et al. (2011), which are explicitly excluded in Eq. (12). Instead, we can rationalize the continued non-trivial evolution by observing that the fidelity of the obtained states at time $t=t_{\mathrm{f}}$ is not exactly one, or equivalently the loss as defined in Eq. $\eqref{eq:lost_fidelity}$ is non-zero, cf. Fig. 4 and discussion. Since in the adiabatic (zero-loss) limit the system should have reached a state equivalent to the initial state but with the Majoranas exchanged, this implies that instead, the system no longer resides in an eigenstate of its Hamiltonian for $t\geqslant t_{\mathrm{f}}$. Such non-adiabatic spectral losses hence correspond to a leakage of wave function weight out of the Majorana, and into the excited-state sector. Consequently, the state at $t_{\mathrm{F}}\geqslant t_{\mathrm{f}}$ – when projected to the (initial) qubit Hilbert space $\{|{0}\rangle,|{1}\rangle\}$ – will “dephase” with time. This finding emphasizes the danger of non-adiabatic errors and excitations even in the most basic braiding experiments, and cannot be alleviated by just increasing the system size Pedrocchi and DiVincenzo (2015); Knapp et al. (2016). Since there is only a weak dependence on $\varphi$ as shown in Fig. 2, we will constrain ourselves to the symmetric case $\varphi=\pi/3$ in the following. Note that for the finite system size considered, this corresponds to the biggest time-evolving protective gap among all other geometric configurations of Y-junctions.

Braiding errors.

There are several sources of errors that destroy the perfect quantization of the braiding phase during the non-adiabatic, unitary time-evolution of the system. One is seen as a reduced fidelity due to the spectral losses, cf. Fig. 2 and the discussion above, which are expected to be larger for smaller protective gap ($\sim\Delta$), lower ramp time $T$, and rougher ramp functions. Further, the condition of exponential localization of the two Majorana zero modes is invalidated when the MBS localization length becomes comparable with the chain size, see Fig. 3. In this case, the hybridization between MBSs is not exponentially small, and the braiding phase encodes dynamic effects stemming from this hybridization. Both the minimal protective gap as well as Majorana localization length depend on the $p$-wave pairing amplitude $\Delta$. Here, $\xi_{\infty}(\Delta)=1/\ln\sqrt{(J+|\Delta|)/(J-|\Delta|)}$ corresponds to the localization length in the infinitely long Kitaev chain with $\mu=0$ where exact results are easily obtained. In Fig. 3a we show $\xi_{\infty}(\Delta)$ and the finite-length $\xi(\Delta,L)$ as a function of $\Delta$. The latter is determined by fitting the weights of the fermionic representation of the MBSs to an exponentially decaying function modulated by Friedel oscillations with period $\pi$ Stoudenmire et al. (2011). $L=11$ and $L=6$ chain sizes correspond to the maximal and minimal length of the topological regions during the braiding in the symmetric Y-junction with equal-size legs $\ell=5$, cf. Fig. 1. Perfectly localized MBSs, $\xi\to 0$, are obtained at $J=\Delta=1$. For $\Delta<1$, MBSs acquire a finite width that increases with decreasing values of $\Delta$, and the localization length $\xi(\Delta)$ diverges for vanishing $\Delta$. We note that the numerical result $\xi(\Delta,L=6)$ starts do deviate from $\xi_{\infty}(\Delta)$ at $\Delta\lesssim 0.5$, indicating that MBSs are not isolated anymore and the independent-particle picture Ivanov (2001); Alicea et al. (2011) breaks down.

In Fig. 4 we show the braiding phase and fidelity losses as a function of the Majorana localization length $\xi(\Delta)$ for the two ramp protocols and various step times $T$. Data in Fig. 4 (bottom row) shows that fidelity losses (10) reduce by orders of magnitude upon increasing the ramping time period $T$ from $250$ to $3000$. In addition, losses can be reduced significantly (by a power of $2$) upon employing the smoother sine-squared ramp function (first derivative is continuous) instead of the simpler guillotine one (first derivative jumps), agreeing with the power-law scaling identified by Knapp et al. Knapp et al. (2016). Similarly, non-adiabatic effects measured directly in the exchange phase (13), see Fig. 4a (top), which can be assessed by the strength of phase-variations as function of the braid parameters (i.e. the sensitivity of the braid phase against small deviations in the protocol), become smaller. In the adiabatic limit for large $T$, these fluctuations become negligible and one obtains smooth curves $\Phi(\xi,\Delta,\dots)$. The braiding phase is stable and closest to $\pi/2$, independent of $T$, for the ideal parameter set $\Delta=J=1$. As expected, for non-ideal values $\Delta<1$, stable braid phases are obtained only for larger times $T$. There are, however, significant deviations from the desired $\pi/2$ braiding phase starting from $\xi>1$ (corresponding to $\Delta<0.7$), which persist largely independent of the protocol time $T$ or smoothness of the ramp function. These can be identified with a break-down of the independent-particle picture of Majorana braiding, where ideal braid statistics do no longer apply. The braid phase in this sense is more sensitive to finite Majorana wave-function overlaps than spectral features like the MBS hybridization, or the predicted localization length in Fig. 3. We emphasize that even though the braid process appears quasi-adiabatic for $T\geqslant 500$ (or $T\geqslant 1500$ for guillotine ramp) and spectral losses even for large $\xi$ are substantially smaller than under quicker braid routines, only the exchange phase itself encodes this dramatic breakdown of braiding. Matching the Majorana wave-function weights or checking for adiabacity alone may thus be too loose a criterion for successful Majorana braid protocols Karzig et al. (2015a); Amorim et al. (2015).

Interactions.

Finally, we investigate the influence of a nearest neighbor interaction $V\neq 0$, see Eqs. (1) and (2), on the braiding operation. In general, interactions in Majorana systems modify the localization length and change the bulk energy gap Gangadharaiah et al. (2011); Sela et al. (2011); Stoudenmire et al. (2011); Thomale et al. (2013). Large values of interaction strength $V$ of either sign eliminate the topological phase. Here, we only consider values $V$ for which the topological phase is not destroyed at $\mu=0$, and where $\mu=4$ is still large enough to drive the corresponding junction leg into the trivial phase. In Fig. 5, we plot the braiding phase $\Phi(V;\Delta,\varphi,T)$ and loss of fidelity $w^{e/o}_{\mathrm{loss}}(V;\Delta,\varphi,T)$ as a function of $V$. We consider a symmetric Y-junction with equal-size legs $\ell=5$, now with pairing strengths $\Delta=1.0,~{}0.8$, and $1.1$. We only show results for the sine-squared ramp function with $T=750$, which in the non-interacting case was found to be sufficient to ensure close-to-adiabatic conditions, cf. Fig. 4. The Majorana localization length $\xi(V;\Delta,L)$ is depicted for Kitaev chains of up to $L=22$ sites, and obtained by fitting the weights of the free fermionic representation of the MBS Stoudenmire et al. (2011). For the ideal parameter set $\Delta=J=1$, the Majorana localization length increases with either increasing repulsive or attractive interaction $V$. This effect turns out to be larger for an attractive potential, where $\xi>1$ for $|V|>0.35$, and the braiding phase starts to deviate significantly from $\pi/2$. For $\Delta=0.8$, the minimum of $\xi(V;\Delta=0.8,L)$ is shifted to $V\approx-0.35$, where for the given $T$ one also observes a stabilization of the correct braiding phase. For large $\Delta=1.1$, the situation is reversed since then repulsive interactions lead to a reduction of the effective superconducting gap, pushing the system closer to the ideal point $\Delta_{\mathrm{eff}}=J$. The fidelity losses mainly depend on the braid velocity $\sim\ell/T$, and hence remain of similar magnitude for moderate $V$. As a weak trend, we find that small attractive or repulsive interactions play a stabilizing role for the braiding operation, depending on the relative size of the bare couplings $J$ and $\Delta$. Again, these effects appear not to be systematically encoded in spectral losses alone, which one can investigate by tracking the Majorana wave-function weights Amorim et al. (2015), but rather follow directly the behavior of the MBS localization length $\xi$, cf. both Figs. 4 and 5.

Discussion and conclusions.

Previous studies of Majorana braiding have predominantly focused on effective models of the low-energy (Majorana) sector Alicea et al. (2011); Clarke et al. (2011), introducing non-adiabatic effects and other sources of error on an effective but not microscopic level Cheng et al. (2011); Karzig et al. (2013); Scheurer and Shnirman (2013). Alternatively, the fidelity of the dynamically evolved state was considered with respect to some ideal, adiabatic reference state Karzig et al. (2015a, b); Amorim et al. (2015). Resulting estimates for the minimal protocol time for near-adiabatic time-evolution in Majorana systems were obtained by many groups Karzig et al. (2013); Scheurer and Shnirman (2013); Karzig et al. (2015a, b); Amorim et al. (2015). As an overarching implication, until recently, it was assumed that by just executing the braiding manipulations sufficiently slowly, one could recover the ideal operation completely. This result, however, disregards the generation and propagation of non-adiabatic excitations Knapp et al. (2016), which in turn directly affect the geometric braiding phase Pedrocchi and DiVincenzo (2015).

In our work, instead of further extending previous studies to more complex setups and protocols Karzig et al. (2015a, b); Aasen et al. (2016), we have revisited the most elementary case of a single MBS braiding operation in a Y-junction of three Kitaev chains Alicea et al. (2011). We considered the full Hilbert space, including the quasiparticle background beyond the low-energy Majorana sector, and investigated several error sources that destroy the perfect braiding phase during the non-adiabatic, unitary time-evolution of the system. By adopting the Bargmann vertex formula for the geometric phase, customarily used in the characterization of geometric quantum gate operations Mukunda and Simon (1993), we identified explicitly the MBS braiding phase, and analyzed its behavior during the braid. We then accounted for both spectral losses (through the fidelity) and direct braiding phase errors (by calculating the geometric phase evolution explicitly) in a numerically exact framework, in order to understand in detail how non-adiabacity, the Majorana mode localization, and interactions affect the braiding operation.

We find that due to the induced non-adiabatic leakage of wavefunction weight out of the ground state into excited states of the final Hamiltonian, the acquired braiding phase does not remain constant and evolves linearly in time even after the braiding process is finished, cf. Fig. 2. Consequently, the final (qubit) state will “dephase” with time, even in the zero-temperature case considered here. While spectral losses can be reduced by orders of magnitude through increasing the ramping time period $T$ or considering smoother ramp functions, cf. Fig. 4, it is clear that in extended protocols such as necessary for serious quantum computations, non-adiabatic errors will accumulate. Our results, based on a numerically exact treatment of the full Hilbert space, thus support and extend the findings of Pedrocchi and DiVincenzo Pedrocchi and DiVincenzo (2015) and Knapp et al. Knapp et al. (2016) that employed an effective description of the low-energy Majorana sector (incl. quasiparticle background via baths/dissipation).

By varying the superconducting wire bulk gap, we additionally resolved how the localization length of individual Majorana zero modes affects the braiding operation, beyond simple dynamical errors induced through a finite MBS hybridization Cheng et al. (2011); Amorim et al. (2015). For the studied Y-junction configuration in Fig. 1, a MBS localization length $\xi>1$ already causes significant deviations from the perfect braiding phase $\pi/2$, since the independent-particle picture of Majorana braiding Ivanov (2001); Kitaev (2001); Alicea et al. (2011) breaks down. We emphasize that this finite-size breakdown occurs independent of adiabacity, cf. Fig. 4. It implies that measures of adiabacity – such as the fidelity or Majorana wave function weights alone – may suggest overoptimistic efficiencies of the braiding operation. At the same time, the behavior in Fig. 4 shows how an extraction of the geometric phase (12) can yield valuable information about the anyonic statistics encoded in the many-body states of a complex system. Finally, nearest neighbor interaction in the wires is found to predominantly affect the braiding operation through its impact on the wire bulk gap. Depending on the non-interacting gap versus hopping strength, for weak attractive (repulsive) interaction, the braiding operation is more stable due to a change of the effective superconducting gap and the reduction of the MBS localization length, cf. Refs. Gangadharaiah et al. (2011); Stoudenmire et al. (2011); Sela et al. (2011).

In summary, finite operation times and non-ideal parameter settings in (simple) Majorana braiding schemes pose serious constraints on the accuracy of braiding operations. We found that even for the ideal case of a closed, zero-temperature Kitaev chain system, non-adiabatic errors and finite Majorana hybridization can become bottlenecks for the fidelity and feasibility of (extended) Majorana braiding routines towards quantum computation. An accumulation of non-adiabatic errors has to be avoided, likely leading to unfavorable finite-time scaling if no explicit quasiparticle relaxation mechanism is included. The optimal engineering of dissipation out of the excited-state sector, such as trapping of quasiparticles, as well as measurement-based (dissipative) topological quantum computation schemes Plugge et al. (2017); Vijay and Fu (2016a); Karzig et al. (2016) may thus behave favorably in comparison to conventional adiabatic-manipulation schemes Alicea et al. (2011); Hyart et al. (2013); Aasen et al. (2016), cf. Ref. Knapp et al. (2016). Similarly, finite Majorana overlaps have to be avoided to high accuracy, since these directly affect the braiding statistics encoded in the low-energy sector of the system – in fact, the braid statistics become “non-Majorana”, cf. Fig. 4 and discussion. As an order of magnitude comparison, we note that current-date nanowire and iron-chain architectures Albrecht et al. (2016); Deng et al. (2016); Zhang et al (2016); Pawlak et al. (2016); Feldman et al. (2017) find localization vs. device lengths of $\xi/L\simeq 1...5$, where our numerics indicate that assuming individual Majoranas (in terms of their braiding statistics) might be too optimistic. Advancing towards longer nanowire devices or a reduction of the MBS localization length, also for the sake of measurement-based qubit experiments Plugge et al. (2017); Vijay and Fu (2016a); Karzig et al. (2016), is thus clearly desirable.

Finally, an analysis of braiding and computation protocols in extended Majorana wire networks Landau et al. (2016); Plugge et al. (2016, 2017); Vijay and Fu (2016a); Karzig et al. (2016) will be the subject to future work, using more sophisticated techniques applicable to larger system sizes. Clearly, these measurement-based schemes for topological quantum computation should face the same detailed numerical investigation and scrutiny as the direct braiding and adiabatic-manipulation schemes.

Acknowledgements.

We thank A. Akhmerov, R. Egger, S. Frolov, and R. Lutchyn for fruitful discussions. This work was supported by DFG-SFB 1170 and ERC-StG-Thomale-TOPOLECTRICS-336012. M.S. acknowledges support by the Rustaveli national science foundation through the grand no. FR/265/6-100/14.