Braiding Majorana corner modes in a second-order topological superconductor

TL;DR

This paper proposes a device using a second-order topological superconductor to host Majorana corner modes, demonstrating their braiding and non-Abelian statistics, advancing potential quantum computing applications.

Contribution

It introduces a novel device concept for braiding Majorana corner modes in a second-order topological superconductor, including a design for experimental realization.

Findings

Successful braiding of Majorana modes with a non-trivial statistical phase π

Demonstration of control over Majorana modes via magnetic field rotation and pairing tuning

Confirmation of Majorana nature through fractional statistics and operator algebra

Abstract

We propose the concept of a device based on a square-shaped sample of a two-dimensional second-order topological helical superconductor which hosts two zero-dimensional Majorana quasiparticles at the corners. The two zero-energy modes rely on particle-hole symmetry (PHS) and their spacial position can be shifted by rotating an in-plane magnetic field and tuning proximity-induced spin-singlet pairing. We consider an adiabatic cycle performed on the degenerate ground-state manifold and show that it realizes the braiding of the two modes whereby they accumulate a non-trivial statistical phase $\pi$ within one cycle. Alongside with the PHS-ensured operator algebra, the fractional statistics confirms the Majorana nature of the zero-energy excitations. A schematic design for a possible experimental implementation of such a device is presented, which could be a step towards realizing non-Abelian braiding.

Tables

Table 1. Table 1: Symmetries ( 𝒮 𝒮 \mathcal{S} ) of the bare Hamiltonian ℋ ( 𝟎 , 𝟎 ) subscript ℋ 0 0 \mathcal{H}_{(\bm{0},\bm{0})} , on lines: particle-hole antiunitary antisymmetry ( 𝒫 𝒫 \mathcal{P} ), time-reversal antiunitary symmetry ( 𝒯 𝒯 \mathcal{T} ), chiral antisymmetry ( 𝒞 = 𝒫 ​ 𝒯 𝒞 𝒫 𝒯 \mathcal{C}=\mathcal{PT} ), inversion symmetry ( ℐ ℐ \mathcal{I} ), two-fold rotation (around j 𝑗 j axis, ℛ 2 j superscript subscript ℛ 2 𝑗 \mathcal{R}_{2}^{j} ), mirror reflection ( ℳ j subscript ℳ 𝑗 \mathcal{M}_{j} maps x j ↦ − x j maps-to subscript 𝑥 𝑗 subscript 𝑥 𝑗 x_{j}\mapsto-x_{j} ), antiunitary symmetries ( 𝒯 ​ ℛ 2 j 𝒯 superscript subscript ℛ 2 𝑗 \mathcal{TR}_{2}^{j} and 𝒯 ​ ℳ j 𝒯 subscript ℳ 𝑗 \mathcal{TM}_{j} ) and unitary antisymmetries ( 𝒞 ​ ℳ j 𝒞 subscript ℳ 𝑗 \mathcal{CM}_{j} ). For each symmetry, the equation of transformation is given, as well as the corresponding representation U 𝑈 U in Hamiltonian space. In the following six columns it is indicated whether the Hamiltonian ℋ ( 𝒃 , 𝒔 ) subscript ℋ 𝒃 𝒔 \mathcal{H}_{(\bm{b},\bm{s})} still possesses a certain symmetry, for various values of 𝒃 = ( b x , b y ) 𝒃 subscript 𝑏 𝑥 subscript 𝑏 𝑦 \bm{b}=(b_{x},b_{y}) and 𝒔 = ( s x , s y ) 𝒔 subscript 𝑠 𝑥 subscript 𝑠 𝑦 \bm{s}=(s_{x},s_{y}) – depicted as arrows of different colors in the table header. The parameters for three representative steps of the braiding process are considered, t ∈ { 0 , T / 16 , T / 4 } 𝑡 0 𝑇 16 𝑇 4 t\in\{0,T/16,T/4\} . For each step, two cases are compared: [ 𝒃 𝒃 \bm{b} as indicated and 𝒔 = ( 0 , 0 ) 𝒔 0 0 \bm{s}=(0,0) ] in the former column and [ 𝒃 𝒃 \bm{b} as indicated, 𝒔 𝒔 \bm{s} as indicated] in the latter column. For each column, the groups G 𝐺 {G} and N 𝑁 {N} explained in Eq.   11 are specified.

			​​​
	Equation of $𝒮$			$𝒃=(0,1)$	$𝒔=(1,0)$	$𝒃=(-0.5,1)$	$𝒔=(0.5,0)$	$𝒃=(-1,0)$	$𝒔=(0,1)$
$𝒮$	$\mathscr{H}(k_x,k_y)\Lsh$	$U\Lsh$	$G,N$ $\to$	$D_{2h},C_{2h}$	$C_{2v},C_{1h}$	$C_{2h},C_i$	$C_{1h},C_1$	$D_{2h},C_{2h}$	$C_{2v},C_{1h}$
$𝒫$	$-U\mathscr{H}{(-k_x,-k_y)}^{\ast }{U}^{†}$	${\sigma }_3{\tau }_1$		✓	✓	✓	✓	✓	✓
$𝒯$	$+U\mathscr{H}{(-k_x,-k_y)}^{\ast }{U}^{†}$	$i{\sigma }_2{\tau }_0$		$\times$	$\times$	$\times$	$\times$	$\times$	$\times$
$𝒞$	$-U\mathscr{H}(k_x,k_y)U^{†}$	${\sigma }_1{\tau }_1$		$\times$	$\times$	$\times$	$\times$	$\times$	$\times$
$ℐ$	$+U\mathscr{H}(-k_x,-k_y)U^{†}$	${\sigma }_0{\tau }_3$		✓	$\times$	✓	$\times$	✓	$\times$
${ℛ}_2^x$	$+U\mathscr{H}(k_x,-k_y)U^{†}$	${\sigma }_1{\tau }_3$		$\times$	$\times$	$\times$	$\times$	✓	$\times$
${ℛ}_2^y$	$+U\mathscr{H}(-k_x,k_y)U^{†}$	$i{\sigma }_2{\tau }_0$		✓	$\times$	$\times$	$\times$	$\times$	$\times$
${ℛ}_2^z$	$+U\mathscr{H}(-k_x,-k_y)U^{†}$	${\sigma }_3{\tau }_3$		$\times$	$\times$	$\times$	$\times$	$\times$	$\times$
${ℳ}_x$	$+U\mathscr{H}(-k_x,k_y)U^{†}$	${\sigma }_1{\tau }_0$		$\times$	$\times$	$\times$	$\times$	✓	✓
${ℳ}_y$	$+U\mathscr{H}(k_x,-k_y)U^{†}$	$i{\sigma }_2{\tau }_3$		✓	✓	$\times$	$\times$	$\times$	$\times$
${ℳ}_z$	$+U\mathscr{H}(k_x,k_y)U^{†}$	${\sigma }_3{\tau }_0$		$\times$	$\times$	$\times$	$\times$	$\times$	$\times$
$𝒯\cdot$
$\cdot {ℛ}_2^x$		${\sigma }_3{\tau }_3$		✓	✓	$\times$	$\times$	$\times$	$\times$
$\cdot {ℛ}_2^y$		${\sigma }_0{\tau }_0$		$\times$	$\times$	$\times$	$\times$	✓	✓
$\cdot {ℛ}_2^z$		${\sigma }_1{\tau }_3$		✓	$\times$	✓	$\times$	✓	$\times$
$\cdot {ℳ}_x$		${\sigma }_3{\tau }_0$		✓	$\times$	$\times$	$\times$	$\times$	$\times$
$\cdot {ℳ}_y$		${\sigma }_0{\tau }_3$		$\times$	$\times$	$\times$	$\times$	✓	$\times$
$\cdot {ℳ}_z$		${\sigma }_1{\tau }_0$		✓	✓	✓	✓	✓	✓
$𝒞\cdot$
$\cdot {ℳ}_x$		${\sigma }_0{\tau }_1$		✓	$\times$	$\times$	$\times$	$\times$	$\times$
$\cdot {ℳ}_y$		$i{\sigma }_3{\tau }_2$		$\times$	$\times$	$\times$	$\times$	✓	$\times$
$\cdot {ℳ}_z$		$i{\sigma }_2{\tau }_1$		✓	✓	✓	✓	✓	✓

Equations

$$\mathcal{H}(\mathbf{k})=\left(\begin{array}{cc}{\hat{h}}_{\mathbf{k}}& {\hat{\Delta }}_{\mathbf{k}}\\ {\hat{\Delta }}_{\mathbf{k}}^{\ast }& -{\hat{h}}_{\mathbf{k}}^{\ast }\end{array}\right),$$

$$\begin{split}\hat{\Delta}=\left(\begin{matrix}\Delta_{\uparrow\uparrow}&-\Delta_{\uparrow\downarrow}\\ \Delta_{\downarrow\uparrow}&-\Delta_{\downarrow\downarrow}\end{matrix}\right),\quad\begin{array}[]{l}{\Delta}_{\uparrow\uparrow}=t_{0}(\Delta_{\downarrow\downarrow}^{*}=\!\sin k_{x}-i\sin k_{y})\,,\\ \Delta_{\downarrow\uparrow}=\!-\Delta_{\uparrow\downarrow}=\!s_{x}\cos k_{x}+s_{y}\cos k_{y}\,,\end{array}\end{split}$$

$${\hat{h}}_{\mathbf{k}}=t_0(1-\cos k_x-\cos k_y){\sigma }_0+\mathbf{b}\cdot \mathbf{\sigma },$$

$${\mathcal{G}}_{(0,0)}=D_{4h}\times \{1,\mathcal{T},\mathcal{P},\mathcal{C}\},$$

$${\mathcal{G}}_{(0,0)}={\tilde{\mathcal{G}}}_{(0,0)}\times \{1,\mathcal{P}\},$$

$$\mathcal{H}(\mathbf{k})$$

$$\gamma_{\alpha}=\sum_{\bm{r}}\sum_{i=1}^{4}\left[\psi_{\alpha,\bm{r}}\right]^{*}_{i}\left[\hat{\Psi}_{\bm{r}}\right]_{i}\,,\quad\alpha\in\big{\{}\text{A},\text{B}\big{\}}\,,$$

$${\gamma }_{\alpha }={\gamma }_{\alpha }^{†},\{{\gamma }_{\alpha },{\gamma }_{\beta }\}=2{\delta }_{\alpha \beta }.$$

$${\Phi }_{\alpha }={\oint }_0^T{\left[{\mathcal{A}}_t\right]}_{\alpha \alpha }\mathrm{d}t,\alpha \in \{\text{A, B}\},$$

$$\left\{\begin{array}[]{c}{\gamma_{\text{A}}}\quad\to\quad-\;{\gamma_{\text{A}}}\\ {\gamma_{\text{B}}}\quad\to\quad-\;{\gamma_{\text{B}}}\end{array}\right.,$$

$${\tilde{\mathcal{G}}}_{(\mathbf{b},\mathbf{s})}=N_{(\mathbf{b},\mathbf{s})}\oplus \mathcal{T}\left(G_{(\mathbf{b},\mathbf{s})}-N_{(\mathbf{b},\mathbf{s})}\right),$$

$$D_{2h}\;(mmm)\supset\left\{\begin{array}[]{lll}C_{2v}\;(mm2)&\supset C_{1h}\;(m)&\supset C_{1}\;(1)\\ C_{2h}\;(2/m)&\supset C_{1h}\;(m)&\supset C_{1}\;(1)\\ C_{2h}\;(2/m)&\supset C_{i}\;(\bar{1})&\supset C_{1}\;(1)\end{array}\right.$$

$$\left\{\Gamma ,\mathcal{H}(k_{\parallel }=k_0,k_{\perp })\right\}=0,\forall k_{\perp },$$

$$W_{\text{1D}}\left(k_{\parallel }\right)=-\frac{1}{4\pi i}\oint \mathrm{d}{k}_{\perp }\mathrm{Tr}\left\{\Gamma \mathcal{H}(\mathbf{k}{)}^{-1}{\partial }_{k_{\perp }}\mathcal{H}(\mathbf{k})\right\}.$$

$$\begin{array}{r}{U}_{\Gamma }^{†}\Gamma U_{\Gamma }=\mathrm{diag}(-1,-1,1,1)\Rightarrow U_{\Gamma }^{†}\mathcal{H}(\mathbf{k})U_{\Gamma }=\left(\begin{array}{cc}0& f_{\mathbf{k}}\\ f_{\mathbf{k}}^{†}& 0\end{array}\right),\end{array}$$

$$W_{\text{1D}}\left(k_{\parallel }\right)=\frac{1}{2\pi }\oint \mathrm{d}{k}_{\perp }{\partial }_{k_{\perp }}\arg \left(\det f_{\mathbf{k}}\right).$$

$$\begin{array}{rl}{\mathcal{CM}}_z:& \left\{\Gamma ,\mathcal{H}(k_x,k_y)\right\}=0,\\ & \Gamma ={\tau }_1{\sigma }_2,U_{\Gamma }=\frac{1}{\sqrt{2}}\left(\begin{array}{cccc}0& i& 0& -i\\ -i& 0& i& 0\\ 1& 0& 1& 0\\ 0& 1& 0& 1\end{array}\right).\end{array}$$

$$\begin{split}f_{\bm{k}}=&\left(\begin{array}[]{cc}-\xi_{\bm{k}}&i\,\Big{(}b-\Delta_{\bm{k},\uparrow\uparrow}^{*}\Big{)}\\ -i\,\Big{(}b+\Delta_{\bm{k},\uparrow\uparrow}\Big{)}&-\xi_{\bm{k}}\end{array}\right)\\ &\Rightarrow\quad\det f_{\bm{k}}=\xi_{\bm{k}}^{2}+|\Delta_{\bm{k},\uparrow\uparrow}|^{2}-2ib\sin k_{y}-b^{2}\,.\end{split}$$

$$\begin{array}{rl}{\mathcal{CM}}_y:& \left\{\Gamma ,\mathcal{H}(k_x,k_y=0)\right\}=0,\\ & \Gamma ={\tau }_2{\sigma }_3,U_{\Gamma }=\frac{1}{\sqrt{2}}\left(\begin{array}{cccc}i& 0& -i& 0\\ 0& -i& 0& i\\ 1& 0& 1& 0\\ 0& 1& 0& 1\end{array}\right).\end{array}$$

$$\begin{split}\det f_{(k_{x},0)}=&\left|\begin{array}[]{cc}e^{-ik_{x}}&b\\ b&e^{-ik_{x}}\\ \end{array}\right|=e^{-2ik_{x}}-b^{2}\\ &\stackrel{{\scriptstyle\text{\lx@cref{creftype~refnum}{W1Dh}}}}{{\Longrightarrow}}\quad W_{\text{1D}}(k_{y}=0)=-2\,.\end{split}$$

$$\begin{array}{rl}{\mathcal{CM}}_x:& \left\{\Gamma ,\mathcal{H}(k_x=0,k_y)\right\}=0,\\ & \Gamma ={\tau }_1{\sigma }_0,U_{\Gamma }=\frac{1}{\sqrt{2}}\left(\begin{array}{cccc}-1& 0& 1& 0\\ 0& -1& 0& 1\\ 1& 0& 1& 0\\ 0& 1& 0& 1\end{array}\right).\end{array}$$

$$c_{\mathbf{k}\sigma }=\frac{1}{N}\sum _{\mathbf{r}}{e}^{i\mathbf{k}\cdot \mathbf{r}}{c}_{\mathbf{r}\sigma },$$

$${\hat{\Psi }}_{\mathbf{r}}:={\left(\begin{array}{cccc}{c}_{\mathbf{r}\uparrow }& c_{\mathbf{r}\downarrow }& c_{\mathbf{r}\uparrow }^{†}& -c_{\mathbf{r}\downarrow }^{†}\end{array}\right)}^T.$$

$$\mathcal{P}\hat{\Psi }\equiv U_{\mathcal{P}}{\hat{\Psi }}^{\ast }=\hat{\Psi },U_{\mathcal{P}}={\tau }_1{\sigma }_3.$$

$$\mathcal{H}|{\varphi }_n⟩=E_n|{\varphi }_n⟩,$$

$$\begin{split}\ket{\psi_{\text{A}}}&=\ket{\phi_{1}}+\ket{\phi_{2}}\,,\\ \ket{\psi_{\text{B}}}&=i\,\big{(}\ket{\phi_{1}}-\ket{\phi_{2}}\big{)}\,,\\ \end{split}$$

$$\mathcal{P}|{\psi }_{\alpha }⟩={\tilde{U}}_{\mathcal{P}}{|{\psi }_{\alpha }⟩}^{\ast }=|{\psi }_{\alpha }⟩.$$

$$\mathcal{H}\ket{u^{j}_{\bm{k}}}=E_{\bm{k}}^{j}\ket{u^{j}_{\bm{k}}},\qquad U_{\bm{k}}=\left(\begin{array}[]{c|c}\ket{u^{1}_{\bm{k}}}&\ket{u^{2}_{\bm{k}}}\end{array}\right).$$

$$P_{\mathbf{k}}:=\sum _{i=1}^{N_{\text{occ}}}|u_{\mathbf{k}}^i⟩⟨u_{\mathbf{k}}^i|=U_{\mathbf{k}}{U}_{\mathbf{k}}^{†}.$$

$${\Delta }_x=(\Delta ,0),\Delta =\frac{2\pi }{Na}.$$

Full text

Standard $d$-dimensional topological insulators and superconductors have a gapped bulk spectrum and exhibit conducting surface states in $(d-1)$ dimensions Kane and Mele (2005); Bernevig et al. (2006); Fu and Kane (2006); Fu et al. (2007); Fu and Kane (2007); Hasan and Kane (2010); Qi and Zhang (2011); Chiu et al. (2016). The analysis of non-interacting electrons in materials possessing fundamental symmetries – time-reversal ($\mathcal{T}$), particle-hole ($\mathcal{P}$) and/or chiral symmetry ($\mathcal{C}$) – led to the initial classifications of such topological phases Schnyder et al. (2008, 2009); Kitaev et al. (2009). Later on, possible topological phases protected by discrete crystalline symmetries were studied Fu (2011); Hsieh et al. (2012); Slager et al. (2012); Chiu et al. (2013); Morimoto and Furusaki (2013); Shiozaki and Sato (2014); Benalcazar et al. (2014); Alexandradinata et al. (2014a); Shiozaki et al. (2016); Kruthoff et al. (2017); Shiozaki et al. (2017); Thorngren and Else (2018); Khalaf et al. (2018) and, lastly, another class of exotic non-interacting topological phases was discovered, namely the second-order topological insulators (SOTIs) and superconductors (SOTSs) Benalcazar et al. (2017a); Song et al. (2017); Schindler et al. (2018a); Geier et al. (2018); Khalaf (2018); Ahn and Yang (2019); van Miert and Ortix (2018).

Unlike standard gapped topological materials (that may be termed ”first-order”), such systems have insulating $(d-1)$-dimensional surfaces, but host topologically-protected boundary modes in $(d-2)$ dimensions. In particular, for $d=3$ this means the material exhibits one-dimensional (1D) states propagating along its hinges Benalcazar et al. (2017b); Song et al. (2017); Langbehn et al. (2017); van Miert and Ortix (2018); Schindler et al. (2018a); Wang et al. (2018a); Călugăru et al. (2019); Ghorashi et al. (2019); Agarwala et al. (2020) as detected in Bi(111) Schindler et al. (2018b). Accordingly, SOTIs and SOTSs in $d=2$ host localized ($d=0$) corner modes and have been discussed in both static Benalcazar et al. (2014, 2017a, 2017b); Langbehn et al. (2017); Wang et al. (2018a); van Miert and Ortix (2018); Călugăru et al. (2019); Volpez et al. (2019); Agarwala et al. (2020) and Floquet-driven systems Franca et al. (2018); Rodriguez-Vega et al. (2019); Huang and Liu (2020). Second-order topological phases have been experimentally realized in artificial systems, such as mechanical metamaterials Serra-Garcia et al. (2018), microwave Peterson et al. (2018) and topolectrical circuits Imhof et al. (2018), and further implementations in materials have been proposed Schindler et al. (2018a); Ezawa (2018a); Wang et al. (2018a); Hsu et al. (2018).

One of the motivations to seek for topological materials is their potential implementation in devices that would exploit the topology-ensured protection from decoherence. Of particular interest is topological quantum computation, considered the cornerstone of quantum technology, in which logical operations are done by means of braiding (double exchange) of non-Abelian quasiparticle excitations Kitaev (2001, 2006); Nayak et al. (2008); Alicea (2012); Sarma et al. (2015). Majorana (zero-energy) excitations are known to be the only non-Abelian quasiparticles that can appear in the absence of electron-electron interactions Fidkowski and Kitaev (2010, 2011) and were initially predicted to appear at the edges of a 1D spinless superconductor Kitaev (2001). Almost a decade later, it was demonstrated that this exotic topological phase could be achieved in a semiconductor nanowire with strong spin-orbit coupling placed in an external magnetic field, if superconductivity is induced by proximity effect  Lutchyn et al. (2010); Oreg et al. (2010). Following these theoretical proposals, such devices were experimentally realized by several groups Mourik et al. (2012); Das et al. (2012); Deng et al. (2012); Finck et al. (2013); Churchill et al. (2013); Deng et al. (2016) and are currently among the most promising building block for topological quantum computing devices. However, despite the existing proposals of braiding Majorana end modes using nanowires Alicea et al. (2011); van Heck et al. (2012), identifying the topological qubit remains extremely challenging due to the engineering complexity of manipulating the topological phase in the wire junctions without losing the superconducting state. This motivates the search for other approaches to realize a Majorana qubit in potentially feasible condensed matter systems Zhu (2018); Yan et al. (2018); Wang et al. (2018b); Liu et al. (2018); Schrade and Fu (2018).

In this paper we propose a device based on a two-dimensional (2D) SOTS which hosts two Majorana modes localized at two corners of a square-shaped crystalline sample (see Fig. 1a). We call them Majorana corner states (MCSs) in the following. Within our model, the corner-localization of the MCSs in the device can be tuned using an in-plane magnetic field and proximity-induced spin-singlet superconductivity. In this setup, the braiding of the Majorana quasiparticles can be realized in a straightforward fashion (schematically illustrated in Fig. 1b), and we show that the associated fractional statistics of the proposed topological Majorana qubit is observed.

SOTS model

In $\bm{k}$-space, our proposed model can be formulated with a minimal four-band Bogoliubov-de Gennes Hamiltonian

[TABLE]

written in the Nambu basis $\hat{\Psi}_{\bm{k}}:=\left(\begin{matrix}c_{\bm{k}\uparrow}&c_{\bm{k}\downarrow}&c_{-\bm{k}\uparrow}^{\dagger}&-c_{-\bm{k}\downarrow}^{\dagger}\end{matrix}\right)^{T}$, with $c^{(\dagger)}_{\bm{k}\,s}$ being the annihilation (creation) operator for an electron with a 2D momentum $\bm{k}$ and spin $s\in\{\uparrow,\downarrow\}$. The gap function considered is (with $\bm{k}$-labels suppressed)

[TABLE]

while the diagonal blocks are given by

[TABLE]

where $\bm{\sigma}=(\sigma_{1},\sigma_{2},\sigma_{3})$ and $\sigma_{0}$ are the three Pauli matrices and identity in spin space, respectively. The vectors $\bm{b}\equiv(b_{x},b_{y})$ and $\bm{s}\equiv(s_{x},s_{y})$ parametrize the in-plane magnetic field and the spin-singlet pairing amplitudes, respectively, with magnitudes smaller than $t_{0}$, taken as unit in the following.

For $\bm{b}=\bm{s}=\bm{0}$, the Hamiltonian $\mathcal{H}_{(\bm{0},\bm{0})}$ describes electrons with nearest-neighbor hopping on a square lattice, subject to $p$-wave pairing (for clarity, we use here the vectors $\bm{b}$ and $\bm{s}$ as subscripts to refer to a specific parameter configuration). $\mathcal{H}_{(\bm{0},\bm{0})}$ possesses all three non-crystalline symmetries ($\mathcal{T}$, $\mathcal{P}$ and $\mathcal{C}$) and falls into the symmetry class DIII Schnyder et al. (2008, 2009); Kitaev et al. (2009), supporting helical Majorana edge states (Fig. 2a). Including crystalline symmetries, the total point symmetry group is

[TABLE]

where $\mathbb{1}$ represents the identity. The main symmetries elements of $\mathcal{G}_{(\bm{0},\bm{0})}$ are listed in Table 1, alongside with the corresponding representations. Since $\mathcal{P}$ is always present in our model, $\mathcal{G}_{(\bm{0},\bm{0})}$ can be rewritten as:

[TABLE]

where the magnetic group for $\bm{b}=\bm{s}=0$ reads $\tilde{\mathcal{G}}_{(\bm{0},\bm{0})}=D_{4h}\oplus\mathcal{T}D_{4h}$ (or $4/mmm1^{\prime}$). The magnetic group is reduced for other values of the parameters $\bm{b}$ and $\bm{s}$, as explained below, and detailed in Appendix A.

A finite in-plane magnetic field $\bm{b}$ conserves $\mathcal{P}$ for any direction, but breaks $\mathcal{T}$. The reduced magnetic group is $\tilde{\mathcal{G}}_{(\bm{b},\bm{0})}=C_{2h}\oplus\mathcal{T}(D_{2h}-C_{2h})$ 111Here, the set $(D_{2h}-C_{2h})$ represents the relative complement of $C_{2h}$ with respect to $D_{2h}$. It contains the elements of $D_{2h}$ which are not part of $C_{2h}$ and it is not a group by itself. if $\bm{b}$ is parallel to $\hat{\bm{x}}$, $\hat{\bm{y}}$ or $\hat{\bm{x}}\pm\hat{\bm{y}}$, and $C_{i}\oplus\mathcal{T}(C_{2h}-C_{i})$ otherwise (i.e. $m^{\prime}m^{\prime}m$ and $2^{\prime}/m^{\prime}$, respectively). The magnetic field hybridizes the helical boundary modes unless they are protected by symmetry: as illustrated in Fig. 2b, for $\bm{b}\parallel\hat{\bm{x}}$ only the $x$ edges remain metallic protected by $\mathcal{CM}_{y}\in\mathcal{G}_{(\hat{\bm{x}},\bm{0})}$, a product of $\mathcal{C}$ and the mirror reflection $\mathcal{M}_{y}$ which maps $(x,y,z)\mapsto(x,-y,z)$. Analogously, for $\bm{b}\parallel\hat{\bm{y}}$ the gapless states at the $y$ edges are protected by $\mathcal{CM}_{x}$ (see Appendix B for a proof). Both edges are gapped for generic orientations $\bm{b}\nparallel\hat{\bm{x}},\hat{\bm{y}}$, but, if $\bm{b}\parallel(\hat{\bm{x}}\pm\hat{\bm{y}})$, the mirror-symmetric corners bisected by $\bm{b}$ will host zero-energy modes Geier et al. (2018); Ezawa (2018b).

The spin-singlet pairing controlled by $\bm{s}$ ($|\bm{s}|\neq 0$) reduces the symmetry further and in the generic case $\tilde{\mathcal{G}}_{(\bm{b},\bm{s})}=C_{1}\oplus\mathcal{T}(C_{s}-C_{1})$ (or $m^{\prime}$). Thus, the (minimal) total group of the total Hamiltonian is $\mathcal{G}_{(\bm{b},\bm{s})}=\{\mathbb{1},\mathcal{TM}_{z}\}\times\{\mathbb{1},\mathcal{P}\}$, including two particularly important elements, $\mathcal{P}$ and the effective chiral symmetry $\widetilde{\mathcal{C}}=\mathcal{CM}_{z}=\mathcal{PTM}_{z}$, which transform the Hamiltonian in the following way:

[TABLE]

where the Pauli matrix $\tau_{1}$ acts on the particle-hole space.

Majorana corner states

On a finite square-shaped lattice, the $\mathcal{P}$-symmetric Hamiltonian hosts two degenerate zero-energy excitations $\ket{\psi_{\alpha}}$, with $\alpha\in\{\mathrm{A},\mathrm{B}\}$, whose features in the second-order topological phase are displayed in Fig. 1a for $\bm{b}=(-0.3,0)$ and $\bm{s}=(0,0.3)$. In this regime, although the magnetic field $\bm{b}\parallel\hat{\bm{x}}$ would allow for $\mathcal{CM}_{y}$-protected boundary modes at the $x$ edges, the spin-singlet pairing controlled by $\bm{s}\parallel\hat{\bm{y}}$ gaps the $x$ edges by breaking this symmetry (Fig. 2c) and enforces the zero-energy modes to localize at two adjacent corners of the device (as shown in the inset of Fig. 1a). More precisely, the pairing term $\Delta_{\bm{k},\downarrow\uparrow}$ reduces $\tilde{\mathcal{G}}_{(\hat{\bm{x}},\bm{0})}=C_{2h}\oplus\mathcal{T}(D_{2h}-C_{2h})$ to $\tilde{\mathcal{G}}_{(\hat{\bm{x}},\hat{\bm{y}})}=C_{s}\oplus\mathcal{T}(C_{2v}-C_{s})$ (i.e. $m^{\prime}m^{\prime}m$ to $m^{\prime}m2^{\prime}$).

The operators corresponding to the MCSs can be written as

[TABLE]

where $[\hat{\Psi}_{\bm{r}}]_{i}$ is the $i$-th component of the Fourier-transformed Nambu basis for the site at position $\bm{r}$ and $\left[\psi_{\alpha,\bm{r}}\right]_{i}$ the $i$-th component of the corresponding sector $\ket{\psi_{\alpha,\bm{r}}}$ of the MCS $\ket{\psi_{\alpha}}$. Because $\ket{\psi_{\alpha}}$ are self-conjugate under $\mathcal{P}$ ($\mathcal{P}\ket{\psi_{\alpha,\bm{r}}}\equiv U_{\mathcal{P}}\ket{\psi_{\alpha,\bm{r}}}^{*}=\ket{\psi_{\alpha,\bm{r}}}$), these operators satisfy the following relations (for details see Appendix C)

[TABLE]

The energies of the midgap states go to zero exponentially with the system size (data shown in Appendix C), are robust against $\mathcal{P}$-preserving disorder (tested numerically, data not shown) and are well separated from the bulk states by an energy gap (emphasized in Fig. 1a). The localization of the boundary modes is strongly dependent on the Hamiltonian parameters $\bm{b}$ and $\bm{s}$, whereby their energy remains zero for all the considered configurations, in which the device is in a SOTS phase. This degenerate ground state manifold is crucial to realize the braiding of the topological MCSs.

Computing topological invariants

We can assess the topological features of our device using the Wilson loop operator (WLO). Several technical aspects are explained in Appendix D.

The WLO can be adiabatically connected to the Hamiltonian at the boundary Blount (1962); King-Smith and Vanderbilt (1993); Fidkowski et al. (2011); Taherinejad et al. (2014) and can thus characterize the surface topology of materials, such as $\mathbb{Z}_{2}$ topological (crystalline or not) insulators in 2D and 3D Fu and Kane (2006); Soluyanov and Vanderbilt (2011a, b); Yu et al. (2011); Taherinejad et al. (2014); Alexandradinata et al. (2014b); Alexandradinata and Bernevig (2016); Gresch et al. (2017) or $\mathbb{Z}$ Chern insulators Coh and Vanderbilt (2009); Qi (2011); Alexandradinata et al. (2014b); Gresch et al. (2017).

Moreover, the nested WLO introduced in Ref. Benalcazar et al. (2017a) has been used to address the physics at ”the boundary of the boundary” Schindler et al. (2018a); Franca et al. (2018); Xie et al. (2018). In our 2D model, using the nested WLOs we obtain two $\mathbb{Z}_{2}$ corner topological invariants $\phi_{xy},\phi_{yx}\in\{0,\pi\}$, which are quantized by the effective chiral symmetry introduced in Eq. (6b) (confirmed by numerical simulations). The device is in a SOTS phase whenever either of $\{\phi_{xy},\phi_{yx}\}$ is equal to $\pi$; for example, in the case $\bm{b}=(-0.3,0)$ and $\bm{s}=(0,0.3)$, shown in Fig. 1a, we have $(\phi_{xy},\phi_{yx})=(\pi,0)$.

Braiding process

We perform the adiabatic braiding of the MCSs by rotating the magnetic field ($\bm{b}$) and adjusting the spin-singlet pairing ($\bm{s}$). Let $T$ be the time period of the cyclic process and $t\in[0,T]$ the parameter which controls the process, such that we have periodic boundary conditions $\mathcal{H}^{(t=0)}=\mathcal{H}^{(t=T)}$. The evolution of the parameters $\bm{b}(t)$ and $\bm{s}(t)$ over a full period is shown in Fig. 3, where the alternations of the second-order topological states characterized by $(\phi_{xy},\phi_{yx})=(\pi,0)$ and $(0,\pi)$ are also highlighted. Along this path in parameter space, the bulk and the edges remain gapped, while the energy of the degenerate MCSs is unaffected – they remain well-separated from the bulk states by an energy gap (feature displayed in the inset of Fig. 3). The MCSs thus do not mix with the bulk states during the cyclic process, which we illustrate in Fig. 4 by computing the spatial probability density of $\ket{\,\psi_{\alpha}(t)\,}$ in the proposed device for several $t$ ($\alpha\in\{\text{A},\text{B}\}$).

Statistical phase

A fermionic many-body wavefunction picks a factor $(-1)$ upon exchanging two fermions and it reverts to the initial wavefunction if the exchange of any other two is performed. In contradistinction, in a system with $2M$ Majorana quasiparticles ($M>1$) two subsequent exchanges (say, $\gamma_{1}\leftrightarrow\gamma_{2}$ and $\gamma_{2}\leftrightarrow\gamma_{3}$) generally give a different state (non-Abelian braid group) Kitaev (2006); Nayak et al. (2008); Alicea (2012); Sarma et al. (2015). In particular, if the system contains a single pair of Majoranas ($M=1$), such a double exchange (braiding) reveals the fractional statistics of the two quasiparticles: each acquires a non-trivial statistical phase $\Phi=\pi$.

The device we propose here realizes the braiding of two midgap states $\psi_{\alpha}$ [with corresponding operators $\gamma_{\alpha}$, $\alpha\in\{\text{A},\text{B}\}$, see Eq. (7)] and we show below their statistics is indeed that of Majorana quasiparticles (fractional, $M=1$). We note, however, that in the proposed setup the single exchange $\text{A}\leftrightarrow\text{B}$ is not feasible with the parameter cycle considered in Fig. 3. Although at $t=T/2$ the centers of the wavefunctions $\psi_{\text{A}}$, $\psi_{\text{B}}$ are swapped (see Fig. 4), the Hamiltonian after a half-period is different from the initial one. Consequently, a statistical phase between $\psi_{\text{A}}$ ($\psi_{\text{B}}$) at $t=0$ and $\psi_{\text{B}}$ ($\psi_{\text{A}}$) at $t=T/2$ would not be well-defined nor robust.

To express the braiding phase accumulated by each MCS within one full cycle (i.e. after a double exchange), we use the Berry phase defined as

[TABLE]

where the Berry connection matrix has elements $\left[\mathcal{A}_{\;t}\right]_{\alpha\beta}=i\bra{\,\psi_{\alpha}(t)\,}\partial_{t}\ket{\,\psi_{\beta}(t)\,}$ (more details are given in Appendix E). The quantity in Eq. (9) is gauge-invariant because (i) the Hamiltonian parameter configurations at $t\in\{0,T\}$ are identical and (ii) the matrix $\mathcal{A}_{\;t}$ is diagonal $\forall t\in[0,T]$, as the two MCSs are kept far apart from each other and do not mix during the course of the adiabatic cycle. The practically null spatial overlap of the two states A and B was checked numerically and is illustrated in Fig. 1b and Fig. 4.

With Eq. (9), we obtain non-trivial statistical phases $\Phi_{\text{A}}=\Phi_{\text{B}}=\pi$, which means

[TABLE]

and, together with the operator algebra (8), proves the Majorana nature of the quasiparticles. As confirmed numerically, the quantization of the phases is guaranteed by $\mathcal{P}$, Eq. (6a).

However, the relative phase described here would be an inaccessible experimental quantity within the proposed protocol, since the two MCSs necessarily have different fermion parities. Nonetheless, we believe it might be possible to extend our scheme as to construct a system hosting $M>1$ Majorana pairs (for example, by distributing $M$ devices as a chain or a two-dimensional array). In such a setup, we would expect that the fractional statistics (10) of our MCSs ($M=1$) would generalize to a non-Abelian braid group.

Schematic device design

Finally, we briefly describe a possible structure for the proposed device that would consist of three square-shaped thin-film superconductors stacked as in Fig. 5a. The superconducting state in our model ($\hat{\Delta}_{\bm{k}}$) can be realized if we induce spin-singlet pairing in the $p$-wave film through proximity effect, by placing it between an extended $s$-wave superconductor and a $d_{x^{2}-y^{2}}$-wave superconductor, with even-parity gap functions $\psi_{s/d}=\eta_{s/d}(\cos k_{x}\pm\cos k_{y})$. In terms of the $\bm{s}$-parameter, the gap magnitudes can be written as $\eta_{s/d}=(s_{x}\pm s_{y})/2$. The relative phases between the superconductors can be fixed by Josephson contacts, as explained in the following.

First, the order parameter $\eta_{p}$ corresponding to the $p$-wave couples only in second order to the spin-singlet ones, which locks the $p$-$s$ and $p$-$d$ relative phases to either $\pi/2$ or $-\pi/2$ (in order to minimize the Josephson energy, in lowest order being typically $E_{\text{J}}\propto+\cos 2\phi_{p-s/d}$). Indeed, first-order coupling terms as $(\eta_{s/d}^{*}\,\eta_{p}+\text{c.c.})$ are assumed to be absent or extremely weak, due to the mismatch of parity. Note that the $p$-wave superconducting state, $d$-vector $\bm{d}_{\bm{k}}=\hat{\bm{x}}\sin k_{y}-\hat{\bm{y}}\sin k_{x}$, belongs to irreducible representation $A_{2u}$, while the order parameters $\psi_{s}$ and $\psi_{d}$ belong to the even-parity representations $A_{1g}$ and $B_{1g}$, respectively. Thus, the phase coupling between the odd-parity and even-parity order parameter would be dominated by the second order coupling of the form $(\eta_{s/d}^{*2}\,\eta_{p}^{2}+\text{c.c.})$. The appropriate phase factor $e^{i\pi/2}$ is already considered in Eq. (1), where the configuration ”$\eta_{s}+\eta_{d}+i\eta_{p}$” is presented.

Second, the relative phase between the spin-singlet condensates can be controlled using alternatively two direct Josephson contacts (see Fig. 5a). The gap functions $\psi_{s}$ and $\psi_{d}$ exhibit in the $xy$ plane positive and negative sectors, as schematically represented in Fig. 5b, and the Josephson contacts are set up such that they link equal- or different-sign sectors, respectively. Switching the active contact increases by $\pi$ not only the relative $s$-$d$ phase, but consequently also that singlet-triplet phase, where the Josephson coupling is weaker. As emphasized in Fig. 5c, during the proposed adiabatic cycle such an operation changes the $p$-$s$ ($p$-$d$) phase by $\pi$ when $|\eta_{s}|<(>)|\eta_{d}|$.

Manipulating the spin-singlet gap magnitudes $\eta_{s/d}$ during the proposed cycle remains a challenging task. At a conceptual level, two strategies could be taken under consideration: using tunable tunneling Josephson contacts between the $p$-wave and the spin-singlet superconductors or modifying the order parameters amplitudes, $\eta_{s}$ and $\eta_{d}$ individually by local heating or other means. In fact, small changes to the proposed functions $\eta_{s/d}(t)$, that might come within the experimental implementation, are not expected to change the topological features of the cycle. The latter argument applies as well to the in-plane magnetic field $\bm{b}(t)$. Moreover, we assume this to be weak enough such that in the superconducting thin films the interference with the Josephson effect as well as the depairing effects are negligible.

*Note. The physical system considered in the previous arXiv version (1904.07822v2) differs from the current one, however, their topological and symmetry features are closely related. *

Acknowledgements.

T.E.P. would like to thank J.L. Lado, T. Kawakami and K. Viebahn for helpful discussions. T.E.P. and M.S. are grateful for the financial support from the Swiss National Science Foundation (SNSF) through Division II (No. 163186 and 184739). A.A.S. acknowledges the support of Microsoft Research, SNSF NCCR MARVEL and QSIT programs, and the SNSF Professorship grant.

Appendix A Symmetries of the model Hamiltonian

The group $\tilde{\mathcal{G}}_{(\bm{0},\bm{0})}$ of Eq. 5 is a Heesch-Schubnikov magnetic group of type II; a magnetic field ($|\bm{b}|\neq 0$) breaks $\mathcal{T}$ and reduces it to a magnetic group of type III. In the general case (possibly $|\bm{s}|\neq 0$), we can write

[TABLE]

where ${G}_{(\bm{b},\bm{s})}$ is a subgroup of $D_{2h}$ and ${N}_{(\bm{b},\bm{s})}$ a halving (and normal) subgroup of ${G}_{(\bm{b},\bm{s})}$. The following relations between the subgroups of $D_{2h}$ hold:

[TABLE]

The groups ${G}_{(\bm{b},\bm{s})}$ and ${N}_{(\bm{b},\bm{s})}$ are given explicitly in Table 1 for certain values of $\bm{b}$ and $\bm{s}$; in those six cases, the magnetic point groups $\tilde{\mathcal{G}}_{(\bm{b},\bm{s})}$ in international notation are $m^{\prime}m^{\prime}m$,   $m^{\prime}m2^{\prime}$,   $2^{\prime}/m^{\prime}$,   $m^{\prime}$,   $m^{\prime}m^{\prime}m$,   $m^{\prime}m2^{\prime}$. In the Table, we denote by $\{\tau_{j}\}$, with $j\in\{0,1,2,3\}$, the set of Pauli matrices acting in particle-hole space.

Appendix B Edge states and their symmetry protection

The two counter-propagating edge modes of opposite spins hosted by the system described by $\mathcal{H}_{(\bm{0},\bm{0})}$ (helical state, depicted in Fig. 2a) continue to exist on certain edges even when the magnetic field is switched on, while they aquire a gap on other edges. We prove their symmetry protection by introducing a one-dimensional (1D) winding number.

Consider a 2D system with a finite number of unit cells in $k_{\perp}$ direction, but with periodic boundary conditions in the other direction $k_{\parallel}$. Provided the bulk Hamiltonian $\mathcal{H}(\bm{k})\equiv\mathcal{H}(k_{x},k_{y})$ has an effective chiral symmetry (with representation $\Gamma$) for a certain $k_{\parallel}$,

[TABLE]

then the topology of the point $k_{0}$ is characterized by the following 1D winding number Sato et al. (2011):

[TABLE]

The chiral symmetry (13) implies that the matrix diagonalizing $\Gamma$ off-diagonalizes the Hamiltonian:

[TABLE]

such that Eq. 14 can be simplified:

[TABLE]

This expresses the number of times the phase of $\det f_{\bm{k}}$ winds when the parameter $k_{\perp}$ completes a full period.

System infinite in the $y$-direction. We consider the case with magnetic field $\bm{b}=(b,0)\parallel\hat{\bm{x}}$, for which there are boundary modes only at the $x$ edges, dispersing as $\pm k_{y}$ at low energies (see Fig. 6a-b). For our model, the chiral symmetry $\mathcal{CM}_{z}$ is present at any momentum $\bm{k}$:

[TABLE]

Applying $U_{\Gamma}$ of Eq. 17 in Eq. 15 and using ${\Delta}_{\bm{k},\uparrow\uparrow}=\Delta_{\bm{k},\downarrow\downarrow}^{*}$, we obtain the off-diagonal block:

[TABLE]

According to Eq. 16, the winding number at $k_{y}=0$ is trivial, since $\det f_{(k_{x},0)}=1-b^{2}$ is a constant. Therefore, the boundary modes are not protected by $\mathcal{CM}_{z}$.

We now consider the effective chiral symmetry (see Fig. 2b)

[TABLE]

This gives a nontrivial winding number as long as $|\bm{b}|<1$, which means $\mathcal{CM}_{y}$ protects the boundary modes:

[TABLE]

The 1D (edge) system has a $\mathbb{Z}$ topological classification Schnyder et al. (2008); Geier et al. (2018) and the topological invariant is precisely the computed winding number, which is equal to the number of states in one direction.

System infinite in the $x$-direction. The calculations are similar for perpendicular magnetic field $\bm{b}=(0,b)\parallel\hat{\bm{y}}$, for which we have edge states at the $y$ edges dispersing as $\pm k_{x}$ at low energies. Using the symmetry $\mathcal{CM}_{z}$, at $k_{x}=0$ we obtain $\det f_{(0,k_{y})}=1-b^{2}+2ib\sin k_{y}$, which gives as well a trivial winding number. The edge states are protected by $\mathcal{CM}_{x}$ (see Fig. 2a):

[TABLE]

This gives $\det f_{(0,k_{y})}=e^{2ik_{y}}-b^{2}$ and thus winding number $W_{\text{1D}}(k_{x}=0)=2$ if $|\bm{b}|<1$.

Appendix C Majorana corner states

Consider a squared-shaped sample of a square lattice with $N\times N$ sites. With the Fourier transform ($\bm{r}$ goes over the positions of all unit cells)

[TABLE]

we obtain the superconducting Hamiltonian $\mathcal{H}$ in the basis

[TABLE]

The Nambu basis $\hat{\Psi}$ ($\hat{\Psi}_{\bm{k}}$ or $\hat{\Psi}_{\bm{r}}$) is $\mathcal{P}$-symmetric, satisfying

[TABLE]

We diagonalize the real-space $(4N^{2})\times(4N^{2})$ Hamiltonian,

[TABLE]

and plot the spectrum $\{E_{n}\}$ in Fig. 1a. The system exhibits two zero-energy modes, which we denote $\ket{\phi_{1}}$ and $\ket{\phi_{2}}$. Their energy goes to zero exponentially with increasing $N$ (see Fig. 7a), scaling like $\exp(-N/N_{0})$, with $N_{0}$ a constant. Each of the zero modes is localized on two adjacent corners (Fig. 7b-c). However, the following two orthogonal linear combinations of the degenerate fermionic states $\{\ket{\phi_{1}},\ket{\phi_{2}}\}$

[TABLE]

are localized on separate corners (illustrated in Fig. 7d-e). Within our model, it is always possible to build such states $\{\ket{\psi_{\alpha}}\}$, with $\alpha\in\{\mathrm{A},\mathrm{B}\}$, which satisfy

[TABLE]

The Majorana corner states (MCSs) $\ket{\psi_{\alpha}}$ are thus self-conjugate under $\mathcal{P}$; together with Eq. 24, this property guarantees the operator algebra is that of Majorana excitations Elliott and Franz (2015) presented in Eq. 8.

Appendix D Wilson loop operators

We used the nested Wilson loop operator (WLO) to calculate the second-order topological invariants $\phi_{xy}$ and $\phi_{yx}$. This technique was introduced and explained in detail in Refs. Benalcazar et al. (2017a, b); here we only outline the essential steps.

The first-order WLO. First we diagonalize the four-band Bloch Hamiltonian and collect the $N_{\text{occ}}=2$ occupied states in a matrix:

[TABLE]

In this way, the projector to occupied space reads

[TABLE]

For the $N\times N$ square lattice with lattice constant $a$, we analyze the adiabatic evolution of a state $\ket{u^{m}_{\bm{k}}}$ along a closed path in reciprocal space, with $k_{y}$ constant, into a final state $\ket{u^{n}_{\bm{k}+(2\pi/a,0)}}$. We discretize the path into $N$ points $\bm{k}_{j}=\bm{k}+j\Delta_{x}$, where

[TABLE]

If we iteratively project the initial state until the final state is reached, the amplitude of this parallel transport process can be written as:

[TABLE]

where the symbol $\overline{\prod}$ stands for ordered product. The unitary matrix $\mathcal{W}_{x,\bm{k}}$ represents the WLO in the $x$-direction; introducing here Eq. 29, it reads

[TABLE]

where $\mathcal{A}^{x}_{\bm{k}}$ represents the $x$ component of the Berry connection constructed with the occupied Bloch eigenfunctions, \Big{(}\mathcal{A}^{x}_{\bm{k}}\Big{)}_{nm}=i\bra{u^{n}_{\bm{k}}}\partial_{k_{x}}\ket{u^{m}_{\bm{k}}}. We define the overlap matrices $O_{x,\bm{k}}$ with elements

[TABLE]

The last equality of Eq. 32 represents the formal definition of the WLO, although in practice the first equality therein is used (with a finite number $N$ of $\bm{k}$ points). The operator $\mathcal{W}_{x,\bm{k}}$ characterizes the $x$ boundary Coh and Vanderbilt (2009); Yu et al. (2011); Qi (2011) and, for this reason, we may refer to it as ”the edge Hamiltonian” 222The Hermitian edge Hamiltonian $H_{\text{ edge}}$ could be related to the unitary WLO through the mapping $\mathcal{W}_{x,\bm{k}}=\exp\left(iH_{\text{ edge}}\right)$.. This becomes clear if one proceeds to obtaining the spectrum of the WLO by solving

[TABLE]

The eigenvalues $\lambda_{x}^{n}$ are related to the surface charge polarization Blount (1962); King-Smith and Vanderbilt (1993), since their phases (non-Abelian Berry phases, $\phi_{x}^{n}$) are proportional to the Wannier charge centers (WCCs) $\bar{x}_{n}$:

[TABLE]

Here, the quantity $\bar{x}_{n}$ is defined $\mod\,a$, as the center of the $n$-th hybrid Wannier function, which is maximally localized in $x$-direction and Bloch-like in $k_{y}$-direction Soluyanov and Vanderbilt (2011a). Typical shapes for the flow of WCCs $\bar{x}_{n}$ in our model are displayed in Fig. 8 ($n$ goes from $1$ to $N_{\text{occ}}=2$).

Symmetry protection of the WCC crossings. Let us assume first $\bm{s}=(0,0)$ and $\bm{b}=(-0.3,0)$, for which the WCCs $\bar{x}_{j}(k_{y})$ anticross, as shown in Fig. 8c. The discussion of Eq. 19 proves that the edge states are protected by $\mathcal{CM}_{y}$ symmetry, which arise from the crystalline symmetry $\mathcal{R}_{2x}$ (multiplied with the effective 2D chiral symmetry $\mathcal{CM}_{z}$, see Table 1). This is the only crystalline symmetry of the 1D system, however the bulk (and therefore the WCCs, too) possesses two additional symmetries: mirror reflection $\mathcal{M}_{x}$ and inversion $\mathcal{I}$. The anticrossing of the WCCs is, in fact, owed to inversion symmetry $\mathcal{I}$, and this feature survives even when $\bm{b}\nparallel\hat{\bm{x}}$, since the magnetic field pseudovector $\bm{b}$ is invariant to inversion. If we have additionally $\bm{s}=(0,0.3)$ [while $\bm{b}=(-0.3,0)$], then only $\mathcal{M}_{x}$ survives in the bulk and thus both the 1D ribbon and the WCCs (Fig. 8a) are gapped.

We remark that $\bar{x}(k_{y})$ in Fig. 8a are fully gapped when $\bm{s}=(0,s_{y})$, but this spin-singlet pairing only removes the anticrossing $\bar{y}(\pi)$. The other crossing, $\bar{y}(0)$, is not protected by any symmetry and can be annihilated by additionally considering an infinitesimal $s$-wave pairing $s_{0}\ll|s_{y}|$ in $\Delta_{\bm{k},\downarrow\uparrow}$.

The nested WLO. If one can separate two sets of Wannier bands which do not touch over the entire range $k_{y}\in[0,2\pi/a]$, one set can be conventionally called the ”occupied” Wannier subspace and can be topologically characterized separately. The nested WLO, introduced in Ref. Benalcazar et al. (2017a) and presented below in Eq. 38, can be computed using the ”occupied” kets $\ket{\bar{x}_{n}(\bm{k})}$ (only one in our model), by first constructing the linear combinations

[TABLE]

The overlap matrices are, similarly to Eq. 33,

[TABLE]

where $\Delta_{y}=(0,\Delta)$ and the new Berry connection $\tilde{\mathcal{A}}^{xy}_{\bm{k}}$ has elements $[\tilde{\mathcal{A}}^{xy}_{\bm{k}}]_{nm}=i\bra{w^{n}_{x,\bm{k}}}\partial_{k_{y}}\ket{w^{m}_{x,\bm{k}}}$. The corresponding WLO follows Eq. 32,

[TABLE]

and may be regarded as the WLO of the edge Hamiltonian $\mathcal{W}_{x,\bm{k}}$. The second-order topological invariant reads

[TABLE]

In our model, $N_{\text{occ}}=2$ and $\widetilde{\mathcal{W}}_{xy,k_{x}}$ is, in fact, a complex number (not a matrix), such that no ”$\det$” is needed and the topological invariant is encoded by a single Wannier band.

The other topological invariant, $\phi_{yx}$, is computed in a similar fashion, by formally swapping the labels $x\leftrightarrow y$. It is worth mentioning that $\phi_{xy}$ and $\phi_{yx}$ are related to the the corner polarization Benalcazar et al. (2017a, b) (up to a factor $2\pi/e$, where $e$ is the electron charge).

Appendix E The statistical phase

We computed the statistical phase as a Berry phase in $t$-space – see Eq. 32 for the Brillouin zone WLO. We evaluated the corresponding WLO by making a discretization of the $t$-space in $n$ steps:

[TABLE]

Here $\alpha\in\{\text{A},\text{B}\}$ refers to either of the two MCSs and the factor $\rho_{N}(n)$ – independent of $\alpha$, but depending on the system size $N$ and on the density of the time mesh $n$ – has been introduced to account for the non-unitarity of the Wilson loop product on the right-hand side (note that $N$ is finite, thermodynamic limit ${N\to\infty}$ is not assumed). The periodic boundary condition $\mathcal{H}^{(t=0)}=\mathcal{H}^{(t=T)}$ implies $\ket{\,\psi_{\alpha}(t_{0})\,}\equiv\ket{\,\psi_{\alpha}(t_{n})\,}$.

In Fig. 9  we plot $\rho_{N}(n)$ for several system sizes. Not surprisingly, it increases when the number of steps $n$ increases and it would reach unity if the braiding speed went to zero (perfect adiabaticity, $\lim_{n\to\infty}\rho_{N}(n)=1$). On the other hand, for a certain $n$, it decreases when the system size increases; physically, this means the speed of the braiding process should be lower for larger systems because the phase space where the MCSs live is larger.