Exact zero modes in twisted Kitaev chains

TL;DR

This paper demonstrates the existence and properties of exact zero modes in twisted Kitaev chains with generalized boundary conditions, including their localization, conditions for occurrence, and robustness against disorder and interactions.

Contribution

It provides a comprehensive analytical study of zero modes in twisted Kitaev chains, including conditions for their existence and their robustness, extending understanding beyond standard boundary conditions.

Findings

Exact zero modes exist in large topological Kitaev chains with twisted boundaries.

Zero modes are localized and can be rigorously characterized for specific boundary phase parameters.

Zero modes persist even with disorder and interactions, indicating robustness.

Abstract

We study the Kitaev chain under generalized twisted boundary conditions, for which both the amplitudes and the phases of the boundary couplings can be tuned at will. We explicitly show the presence of exact zero modes for large chains belonging to the topological phase in the most general case, in spite of the absence of "edges" in the system. For specific values of the phase parameters, we rigorously obtain the condition for the presence of the exact zero modes in finite chains, and show that the zero modes obtained are indeed localized. The full spectrum of the twisted chains with zero chemical potential is analytically presented. Finally, we demonstrate the persistence of zero modes (level crossing) even in the presence of disorder or interactions.

Equations

$$H_{\mathrm{bulk}}$$

$$-\sum _{j=1}^L{\mu }_j\left(c_j^{†}{c}_j-\frac{1}{2}\right),$$

$$H_{\mathrm{boundary}}=\mathsf{b}\left[-t(e^{\text{i}{\varphi }_1}{c}_L^{†}{c}_1+h.c.)+\Delta (e^{\text{i}{\varphi }_2}{c}_L{c}_1+h.c.)\right],$$

$$P{\mid }_{(0,0)}\times P{\mid }_{(\pi ,\pi )}=\{\begin{array}{ll}1& \left|\mu /2t\right|>1,\\ -1& \left|\mu /2t\right|<1,\end{array}$$

$$a_j:=c_j+c_j^{†},b_j:=(c_j-c_j^{†})/\text{i}.$$

$$H_{\mathrm{bulk}}$$

$$-\frac{\text{i}}{2}\sum _{j=1}^L{\mu }_j{a}_j{b}_j,$$

$$H_{\mathrm{boundary}}$$

$$=\frac{\text{i}\mathsf{b}}{2}\left[(t\sin {\varphi }_1-\Delta \sin {\varphi }_2)a_1{a}_L+(t\sin {\varphi }_1+\Delta \sin {\varphi }_2)b_1{b}_L\right]$$

$$-\frac{\text{i}\mathsf{b}}{2}\left[(t\cos {\varphi }_1+\Delta \cos {\varphi }_2)a_1{b}_L-(t\cos {\varphi }_1-\Delta \cos {\varphi }_2)b_1{a}_L\right].$$

$$H=\text{i}t\sum _{j=2}^{L-1}{b}_j{a}_{j+1}+\text{i}t{b}_1\left(a_2+b_L\right).$$

$$H=\frac{\text{i}}{4}\sum _{i,j}{d}_i{\left[M_L\right]}_{ij}{d}_j$$

$$M_L=\left(\begin{array}{ccccc}\mathsf{a}\hat{m}& {\hat{t}}_0& & & \mathsf{b}{\hat{t}}_1\\ -{{\hat{t}}_0}^T& \hat{m}& & & \\ & \ddots & \ddots & \ddots & \\ & & & \hat{m}& {\hat{t}}_0\\ -\mathsf{b}{{\hat{t}}_1}^T& & & -{{\hat{t}}_0}^T& \mathsf{a}\hat{m}\end{array}\right),$$

$$\hat{m}$$

$${\hat{t}}_0$$

$${\hat{t}}_1$$

$$\begin{array}{rl}\mathrm{Pf}{M}_L=& \mathrm{Pf}{\tilde{M}}_L^{(\mathrm{open})}+{\mathsf{b}}^2\left(t^2-{\Delta }^2\right)\mathrm{Pf}{\tilde{M}}_{L-2}^{(\mathrm{open})}\\ & -\mathsf{b}\left(t\cos {\varphi }_1+\Delta \cos {\varphi }_2\right){\left(t+\Delta \right)}^{L-1}+\mathsf{b}\left(t\cos {\varphi }_1-\Delta \cos {\varphi }_2\right){\left(-(t-\Delta )\right)}^{L-1},\end{array}$$

$$\begin{array}{r}\mathsf{b}\left(t\cos {\varphi }_1+\Delta \cos {\varphi }_2\right)-\mathsf{b}\left(t\cos {\varphi }_1-\Delta \cos {\varphi }_2\right){\left(-\frac{t-\Delta }{t+\Delta }\right)}^{L-1}=\left(t+\Delta \right)\frac{\mathrm{Pf}{\tilde{M}}_L^{(\mathrm{open})}}{{\left(t+\Delta \right)}^L}+{\mathsf{b}}^2\left(t-\Delta \right)\frac{\mathrm{Pf}{\tilde{M}}_{L-2}^{(\mathrm{open})}}{{\left(t+\Delta \right)}^{L-2}}.\end{array}$$

$$\mathrm{Pf}{\tilde{M}}_L^{\left(\mathrm{open}\right)}$$

$$=-\mu \mathrm{Pf}{\tilde{M}}_{L-1}^{\left(\mathrm{open}\right)}-\left(t^2-{\Delta }^2\right)\mathrm{Pf}{\tilde{M}}_{L-2}^{\left(\mathrm{open}\right)},$$

$$\mathrm{Pf}{\tilde{M}}_L^{\left(\mathrm{open}\right)}\sim {\Lambda }^L,$$

$$\Lambda =-\frac{\mu }{2}\left(1+\sqrt{1-\frac{4\left(t^2-{\Delta }^2\right)}{{\mu }^2}}\right)$$

$$\left|\frac{\Lambda }{t+\Delta }\right|<1,$$

$$\left|\frac{\mu }{2t}\right|<1.$$

$$\mathsf{b}{\left(t+\Delta \right)}^{L-1}\left[\left(t\cos {\varphi }_1+\Delta \cos {\varphi }_2\right)-\epsilon \left(L\right)\right]=0,$$

$$\epsilon (L)=O\left({\left(\frac{t-\Delta }{t+\Delta }\right)}^{L-1}\right)+O\left({\left(\frac{\Lambda }{t+\Delta }\right)}^{L-1}\right)$$

$$\mathrm{Pf}{M}_L^{\left(\mathrm{open}\right)}=-\mu \mathrm{Pf}{M}_{L-1}^{\left(\mathrm{open}\right)}-{\tau }^2\mathrm{Pf}{M}_{L-2}^{\left(\mathrm{open}\right)}$$

$$\mathrm{Pf}{M}_1^{\left(\mathrm{open}\right)}=-\mu ,\mathrm{Pf}{M}_2^{\left(\mathrm{open}\right)}={\mu }^2-{\tau }^2,$$

$$\mathrm{Pf}{M}_L^{\left(\mathrm{open}\right)}={\tau }^L\cdot U_L\left(-\frac{\mu }{2\tau }\right).$$

$$\frac{\mu }{t}=2\sqrt{1-{\left(\frac{\Delta }{t}\right)}^2}\cos \frac{k\pi }{L+1}\left(k=1,2,\cdots ,L\right)$$

$$\mathrm{Pf}{M}_L$$

$$U_L\left(z\right)-\frac{2\mathsf{a}\left(1-\mathsf{a}\right)+{\mathsf{b}}^2}{{\mathsf{a}}^2}{U}_{L-2}\left(z\right)+{\left(\frac{\mathsf{a}-1}{\mathsf{a}}\right)}^2{U}_{L-4}\left(z\right)=0.$$

$$\left({\cos }^2\theta -\frac{2\mathsf{a}+{\mathsf{b}}^2}{{\mathsf{a}}^2}\right)\sin \left[\left(L-1\right)\theta \right]+\frac{1-2\mathsf{a}}{4\mathsf{a}}\sin \left[\left(L-3\right)\theta \right]=0.$$

$$2\cos L\theta =0,$$

Full text

I Introduction

Majorana zero modes have played an important role in condensed matter physics in recent years Wilczek09 ; Alicea12 ; LeijnseFlensberg12 ; ElliottFranz15 . The Majorana modes are the same as their own antimodes by definition, and they have been anticipated to appear as zero-energy bound states. There have been considerable efforts towards the realization of Majorana modes in condensed matter settings Mourik-12 ; Rokhinson-12 ; Das-12 ; DengYu-12 ; Finck-13 ; Churchill-13 ; Nadj-Perge-14 ; Higginbotham-15 . The existence of Majorana zero modes is of special interest because it can be applied to the physical construction of qubits for topological quantum computing Ivanov-01 ; Nayak-08 ; Pachos-12 . From an experimental point of view, it is essential to investigate the effects of disorder Brouwer-11 ; Lobos-12 ; DeGottardi-14 ; Altland-14 ; Crepin-14 and interactions FidkowskiKitaev10 ; FidkowskiKitaev11 ; Turner11 ; Gurarie11 ; Gangadharaiah-11 ; Stoudenmire-11 ; Hassler-12 ; Thomale-13 ; Katsura_Int_Majorana ; Rahmani-15-L ; Rahmani-15-B ; Gergs-16 . Furthermore, various theoretical aspects have been revealed, including the connection with supersymmetry Grover-14 ; Jian-15 ; Rahmani-15-L ; Rahmani-15-B ; Hsieh-16 , the generalization to parafermion modes Fendley12 ; Clarke-14 ; KlinovajaLoss14 ; Mong-14 ; Jermyn-14 ; Aris-15 ; Alicea-16 ; Iemini-16 , and the construction of topologically invariant defects Aasen-16 .

The emergence of Majorana zero modes in condensed matter systems was first proposed by Kitaev Kitaev01 . The Kitaev chain is a one-dimensional lattice model that describes a spin-polarized $p$-wave superconductor with open boundaries. This model possesses a topological phase with two-fold degenerate ground states which are robust against local perturbations that preserve the fermion parity symmetry. The origin of the ground-state degeneracy in the topological phase is the presence of Majorana zero modes. These zero-energy modes, often called strong zero modes, commute with the Hamiltonian and anti-commute with the fermionic parity, and are localized near the boundaries. However, the existence of a topological phase can be manifested by only a minimal degeneracy, the ground-state degeneracy, say, leading to the notion of a ‘weak zero mode Alicea-16 ’ that only commutes with the projected Hamiltonian onto the low-energy manifold. Moreover, the celebrated bulk-boundary correspondence is one of the most crucial properties for the topological phases of matter as a general rule, which topological insulators also exhibit HasanKane10 ; QiZhang11 ; shortcourse topo.ins. . Here, the fundamental question that should be addressed is whether Majorana zero modes (either strong or weak) can persist for a system without boundaries. Naively, one might think that there will be no zero-mode since there is no “edge” in the system.

In this paper, we will answer the above question by investigating the Kitaev chains with generalized twisted boundary conditions (TBCs), i.e., we arbitrarily control the amplitudes and the phases of the couplings on the boundaries. The usual periodic boundary condition (PBC) and the anti-periodic boundary condition (APBC) are included in the TBC as limiting cases. In practice, such a boundary condition can be realized by magnetic fluxes and Josephson junctions Kitaev01 ; Kwon04 ; Fu-09 ; Jiang-11 ; San-Jose12 ; Rokhinson-12 ; Deng-12 ; Deng-13 ; Beenakker-13-Josephson ; Hansen-16 ; Marra-16 ; Alase-16 ; Dmytruk-16 ; Giuliano-16 ; Cobanera-16 . Our work is motivated by the observation that the fermionic parities in the ground states of the Kitaev chain with the PBC and the APBC have opposite signs in the topological phase, hence indicating a level crossing when one continuously changes the parameters so as to connect the PBC with the APBC. When the open Kitaev chain resides in the topological phase, we show that for sufficiently large chains, Majorana zero modes do appear when the phase parameters at the boundary are tuned to specific values, in spite of the absence of the edges def_Majorana_zero . Then, we obtain the condition for the presence of Majorana zero modes in finite chains, as well as explicitly determine the spatial profile of the zero modes. We note that the emergence of Majorana zero modes (protected level crossing) in the fractional Josephson effect has been well-established since the seminal work of Kitaev Kitaev01 ; Jiang-11 ; Alicea12 ; Alase-16 ; Dmytruk-16 ; Giuliano-16 . However, we emphasize that the boundary couplings studied in this paper are more general than those in previous work, as their amplitudes and phases can be arbitrary. More importantly, our derivation of the condition for the Majorana zero modes does not require the assumption that the Majorana edge modes that appear in the absence of the boundary term do not hybridize much with the bulk states, which was implicitly assumed from the outset in most previous work.

We also investigate the robustness of the zero modes against disorder and interactions. We show for a particular set of parameters that the Majorana zero operators that commute with the Hamiltonian exist even in the presence of spatially varying couplings. The level crossing signaling a topological order is found to be robust against nearest-neighbor interactions, given that the bulk parameters are those of the interacting Kitaev chain in the topological phase. This means that the zero modes (at least in the weak sense Alicea-16 ) survive even in the presence of interactions.

The paper is organized as follows. In Sec. II, we introduce the model, and point out that there should be Majorana zero modes when we smoothly change the boundary conditions between the PBC and the APBC. In Sec. III, we directly compute the Pfaffian of the Hamiltonian and show the precise condition for the existence of zero modes. In Sec. IV, we present explicit forms of the zero modes and some other properties, including the full spectrum of the chain at zero chemical potential. In Sec. V, we investigate the effects of nearest-neighbor interactions. We numerically demonstrate the presence of the zero modes in the weak sense and determine the condition for the emergence of the zero modes. We conclude this paper in Sec. VI. In Appendix A, we calculate in detail the fermionic parity of the ground state of the Kitaev chain for the PBC and the APBC. In Appendix B, we present an explicit expression for the matrix which in general determines whether zero modes exist or not. In Appendix C, we give a detailed procedure for calculating all the eigen-energies of the chain with zero chemical potential. In Appendix D, the presence of the zero modes in interacting chains is confirmed for a solvable (frustration-free) modelKatsura_Int_Majorana .

II Model and Phases

We consider a system of spinless fermions on a chain of length $L$. For each site $j=1,2,\cdots,L$, we denote by $c^{\dagger}_{j}$ and $c_{j}$ the creation and the annihilation operators, respectively. We impose a twisted boundary condition, in which the amplitudes and the phases of the parameters on the boundaries can be tuned arbitrarily.

The Hamiltonian in question consists of two parts: $H=H_{\mathrm{bulk}}+H_{\mathrm{boundary}}$. The Hamiltonian for the bulk, $H_{\mathrm{bulk}}$, is described by

[TABLE]

where $t$ is the hopping amplitude and $\Delta$ is the $p$-wave pairing gap, both of which can be assumed to be nonnegative without loss of generality. Here, $\mu_{j}$ is the on-site (chemical) potential, and we set those on the boundaries as $\mu_{1}=\mu_{L}={\sf a}\mu$ and those in the bulk as $\mu_{j}=\mu$ ($j=2,3,\cdots,L-1$), where ${\sf a}\geq 0$ is a constant. In the case of ${\sf a}=1$, $H_{\mathrm{bulk}}$ reduces to the Kitaev’s $p$-wave superconductor model with open boundaries Kitaev01 , in which Majorana edge zero modes occur provided that it is in the topological phase $\left|\mu/2t\right|\leq 1$. The boundary Hamiltonian is given by

[TABLE]

where $\phi_{1},\phi_{2}\in\left[0,2\pi\right)$ are two independent phases that define the twisted boundaries, and ${\sf b}\geq 0$ is a constant. The TBC reduces to the open boundary condition (OBC) when ${\sf b}=0$. In the case of ${\sf a}={\sf b}=1$, the TBC boils down to the PBC for $\left(\phi_{1},\phi_{2}\right)=(0,0)$, or the APBC for $\left(\phi_{1},\phi_{2}\right)=(\pi,\pi)$.

Although the Hamiltonian $H$ does not conserve the total fermion number $F:=\sum_{j=1}^{L}c_{j}^{\dagger}c_{j}$, the parity of the fermion number, i.e., the fermion number modulo $2$, is conserved since $H$ commutes with $P:=(-1)^{F}$. Besides, time-reversal symmetry, i.e. invariance under complex conjugation, is respected if the chain has the periodic or anti-periodic boundaries.

It was pointed out in Ref. Greiter-14, that the fermionic parity of the ground state of $H$ with the PBC is odd in the topological phase. In order to see how the PBC is connected to the APBC through varying the twist parameters $(\phi_{1},\phi_{2})$, we diagonalize $H$ by the usual Fourier transform followed by a Bogoliubov transformation under the PBC and the APBC (see Appendix A). The fermionic parity $P\mid_{(0,0)}$ or $P\mid_{(\pi,\pi)}$ in the ground state for the PBC or the APBC is summarized as follows:

Even $L$

[TABLE] 2. 2.

Odd $L$

[TABLE]

Thus, we always have

[TABLE]

regardless of whether $L$ is even or odd. The fact that $P\mid_{(0,0)}$ and $P\mid_{(\pi,\pi)}$ have opposite signs in the topological phase $\left\lvert\mu/2t\right\rvert<1$ indicates that, if we consider an evolution of parameters along a continuous path $\mathbf{\Phi}\left(s\right)=(\phi_{1}(s),\phi_{2}(s))\in\mathbb{R}^{2}$ with $s\in\left[0,1\right]$, such that $\mathbf{\Phi}\left(0\right)=\left(0,0\right)$ and $\mathbf{\Phi}\left(1\right)=\left(\pi,\pi\right)$, the ground state must be degenerate at some particular $s\in\left[0,1\right]$ since $P$ cannot be changed discontinuously without gap closing. As we will show below, the degeneracy of ground states actually indicates the existence of Majorana zero modes Beenakker-13 ; Ortiz-14 ; Hegde-16 .

For each site $j$, we define the Majorana fermions by

[TABLE]

One can easily see that they satisfy the canonical anti-commutation relations for Majorana fermions. When written in terms of $a_{j}$ and $b_{j}$, the bulk Hamiltonian in Eq. (1) becomes

[TABLE]

and the Hamiltonian on the boundaries in Eq. (2) becomes

[TABLE]

We use in the following analysis the Majorana representation $a_{j}$ and $b_{j}$, instead of the ordinary fermionic representation $c_{j}$.

III Emergence of zero modes

In the previous section, it was demonstrated that the degeneracy of the ground states should happen between the PBC and the APBC, which implies the appearance of Majorana zero modes. In this section, we show that the phase parameters $\left(\phi_{1},\phi_{2}\right)$ can be tuned so that the Majorana zero modes appear if and only if the system belongs to the topological phase. Using the expansion formula for Pfaffians, we explicitly calculate the parameter conditions for the existence of the zero modes for large chains. Note that Nava et al. Giuliano-16 studied the Hamiltonian $H=H_{\rm bulk}+H_{\rm boundary}$ with $t=\Delta$, ${\sf a}=1$ in $H_{\rm bulk}$ and $\Delta=0$ in $H_{\rm boundary}$, and obtained the condition for the presence of Majorana zero modes using a different approach. Our method applies to the entire parameter region of the model and generalizes their results.

In order to see how the Majorana zero modes emerge from the twisted boundary conditions, let us first consider the simplest case where ${\sf a}={\sf b}=1,~{}t=\Delta$, $\mu=0$, and $\phi_{1}=\phi_{2}=\pi/2$. In this case, the Hamiltonian is represented schematically in Fig.1, and reduces to

[TABLE]

We see that the Majorana operator $a_{1}$ does not enter the Hamiltonian and thus corresponds to a Majorana zero mode. The other “edge” mode that commutes with the Hamiltonian is $\left(b_{L}-a_{2}\right)/\sqrt{2}$. We thus have a chain supporting zero modes even though there is no edge in the system.

In general, the Hamiltonian given by Eqs. (5) and (6) can be written in a quadratic form of Majorana fermions as

[TABLE]

where $d_{2i-1}:=a_{i},~{}d_{2i}:=b_{i}$ ($i=1,2,\cdots,L$). The ground state degeneracy at some points on $\mathbf{\Phi}\left(s\right)$, mentioned in the preceding section, implies the vanishing of $\det M_{L}$ at these points Remark_InfiniteSystem . The $2L\times 2L$ real skew symmetric matrix $M_{L}$ can be expressed as

[TABLE]

where the empty entries are zero and

[TABLE]

Using the expansion formula (see Sec. 2.8 of Ref. Laplace expansion, ), the Pfaffian of $M_{L}$ can be obtained as

[TABLE]

where $\tilde{M}^{\left(\mathrm{open}\right)}_{L}:=M_{L}|_{{\sf b}=0}$ is the matrix representation for the chain with open boundaries (but with arbitrary ${\sf a}$). Since the determinant of $M_{L}$ is related to $\mathrm{Pf}~{}M_{L}$ via $\det M_{L}=\left(\mathrm{Pf}~{}M_{L}\right)^{2}$, we impose $\mathrm{Pf}~{}M_{L}=0$, yielding

[TABLE]

When $L$ is large, the phase parameters $\left(\phi_{1},\phi_{2}\right)$ that exactly satisfy $\mathrm{Pf}~{}M_{L}=0$ exist only when the right hand side of Eq. (14) converges. Note that $|(t-\Delta)/(t+\Delta)|\leq 1$, because both $t$ and $\Delta$ are nonnegative. From the recurrence relation of $\mathrm{Pf}~{}M^{\mathrm{\left(open\right)}}_{L}$,

[TABLE]

the asymptotic form of $\mathrm{Pf}~{}\tilde{M}_{L}^{\left(\mathrm{open}\right)}$ for large $L$ is obtained as

[TABLE]

where

[TABLE]

is the solution of the characteristic equation of the recurrence relation Eq. (15) with the largest absolute value. Then, the condition that the right hand side of Eq. (14) converges for large $L$ is given by

[TABLE]

which just reduces to the condition for the topological phase:

[TABLE]

In this case, Eq. (14) becomes

[TABLE]

where

[TABLE]

is a small term for sufficiently large $L$. Thus, Eq. (20) gives an $L$-dependent solution $\left(\phi_{1}\left(L\right),\phi_{2}\left(L\right)\right)$ of $\mathrm{Pf}~{}M_{L}=0$.

In the case of $\phi_{1}=\phi_{2}$, Eq. (20) reduces to $\phi_{1}=\phi_{2}\simeq\pi/2$, $3\pi/2$ [Remark_Duality_Defect, ]. Since $\left(\phi_{1},\phi_{2}\right)=(0,0)$ and $(\pi,\pi)$ are separated by the curve described by Eq. (20) on a $\phi_{1}$-$\phi_{2}$ plane (see Fig.2), an arbitrary path ${\bm{\Phi}}(s)$ from $(0,0)$ to $(\pi,\pi)$ must intersect the curve when the system belongs to the topological phase. Note that Eq. (20) applies to arbitrary boundary conditions, while Eq. (3) is valid only for the PBC and the APBC. Therefore, the intersection between ${\bm{\Phi}}(s)$ and the curve Eq. (20) can be viewed as a general condition which characterizes the topological phase.

IV Zero modes for $\phi_{1}=\phi_{2}=\pi/2$

So far we have confirmed the presence of Majorana zero modes in large chains belonging to the topological phase. However, the spatial profile of the zero modes is not obtained and it is unclear whether the zero modes are localized. In this section, we discuss the explicit forms of the zero modes and their properties for finite chains with $\phi_{1}=\phi_{2}=\pi/2$. We first focus on a simple discussion by the Chebyshev polynomials and then determine explicit forms of the zero modes. Moreover, we present the full spectrum of the chain with $\mu=0$, and finally, we demonstrate that the topological order survives even in the presence of disorder.

IV.1 Conditions for the presence of exact zero modes by the Chebyshev polynomials

First of all, we consider the Majorana chain with open boundaries. We here set ${\sf a}=1$ in particular. In this case, the following recurrence relations and initial conditions hold:

[TABLE]

where we define $\tau:=\sqrt{t^{2}-\Delta^{2}}$, which can be either real or imaginary depending on the sign of $t^{2}-\Delta^{2}$. Note that $\mathrm{Pf}~{}M_{1}^{\left(\mathrm{open}\right)},~{}\mathrm{Pf}~{}M_{2}^{\left(\mathrm{open}\right)}$ are defined so that they are compatible with the recurrence relation Eq. (22) and $\mathrm{Pf}~{}M_{3}^{\left(\mathrm{open}\right)},~{}\mathrm{Pf}~{}M_{4}^{\left(\mathrm{open}\right)}$. Because these recurrence relations and initial conditions are the same as those of the Chebyshev polynomials of the second kind $U_{L}\left(z\right)$ Definition_Chebyshev , we could express $\mathrm{Pf}~{}M_{L}^{\left(\mathrm{open}\right)}$ as

[TABLE]

The necessary and sufficient condition for the presence of exact zero modes in finite chains is $\mathrm{Pf}~{}M_{L}^{\left(\mathrm{open}\right)}=0$. Since $U_{L}\left(z\right)$ has zeros only in the interval $z\in\left[-1,1\right]$, exact zero modes for finite $L$ appear for $\mu^{2}<4\left(t^{2}-\Delta^{2}\right)$ and on the curves

[TABLE]

in the $\Delta/t$-$\mu/t$ phase diagram (see Fig. 3), reproducing the results in Ref. Kao-14, . For the remaining regions $\mu^{2}>4\left(t^{2}-\Delta^{2}\right)$, exact zero modes are absent even when the chain belongs to the topological phase $\left|\mu/2t\right|<1$ Kao-14 . The following discussions on the Majorana chains with twisted boundaries become simple because $\mathrm{Pf}~{}M_{L}^{\left(\mathrm{open}\right)}$ can be expressed as the Chebyshev polynomials Wakatsuki-14 .

For a chain with twisted boundaries, especially with $\phi_{1}=\phi_{2}=\pi/2$, we find the following relation between the Pfaffian under the TBC and that under the OBC,

[TABLE]

which leads to the condition that exact zero modes appear in finite chains ($z:=-\mu/2\tau$)

[TABLE]

Letting $z=:\cos\theta$ ($\theta\in\mathbb{C}$), a simple calculation gives

[TABLE]

It seems difficult to solve this equation in general, but some simplification occurs in the following cases:

${\sf a}={\sf b}=1$

In this case, Eq. (28) reduces to

[TABLE]

with solutions

[TABLE] 2. 2.

${\sf a}=1/2$

In this case, the condition Eq. (28) reduces to

[TABLE]

with solutions (for $k=1,2,\cdots,L-2$)

[TABLE]

indicating that exact zero modes appear for all $L$ for the chain with ${\sf a}=1/2$ and $\mu=\pm 2\sqrt{t^{2}-\Delta^{2}}\sqrt{1+{\sf b}^{2}}$. We can also prove the contrary, that is, the necessary condition that exact zero modes appear for all $L$ is ${\sf a}=1/2,~{}\mu=\pm 2\sqrt{t^{2}-\Delta^{2}}\sqrt{1+{\sf b}^{2}}$, except for the trivial case $t=\Delta,~{}\mu=0$. To prove this, we first note that $\mathrm{Pf}~{}M_{L}$ obeys the recurrence relation for $L\geq 5$, which is the same as Eq. (22). We set

[TABLE]

so that they are compatible with the recurrence relation and $\mathrm{Pf}~{}M_{3},~{}\mathrm{Pf}~{}M_{4}$. Then, if the exact zero modes exist for all $L$ with the specific parameters, both $\mathrm{Pf}~{}M_{1}$ and $\mathrm{Pf}~{}M_{2}$ are required to be zero by the Euclidean algorithm, which completes the proof. 3. 3.

$\mu=0$

If $L$ is odd, all of $U_{L}\left(0\right),U_{L-2}\left(0\right)$, and $U_{L-4}\left(0\right)$ are zero and the exact zero mode condition is always satisfied, which leads to ever-presence of the exact zero modes. On the other hand, if $L$ is even, the identities $U_{L}\left(0\right)=-U_{L-2}\left(0\right)=U_{L-4}\left(0\right)$ lead to ${\sf b}^{2}+1=0$, which implies never-presence of the exact zero modes. Actually, the full spectrum of the Hamiltonian can be determined analytically for $\mu=0$, as will be presented in IV.3. 4. 4.

$\mu=\pm 2\tau\left(\neq 0\right)$

This case is of particular importance in the discussion of solvable interacting Majorana chains Katsura_Int_Majorana , where $\mu=\pm 2\tau$ is just the frustration-free condition in the noninteracting limit. In this case, the identities $U_{L}\left(1\right)=L+1$ and $U_{L}\left(-1\right)=\left(-1\right)^{L}\left(L+1\right)$ give

[TABLE]

which leads to the necessary condition ${\sf a}=1/2$ and ${\sf b}=0$ for the chain to have exact zero modes for arbitrary $L$. In other words, exact zero modes do not appear when $\mu=\pm 2\tau$ and at the same time the chain has boundary terms.

IV.2 Explicit forms of zero modes

In this subsection, we show explicit forms of the zero modes. A zero mode $\bm{\Psi_{0}}$ is expressed as $\left(\bm{\Psi_{0}}\right)_{2i-1}=A_{i},~{}\left(\bm{\Psi_{0}}\right)_{2i}=B_{i}$ ($i=1,2,\cdots,L$) and satisfies the following relations for the bulk ($i=1,2,\cdots,L$)

[TABLE]

where we have introduced four virtual variables $A_{0},B_{0},A_{L+1}$, and $B_{L+1}$ that can be incorporated into the following boundary conditions

[TABLE]

where $\tilde{\sf a}:=1-{\sf a}$.

The bulk conditions form second-order linear recurrence equations whose general solutions can be written as $A_{i}=A_{+}\lambda_{+}^{i}+A_{-}\lambda_{-}^{i},~{}B_{i}=B_{+}\lambda_{+}^{L-i+1}+B_{-}\lambda_{-}^{L-i+1}$, with

[TABLE]

Note that the absolute value of $\lambda_{+}$ and $\lambda_{-}$ must be less than $1$ in order for the zero modes to have finite normalization for large $L$, which implies the topological condition Eq. (19). The coefficients $A_{\pm}$ and $B_{\pm}$ are determined by the boundary conditions Eq. (36)

[TABLE]

where $X$ is a $4\times 4$ matrix defined by the parameters (see Appendix B for the specific form of $X$). Therefore, the necessary and sufficient condition for the existence of the zero modes is $\mathrm{det}~{}X=0$. We here remark that this condition can be applied to infinite systems as well as finite systems. In fact, for infinite chains, it is clear that $\det X$ will always become zero when the systems belong to the topological phase Eq. (19). Therefore, we again confirm the presence of the zero modes for large chains in the topological phase. By their explicit forms, it is evident that the zero modes are localized if they exist.

The condition obtained can be simplified in some specific cases. First, for a chain with open boundaries (${\sf a}=1,~{}{\sf b}=0$), $M$ becomes block diagonal and $\det X=0$ reduces to $\left(\lambda_{+}/\lambda_{-}\right)^{L+1}=1$, which recovers Eq. (25) obtained using the Chebyshev polynomials. Next, for a chain with $t=\Delta$, $\det X=0$ is always satisfied, which confirms the existence of exact zero modes for the chain in this case regardless of the other conditions.

Finally, in the case of the chain with $\mu=0$, we have

[TABLE]

The exact zero modes thus appear when the chain length $L$ is odd regardless of other conditions. For even $L$, exact zero modes never appear unless $t=\Delta$. However, as shown above, there always exist the zero modes (though not an exact one) in the thermodynamic limit. All of these observations are consistent with the previous results presented in IV.1. We will discuss the case of $\mu=0$ in more detail in the next subsection.

IV.3 Full spectrum for $\mu=0$

We have discussed the conditions for the existence of zero modes in the case of $\phi_{1}=\phi_{2}=\pi/2$. Besides the zero modes, we are also interested in determining the full spectrum of the chain. In particular, the chain with $\mu=0$ can be solved exactly as follows. We set $J:=t+\Delta,~{}f:=t-\Delta$ and assume that $J\neq 0$ and $f\neq 0$. Then the Hamiltonian with $\mu=0$ reads

[TABLE]

When the energy eigenvalues of $H$ are expressed as $E=\varepsilon/4$, the conditions for the bulk are ($i=1,2,\cdots,L$)

[TABLE]

and the conditions for the boundary are

[TABLE]

where the virtual variables $A_{0},B_{0},A_{L+1},B_{L+1}$ are defined by the bulk conditions Eq. (41). Since the spectrum is chiral, we concentrate on the nonnegative eigenvalues $E\geq 0$ in what follows.

The structure of the Hamiltonian depends on the parity of $L$. When $L$ is even, $H$ forms a nearest-neighboring closed chain (just like a Möbius ring) of length $2L$ with two defects (Fig. 4a) PLA . On the other hand, if $L$ is odd, the Hamiltonian is separated into two decoupled chains with nearest-neighbor hopping, each of which is closed and has one defect bond (Fig. 4b).

For even $L$, the solutions on the foregoing $2L$-chain are obtained by conducting a plane-wave expansion independently in two parts of the Möbius ring and combining them at the two defects. The exact eigenvalues are given by (see Appendix C)

[TABLE]

where the wave numbers $q$ are determined by the following quantization condition ($M:=L/2$, Fig. 5a):

[TABLE]

For $\Delta\neq 0$, the existence of zero modes is only possible for the wavenumber $q$ with nonvanishing imaginary part since the zero energy is out of the energy band. It follows from the reality of $\varepsilon$ that $q$ takes the form $q=m\pi+\text{i}q_{*}$, where $m\in\mathbb{Z}$ and $q_{*}\in\mathbb{R}$. Thus, the dispersion relation reads

[TABLE]

and the quantization condition given by Eq. (44) in the case of large $L$ will be

[TABLE]

After a straightforward calculation, we have an eigenvalue which is exponentially small for large even $L$:

[TABLE]

The analysis for odd $L$ can be performed in a similar way. The dispersion relation is the same as Eq. (43), and the wave numbers $q$ are determined by the following quantization condition ($N:=(L+1)/2$, Fig. 5b):

[TABLE]

The wavenumber determined by the first equality corresponds to the exact zero mode. Here, due to the fact that the spectrum is chiral and the number of eigenvalues is odd, there always exist exact zero modes for odd $L$. In fact, the Majorana zero operators are

[TABLE]

Since we assumed from the outset that $|f/J|<1$, the coefficients fall off exponentially in system size and hence each operator has finite normalization even in the $L\to\infty$ limit, i.e. $(\Psi_{i})^{2}={\rm const.}<\infty$ ($i=1,2$).

When $L$ is large, the quantization condition Eq. (48) except for that of the exact zero mode becomes

[TABLE]

where we have substituted $q=m\pi+\text{i}q_{*}$. Equation (51) does not have a solution because it is incompatible with the necessary condition for the existence of a zero mode, i.e., $m=[1+{\rm sgn}(Jf)]/2$. We therefore conclude that there does exist one zero mode, which is also the exact zero mode, for large odd $L$. It might be intriguing that just one zero mode appears regardless of the parity of $L$ for large chains, in spite of the fact that the exact zero mode exists only when $L$ is odd.

IV.4 Majorana zero operators in inhomogeneous chains with $\mu=0$

The Majorana zero operators localized at the boundary exist even in the presence of couplings varying over space. To see this, we consider the inhomogeneous chain with $\phi_{1}=\phi_{2}=\pi/2$ and $\mu=0$. The Hamiltonian in terms of $a_{j}$ and $b_{j}$ reads

[TABLE]

We here assume that the minimum of $\left\{J_{j}\right\}$ is larger than the maximum of $\left\{f_{j}\right\}$ : $\min_{j}\left|J_{j}\right|>\max_{j}\left|f_{j}\right|$. Then, when we set $K:=\lfloor\left(L+1\right)/2\rfloor$, the Majorana zero operators take the following forms:

[TABLE]

These zero operators remain normalizable even in the $L\to\infty$ limit. For the homogeneous chain with odd $L$, Eqs. (53) and (54) reduce to Eqs. (49) and (50), respectively. When the amplitude on the boundaries vanishes (${\sf b}=0$), the above zero modes boil down to those of the open chain.

In the case of odd $L$, the Majorana zero operators exactly commute with the Hamiltonian:

[TABLE]

On the other hand, in the case of even $L$, the zero operators do not exactly commute with the Hamiltonian. However, their commutators with $H$ are exponentially small in the system size:

[TABLE]

[TABLE]

implying that each commutes with $H$ in the limit $L\to\infty$. Therefore, the Majorana zero operators do exist in the $L\to\infty$ limit, regardless of the parity of $L$. This clearly demonstrates that the topological order survives even in the presence of disorder.

V Effects of nearest-neighbor interactions

In this section, we consider the generalization of $H$ to the case with nearest-neighbor interactions Hassler-12 :

[TABLE]

where $U$ is the strength of the interaction and $n_{j}:=c^{\dagger}_{j}c_{j}$ is the fermion number operator at site $j$. For the non-interacting Kitaev chains, we have demonstrated that the ground-state fermionic parity changes between the PBC and the APBC, when the system resides in the topological phase. Now we remark that a small $U$ that fails to close the many-body gap will not remove the parity switch, which leads to the same conclusion as the free case. To illustrate this heuristic picture, we have performed exact diagonalization of interacting chains. Figure 6 shows the evolution of the spectrum as a function of the phase parameter $\phi:=\phi_{1}=\phi_{2}$. In the topological phase, the level crossing between the two lowest-lying states occurs (see Fig. 6 (a)). On the other hand, this does not happen in the trivial phase (see Fig. 6 (b)). These results imply that the level crossing is a generic feature of topological phases even in the presence of nearest-neighbor interactions. Besides, the parity switches indeed take place for the exactly solvable interacting chains with the fine-tuned parameters satisfying ‘frustration-free conditionKatsura_Int_Majorana ’ (see Appendix D).

For free chains, we have also explicitly determined the parameter conditions for the two-fold degeneracy of the ground states. In particular, the chains with $\phi=\pi/2$ or $3\pi/2$ definitely satisfy the obtained conditions in the infinite-size limit. Now let us see what happens in interacting systems. We can observe the level crossing appears at around $\phi=\pi/2$ in Fig. 6 (a), hence it may be conjectured that the two-fold degeneracy of the ground states at $\phi=\pi/2$ is also the case with interacting chains in the infinite-size limit, irrespective of the details of the other parameters. We have performed exact diagonalization further to confirm the conjecture. Figure 7 provides the difference between the level crossing angle $\phi$ and $\pi/2$ as a function of the system size $L$. Remarkably, $\lvert\phi-\pi/2\rvert$ shows a significant decrease as $L$ is increased. As discussed in Sec. III, $\lvert\phi-\pi/2\rvert$ has exponentially small ($\sim\varepsilon(L)$ in Eq. (21) ) dependence on the system size $L$ in free chains. Fitting the data with an exponential ansatz, $|\phi-\pi/2|\propto\exp(-L/\xi)$, yields the curves shown in Fig. 7. The exponential decrease of $\lvert\phi-\pi/2\rvert$ in $L$ clearly demonstrates that the ground states are two-fold degenerate at $\phi=\pi/2$ in the infinite-size limit even for interacting chain. In addition, the ground-state degeneracy at the crossing point suggests the presence of zero modes at $\phi=\pi/2$. An interesting question is whether the zero modes that map one of the ground state to the other are strong or weak. This could be studied systematically following previous approaches Stoudenmire-11 ; Aris-15 , but we leave it for future work.

VI Conclusion

In this paper, we have studied the Kitaev chains under generalized twisted boundary conditions characterized by the phase parameters $(\phi_{1},\phi_{2})$. We found that the phases $(\phi_{1},\phi_{2})$ can be adjusted so that Majorana zero modes appear as long as the bulk couplings are those of the Kitaev chain in the topological phase. By computing the Pfaffian of the Hamiltonian matrix in the Majorana basis, we rigorously obtained the condition on $(\phi_{1},\phi_{2})$ for the presence of Majorana zero modes. The condition reduces to $\phi_{1}=\phi_{2}=\pi/2$ or $3\pi/2$ in the infinite-size limit when the constraint $\phi_{1}=\phi_{2}$ is imposed.

We then analyzed finite chains at $\phi_{1}=\phi_{2}=\pi/2$ and enumerated conditions on the other parameters under which exact Majorana zero modes exist. A particularly interesting case is $\mu=0$, where the exact zero modes must appear in a chain of odd length, irrespective of the details of the other parameters. The full energy spectrum for this case was analyzed in detail and the explicit expressions for the Majorana zero operators that commute with the Hamiltonian were obtained. The operators obtained are exponentially localized and hence normalizable in the infinite-size limit. We also showed that the presence of Majorana zero operators survive even in the presence of spatially varying couplings, provided that $\phi_{1}=\phi_{2}=\pi/2$ and $\mu=0$.

The robustness of the zero modes at least in the weak sense persists even in the presence of interactions, as demonstrated by our analytical and numerical results. Whether the zero modes at the level crossing point are strong or weak is an intriguing open question. It would also be interesting to see if and how the twisted boundary conditions lead to a level crossing in other systems such as parafermion Fendley12 ; Clarke-14 ; KlinovajaLoss14 ; Mong-14 ; Jermyn-14 ; Aris-15 ; Alicea-16 and XYZ chains Fendley-XYZ , which are not reducible to free fermions.

Acknowledgment

The authors thank Yutaka Akagi and Ken Shiozaki for valuable discussions, and Domenico Giuliano for bringing our attention to Ref. Giuliano-16, . H. K. was supported in part by JSPS KAKENHI Grant No. JP15K17719, No. JP16H00985, and No. JP15K21717. N. W. was supported by a startup fund from the Beijing Institute of Technology.

Appendix A Diagonalization of the Kitaev chain for the PBC and the APBC

The Hamiltonian $H$ with the PBC or the APBC can be easily diagonalized by the Fourier transform. Let $\psi_{k}$ be the Fourier transform of $c_{j}$. The annihilation operator $c_{j}$ is written in terms of $\psi_{k}$ as

[TABLE]

where $\mathcal{K}$ denotes the set of possible wavenumbers which depends on both the parity of $L$ and the boundary condition. Then the Hamiltonian in $k$-space is given by

[TABLE]

where

[TABLE]

For $k\neq 0,\pi$, the ground state of $H_{k}+H_{-k}$ is the vacuum of Bogoliubov quasiparticles:

[TABLE]

where $\lvert 0\rangle_{k}$ is the vacuum state of $\psi_{k},\psi_{-k}$, and $\alpha_{k/-k}$ is the annihilation operator of Bogoliubov quasiparticles, which has odd fermionic parity. Therefore, the ground state at each $k$ has even number of fermions unless $k=0$ or $\pi$. For $k=0,\pi$, the fermionic parity in the ground state of $H_{k}$ is given by

$H_{0}$

[TABLE] 2. 2.

$H_{\pi}$

[TABLE]

Now, it is easy to get the fermionic parity in the ground state and hence Eq. (3).

Appendix B The specific form of the coefficient matrix $X$ in Eq. (38)

The coefficient matrix $X$ in Eq. (38) which determines whether zero modes exist is

[TABLE]

In the case of OBC ($\tilde{\sf a}={\sf b}=0$), $X$ becomes block diagonal and

[TABLE]

Appendix C The detailed calculation of the full spectrum of the chain with $\mu=0$

In this Appendix, we present a detailed exposition of the calculation of the full spectrum of the chain with $\mu=0$. The results depend on the parity of $L$, but the way to get the spectrum is essentially the same. We first conduct a plane-wave expansion in the bulk, and then match the solutions at the boundaries. The calculation of the chain with odd $L$ is easier than that of the chain with even $L$: we have to evaluate the determinant of a $2\times 2$ matrix in the case of odd $L$, but on the other hand, we have to analyze a $4\times 4$ matrix in the case of even $L$.

C.1 Even $L$

The spectrum of the original $2L\times 2L$ Hamiltonian is equivalent to that of the following $2L\times 2L$ matrix:

[TABLE]

Let $\left(C_{1},~{}D_{1},~{}\cdots,C_{M},~{}D_{M},~{}E_{1},~{}F_{1},~{}\cdots,E_{M},~{}F_{M}\right)$ be an eigenvector of $H^{\left(\mathrm{even}\right)}$, where $M:=L/2$. The relation between $\{A_{j}\},~{}\{B_{j}\}$ and $\{C_{j}\},~{}\{D_{j}\},~{}\{E_{j}\},~{}\{F_{j}\}$ is that

[TABLE]

where $j=1,2,\cdots,M$. If we take an ansatz $C_{j}\sim ce^{\text{i}qj},~{}D_{j}\sim de^{\text{i}qj}$, the bulk conditions give

[TABLE]

and the determinant of the coefficient matrix should be zero for the existence of the eigenvectors, which is just the dispersion relation

[TABLE]

Remembering that we focus on the non-negative eigenvalues ($E,\varepsilon\geq 0$),

[TABLE]

Therefore, we can expand $\{C_{j}\},~{}\{D_{j}\}$ as follows

[TABLE]

We can expand $\{E_{j}\},~{}\{F_{j}\}$ in a similar way:

[TABLE]

Replacing $\{A_{j}\}$, $\{B_{j}\}$ in Eq. (42) with $\{C_{j}\}$, $\{D_{j}\}$, $\{E_{j}\}$, $\{F_{j}\}$, we obtain the boundary conditions:

[TABLE]

Substituting the plane wave expansions into the boundary conditions, the consistency condition for the wavenumber $q$ is obtained as

[TABLE]

where $X$ is the following $4\times 4$ matrix:

[TABLE]

The above equation has a nontrivial solution if the determinant of the coefficient matrix vanishes. After some calculations, this condition boils down to Eq. (44).

C.2 Odd $L$

The spectrum of the original $2L\times 2L$ Hamiltonian is equivalent to those of the two independent $L\times L$ matrices, one of which is

[TABLE]

We restrict our attention to the spectrum of this matrix since that of the other one is obtained by exchanging $(J,f)$ for $(-f,-J)$ in this one. Let $\left(C_{1},~{}D_{1},~{}\cdots,C_{N-1},~{}D_{N-1},~{}C_{N}\right)$ be an eigenvector of $H^{\left(\mathrm{odd}\right)}$, where $N:=(L+1)/2$. The relation between $\{A_{j}\},~{}\{B_{j}\}$ and $\{C_{j}\},~{}\{D_{j}\}$ is that

[TABLE]

where $j=1,2,\cdots,N-1$. The bulk conditions are the same as those of the case of even $L$, and the dispersion relation is thus the same:

[TABLE]

We can expand the eigenstates as

[TABLE]

The boundary conditions in terms of $\{C_{j}\},~{}\{D_{j}\}$ are

[TABLE]

Substituting the plane wave expansions into the boundary conditions, the consistency condition for the wavenumber $q$ is obtained as

[TABLE]

This equation has a nontrivial solution if the determinant of the coefficient matrix vanishes. This condition reads

[TABLE]

“$Je^{\text{i}q}+f=0$” means the existence of the exact zero mode. And the condition for the other Hamiltonian is obtained by swapping $(J,f)$ for $(-f,-J)$ in the above condition, which leads to the quantization condition Eq. (48).

Appendix D Parity switches of the ground states for the interacting chains satisfying the frustration-free condition

It was shown for ${\sf a}=1/2$ and ${\sf b}=0$ that the exact ground states of $H_{\rm{int}}|_{{\sf a}=1/2,{\sf b}=0}$ can be obtained, provided that the parameters satisfy the following ‘frustration-free’ condition Katsura_Int_Majorana

[TABLE]

The ground states are found to be two-fold degenerate and have opposite fermionic parities

[TABLE]

where $\alpha:=\sqrt{\cot\left[\arctan\left(2\Delta/\mu^{*}\right)/2\right]}$, $A_{L}^{\rm{(even/odd)}}:=A_{L}^{(+)}\pm A_{L}^{(-)}$ and

[TABLE]

Here, $|\rm{vac}\rangle$ is the vacuum state of $\{c_{j}\}$. The key to obtaining the exact ground states is the fact that the two-site states $e^{\pm\alpha c^{\dagger}_{j}}e^{\pm\alpha c^{\dagger}_{j+1}}|\rm{vac}\rangle$ minimize the following local Hamiltonian $h_{j}$ simultaneously

[TABLE]

In terms of $h_{j}$ ($j=1,2,\cdots L-1$) the total Hamiltonian is expressed as $H_{\rm int}=\sum^{L-1}_{j=1}h_{j}$. We note in passing that the phase diagram of the model has been obtained in Ref. Hassler-12, and the frustration-free line always resides in the topological phase Katsura_Int_Majorana .

One can show that, under the frustration-free condition Eq. (92), $|\Psi^{\rm{(odd)}}_{0}\rangle$ ($|\Psi^{\rm{(even)}}_{0}\rangle$) is the ground state of the interacting Kitaev chain under the PBC (APBC). We note that the same result has also been presented Shiozaki ; Shiozaki2 , but for the reader’s convenience we give an explicit proof of this fact. Since $|\Psi^{\rm{(odd)}}_{0}\rangle$ and $|\Psi^{\rm{(even)}}_{0}\rangle$ are the ground states of the bulk part of the Hamiltonian, it suffices to show that they are ground states of the boundary part $h_{L}$. The proof goes as follows:

PBC: $c_{L+1}=c_{1},c^{\dagger}_{L+1}=c^{\dagger}_{1}$.

We first note that any state of the form $e^{\pm\alpha c^{\dagger}_{L}}e^{\pm\alpha c^{\dagger}_{1}}\cdots|\rm{vac}\rangle$, is a ground states of $h_{L}$, where the part $\cdots$ denotes an arbitrary polynomial in $c^{\dagger}_{2},...,c^{\dagger}_{L-1}$. Consequently, any linear combination of such states is also a ground state of $h_{L}$. Then we rewrite $A_{L}^{(\rm{odd})}$ as

[TABLE]

and find that $\lvert\Psi^{\rm{(odd)}}_{0}\rangle\propto A_{L}^{\rm{(odd)}}\lvert\mathrm{vac}\rangle$ minimizes $h_{L}$. This is the unique ground state of $H_{\rm int}$ at $\mu=\mu^{*}$ with PBC, because the other state $\lvert\Psi^{\rm{(even)}}_{0}\rangle$ does not minimizes $h_{L}$. 2. 2.

APBC: $c_{L+1}=-c_{1},c^{\dagger}_{L+1}=-c^{\dagger}_{1}$.

We first note that any state of the form $e^{\pm\alpha c^{\dagger}_{L}}e^{\mp\alpha c^{\dagger}_{1}}\cdots|\rm{vac}\rangle$, is a ground states of $h_{L}$, where the part $\cdots$ denotes an arbitrary polynomial in $c^{\dagger}_{2},...,c^{\dagger}_{L-1}$. Consequently, any linear combination of such states is also a ground state of $h_{L}$. Then we rewrite $A_{L}^{(\rm{even})}$ as

[TABLE]

and find that $\lvert\Psi^{\rm{(even)}}_{0}\rangle\propto A_{L}^{\rm{(even)}}\lvert\mathrm{vac}\rangle$ minimizes $h_{L}$. This is the unique ground state of $H_{\rm int}$ at $\mu=\mu^{*}$ with APBC, because the other state $\lvert\Psi^{\rm{(odd)}}_{0}\rangle$ does not minimizes $h_{L}$.

The change of the parity between the PBC and the APBC for the interacting Kitaev model again indicates the ground-state degeneracy at some point on a path $(\phi_{1}(s),\phi_{2}(s))$ that connects $H_{\rm{int}}|_{{\sf a}={\sf b}=1;(0,0)}$ and $H_{\rm{int}}|_{{\sf a}={\sf b}=1;(\pi,\pi)}$. Therefore, the level crossing in the spectrum survives in the presence of interactions, under the frustration-free condition.