Multiparticle content of Majorana zero modes in the interacting $p$-wave wire

In the topological phase of $p$-wave superconductors, zero-energy Majorana quasiparticle excitations can be well defined in the presence of local density-density interactions. Here, we examine this phenomenon from the perspective of matrix representations of the commutator $\mathcal{H}=[H,•]$, with the aim of characterizing the multiparticle content of the many-body Majorana mode. To do this we show that, for quadratic fermionic systems, $\mathcal{H}$ can always be decomposed into subblocks that act as multiparticle generalizations of the Bogoliubov–de Gennes/Majorana forms that encode single-particle excitations. In this picture, density-density-like interactions will break this exact excitation-number symmetry, coupling different subblocks and lifting degeneracies so that the eigenoperators of the commutator $\mathcal{H}$ take the form of individual eigenstate transitions $|n\rangle \langle m|$. However, the Majorana mode is special in that zero-energy transitions are not destroyed by local interactions and it becomes possible to define many-body Majoranas as the odd-parity zero-energy solutions of $\mathcal{H}$ that minimize their excitation number. This idea forms the basis for an algorithm which is used to characterize the multiparticle excitation content of the Majorana zero modes of the one-dimensional $p$-wave lattice model. We find that the multiparticle content of the Majorana zero-mode operators is significant even at modest interaction strengths. This has important consequences for the stability of Majorana-based qubits when they are coupled to a heat bath. We will also discuss how these findings differ from previous work regarding the structure of the many-body Majorana operators and point out that this should affect how certain experimental features are interpreted.

Figure 1

The commutator $\mathcal{H}$ for a quadratic fermion model can be decomposed into blocks $A^{\left(s\right)}$ which encode energy transitions involving the same number of fermions $s$. In this example, we show the $2N+1$ blocks due to an $N=4$ fermion model. Quartic interacting terms (dark blue/purple) connect different $A^{\left(s\right)}$ blocks, breaking the symmetry responsible for the excitation-number conservation.

Figure 2

(a) Graphical representation of the block diagonal ${\mathcal{H}}_Q$. (b) Graphical representation of the quartic term ${\mathcal{H}}_I$. (c) The symmetry operator responsible excitation-number conservation of ${\mathcal{H}}_Q$.

Figure 3

In the figure we show the behavior of $|N_3^{\mathrm{\Gamma }}{|}^2$ for a system with $\mathrm{\Delta }=0.9, \mu =1.5$, and $t=1$. For $U=0$, the Majorana mode has a coherence length of $\xi \approx 1.11$ and thus for small values of $U$ permits us to use the full diagonalization method (FD) of Ref. [31]. For this purpose, we use a system size of $N_x=10$. We compare this with results using the commutator method (CM) with $(N_s,N_d)=(7,7)$ and $(7,9)$ for a system of size $N_x=40$. In the inset we plot the $|N_3^{\mathrm{\Gamma }}|$ to emphasize the linear dependence on $U$.

Figure 4

In the figure we plot $|N^{{\mathrm{\Gamma }}_5}{\left(U\right)|}^2$ and $1-|N_1^{gs}{\left(U\right)|}^2$ for a system of $N=40$ with $\mathrm{\Delta }=0.9, \mu =1.5$, and $t=1$. The inset shows that for these parameters $[1-|N_1^{gs}\left(U\right){|}^2]$ and $|N_5^{\mathrm{\Gamma }}{|}^2$ grow as the fourth power of $U$. The quartic dependence of $1-|N_1^{gs}{\left(U\right)|}^2$ would seem to indicate that single-particle contributions should dominate even at higher interaction strengths (note the small scale on the $y$ axis of the main figure). However, as we discuss in Sec. 5, the linear correlators are not a reliable indicator of the $N$-particle content of the quasiparticle itself. Indeed, we note that the $|N_5^{\mathrm{\Gamma }}{|}^2$ shown here is far smaller that $|N_3^{\mathrm{\Gamma }}{|}^2$ shown in Fig. 3 for the same system parameters.

Figure 5

In the figure we plot ${\alpha }_3$, the rate of slope of the $|N_3^{\mathrm{\Gamma }}\left(U\right)|$ as a function of $\mu$. Our results show best convergence for small coherence lengths and values of the $\mu$ that are far from the bottom of the band. We see a clear quadratic dependence centered around $\mu =2t$. The rate of growth of this term is therefore a minimum when we can linearize our dispersion. Fits of quadratic curves that show good overlap with the numerical data are also given. A system size of $N_x=50$, with cutoffs $(N_s,N_d)=(9,9)$, was used for this plot.

Figure 6

In the figure we plot the rate ${\alpha }_3$ as a function of $\mathrm{\Delta }$. The rate of growth clearly increases as $\mathrm{\Delta }$ decreases. The $N$-particle content of the many-body Majorana therefore clearly grows with the increased coherence length associated with a smaller superconducting gap. A system size of $N_x=50$, with cutoffs $(N_s,N_d)=(9,9)$, was used for this plot.

Figure 7

Fermionic doubling: in the $\mathrm{\Sigma }$ basis, the commutator $[H,•]$ can be understood as two disconnected copies of original system.

Figure 8

In the figure we show how $|N_3^{\mathrm{\Gamma }}|$ grows for a fixed value of $\mathrm{\Delta }=0.8$ and different values of $\mu \in [0.7,2.1]$. The values of ${\alpha }_3$ at different parameters represent the slopes of these straight lines. A system size of $N_x=50$, with cutoffs $(N_s,N_d)=(9,9)$, was used for this plot.

Figure 9

In the figure we show how $|N_5^{\mathrm{\Gamma }}{|}^{1/2}$ grows for a fixed value of $\mathrm{\Delta }=0.7$ and different values of $\mu \in [0.7,2.1]$. The values of ${\alpha }_5^{1/2}$ at different parameters represent the slopes of these straight lines. A system size of $N_x=50$, with cutoffs $(N_s,N_d)=(9,9)$, was used for this plot.

Figure 10

In the figure we plot the rate ${\alpha }_5^{1/2}$ as a function of $\mu$. We again see a quadratic dependence about the linearized dispersion point at $\mu =2$. The plot shows some numerical instability. Note that to fit these curves we take the $|N_5{|}^{1/4}|$ and very small numerical errors get magnified to some degree. A system size of $N_x=50$, with cutoffs $(N_s,N_d)=(9,9)$, was used for this plot.

Figure 11

In the figure we plot the rate ${\alpha }_5^{1/2}$ as a function of $\mathrm{\Delta }$. A system size of $N_x=50$, with cutoffs $(N_s,N_d)=(9,9)$, was used for this plot.