A Integrals required for exact Green’s function

# Continuum description for sub-gap states at ferromagnetic and spiral ordered magnetic chains in two-dimensional superconductors

## Abstract

We consider sub-gap bands induced in a two-dimensional superconductor by a densely packed chain of magnetic moments with ferromagnetic or spiral alignments. We show that crucially all wavelengths must be taken into account in the calculation of the sub-gap properties, and that in particular gap closings can occur at high momentum due to a mix of Shiba type and magnetic scattering states. The sub-gap bands are always connected to the bulk bands such that it is impossible to single out a fully one-dimensional subband as required for distinct topological properties. To obtain the latter a finite spacing between the impurities needs to be reintroduced such that by zone folding the bands can disconnect from the continuum.

It has been realized in recent years that combinations of ordinary electronic materials can be tuned to lead to quantum states with novel properties. This provides the exciting possibility of actively designing quantum systems with distinct topological states. Topology adds a twist to systems so that surface states develop at the interface between gapped systems of different topological classes. The statistics of such topological surface states can be distinct from both fermions and bosons, which makes them fundamentally attractive, but also suitable for applications such as topological quantum computing (1). Both aspects are well connected through the recent focus on Majorana surface modes (2); (3); (4); (5).

One approach to tune the topology in a superconductor is by the addition of magnetic impurities. As shown by Yu, Shiba, and Rusinov (YSR) in the 1960s (6); (7); (8); (9) a bound state develops in the superconducting gap at the location of a magnetic impurity. This work has found a renaissance when it was realized that arranging such impurities in chains, such that the bound electrons can hop between different impurity sites, can lead to sub-gap bands of different topological classes upon tuning the impurity strength and the alignment of the magnetic moments (10); (11); (12); (13); (14); (15); (16); (17); (18); (19); (20); (21). These Shiba chains provide a complementary approach for one-dimensional (1D) topological physics to the current vast activities on nanowires with spin-orbit interaction, initiated by (22); (23); (24); (25); (26); (27); (28); (29), on topological insulator-superconductor structures (30), or on Josephson junction arrays (31); (32); (33). Similar physics appears also on the edge of magnetic islands on superconductors (34), or in proximitized nanowires (35); (36); (37); (38); (39); (40) and magnetic chains on substrates (41); (42); (43); (44); (45); (46) with spiral magnetic order.

The tight binding model of Shiba chains provides a transparent approach to the sub-gap bands and, through its resemblance with the Ising/Majorana chain (47), to their topological properties (48). But the solutions usually focus only on the novel properties arising from the hybridization of the YSR states.

In this paper we demonstrate that this focus must generally be widened to fully take into account the magnetic properties that already exist in the normal state and the dimensionality of the substrate. The magnetic scattering strongly modifies the band structure within the superconducting gap and in particular causes further gap closings at momenta larger than the Fermi momentum . Furthermore YSR states in two-dimensional (2D) films have a much larger range (49) than in three dimensions such that a simple hybridization picture of YSR states may become problematic. We therefore approach the 2D superconductor directly from the dense chain limit (see Fig. 1), in which the distance between the magnetic scatterers is small, , and the chain can be modeled as a continuous line of scatterers.

We discuss ferromagnetic and spiral orders but not any mechanism (such as in (50); (51)) causing a specific order. We show that a focus on only the long wavelength behavior, similar to the tight binding Shiba chain, misses essential features in the sub-gap band structure, and that an exact treatment reveals further structure in the bands and further gap closing conditions. While similar features have been reported before (52); (53) here we put them in a fully analytic and very accessible context, and provide a complete explanation of their origin. We also show that the sub-gap bands always continuously vary between a 1D and a 2D character such that the application of topological classifications relying on well defined dimensionalities becomes questionable. Our approach relies on an exact analytic evaluation of the Green’s function that takes into account the broken translational symmetry due to the chain and arbitrary scattering processes on the impurities.

Model. The 2D superconductor is described by the Hamiltonian

 H0=∑k,σϵkc†k,σck,σ+(Δc−k,↓ck,↑+h.c.), (1)

with the electron operators for momentum and spin , the dispersion , assumed to have circular symmetry, with effective mass and Fermi momentum , and the -wave bulk gap . We set throughout. The scattering on the dense chain of classical moments (placed at ) is described by

 Hm=Vm∫dxM(x)⋅S(x,y=0), (2)

where are the spatial coordinates, is the scattering strength, is the electron spin operator, and describes the planar magnetic spiral with period and arbitrary orthogonal unit vectors . We do not impose a self-ordering mechanism that would fix (50); (51); (35); (36); (37); (38); (39); (40) but leave as a tuning parameter. By choosing the spin quantization axis perpendicular to we can absorb the momentum shifts induced by through a gauge transformation (54), , such that , for the chain length, describes a line of ferromagnetic scatterers. This also transforms the kinetic energy, , but leaves the pairing term unchanged.

Green’s functions. Green’s functions conveniently provide a full characterization of the system. Since in the gauge transformed basis we maintain the translational invariance along the direction, we will work in mixed momentum-position space. Due to the spin-flip scattering by we must furthermore consider an extended Nambu-spin basis , with the restriction to avoid double counting of states. The retarded Green’s function of the superconducting substrate then takes the form

 g(ω,k)=∑στσ0ω+τσzϵkσ(1−iη)+στσxΔω2−ϵ2kσ−Δ2+iη, (3)

where is an infinitesimal shift, , and , for the Pauli matrices in Nambu and in spin space (with the unit matrices). To go to space we let and calculate . This can be performed exactly for a quadratic dispersion and yields

 g(ω,kx,y)=∑σ−iπρ4kF√~ω2−~Δ2+iη (4) ×{~ωξστσ0+σ~Δξστσx+(1−iη)[(κ2σ−1)ξσ+χσ]τσz},

where is the 2D density of states at the Fermi energy, , , for , and we have defined , and , with . Since the integration is only over all integrals are convergent and Eq. (4) is the exact partial Fourier transform of Eq. (3). Contrary to the full 2D Fourier transform no cutoff by the Debye frequency is required (12), and Eq. (4) is valid for any , , and (see App. A).

In comparison, the treatment of Shiba chains takes the different approach (12) by proceeding through the full real space dependence of the Green’s (or wave) functions and bringing it into the form of a hopping integral in the limit of . This approach becomes inapplicable for dense chains where . It is nonetheless possible to obtain a natural extension of the latter result that incorporates the same focus on long wavelengths and the physics near the gap center. In the limit of small we can replace the integral by an energy integral with constant density of states rather than a partial density of states at fixed , . For the case this leads to the following often used approximate form (55); (56); (3) (equivalent to expanding Eq. (4) in )

 g≈−iπρ√~ω2−~Δ2+iηcos(y√k2F−k2x)√k2F−k2x(τ0σ0~ω+τxσz~Δ). (5)

We will refer to Eq. (5) as the long wavelength approximation (LWA). The LWA is limited to and becomes singular at , which is a consequence of neglecting the singular behavior of at . The LWA becomes thus inaccurate at larger and, in particular, cannot be extended to .

Scattering on the magnetic line impurity preserves in the gauge transformed basis and can be included in form of a T-matrix,

 Missing or unrecognized delimiter for \bigr (6)

such that the full Green’s function reads . The roots of provide the sub-gap bands (for ) and resonances (for ) in the vicinity of the impurities.

Ferromagnetic interface. The differences between the LWA and the exact result are fully analytically accessible for the ferromagnetic interface (), which allows us to provide full explanations where the LWA results fail and how novel sub-gap structures emerge. There are 6 solutions to that can be reduced to solving a cubic equation (see App. C). All but 2 solutions are irrelevant because they lie either above the gap or are complex. The remaining 2 solutions are related by particle-hole symmetry and take the form

 ~ω =±√−(b+γ+δ0/γ)/3, (7)

where , , , with , , , and .

Fig. 2 shows the comparison between the exact and approximate sub-gap bands. At there is a good agreement and we recover the band structure of the LWA,

 Missing or unrecognized delimiter for \bigl (8)

for and . Equation (8) extends the well known expression for single impurity YSR states (6); (7); (8); (9); (12); (52); (53) into a hybridized band, and we identify this structure as the band formed by the overlapping YSR wave functions in the dense impurity limit.

However, at the exact and LWA band structures completely disagree. Notably, the LWA result is wrong about the gap closing properties at small , and fully misses the double band crossing of the Fermi level at larger and the persistent closure of the gap at for large (as indeed the LWA is limited to ). This novel structure arises entirely from the Zeeman type interaction of the magnetic impurity line and does not require superconductivity. Indeed we can recover it from the poles of the T-matrix in the normal state (or in the limit of vanishingly small ), in which the particle and hole bands induced by the magnetic scatterers become . As is increased these bands form a cross like dispersion around . For the magnetism competes with superconductivity, such that for sufficiently large , at , the magnetism is strong enough to break the superconducting gap and recreate the cross like dispersion around a slightly renormalized momentum .

At the same time, as shown in Fig. 2, superconductivity induces a second gap closing at a smaller momentum . Since such a closing can also be obtained from the LWA (Fig. 2), we identify it as a consequence of the YSR correspondence of Eq. (8). With increasing the crossing moves to and becomes gapped at . The exact sub-gap bands thus mix bands of YSR type states and magnetic bands of normal state origin. In App. B we verify that these results are robust, through a fully self-consistent numerical solution of a 2D tight binding model. We observe similar behavior in the extended Shiba band models for (53).

The double crossing of the Fermi surface at and is a seemingly unique feature for the ferromagnetic interface, and is absent for any spiral interface. We also remark that the YSR type result from the LWA entirely misses the crossing and produces only the features at . In particular the LWA indicates a full reopening of the gap at for (instead of ). The exact solution shows, however, that the gap at closes first for increasing and remains closed even after the reopening of the gap at , with potential implications for changing topological properties. This persistence of sub-gap states should be easily detectable, for instance by STM spectroscopy, with distinct peaks at any of the bands’ extrema and a nonzero density of states through the entire energy range for large .

Comparison to 1D model. The magnetic bands near are no speciality of the dimensionality. Indeed let us consider a purely 1D system, a chain with induced superconductivity and magnetic moments. The diagonalization of the Hamiltonian or, equivalently, the determination of the poles of the analogue of Eq. (6), is straightforward (see App. F) and produces the same two sub-gap features, but with gap closings at different strengths of : for at and for at . The role of substrate dimensionality per se is not surprising, but we see that it has substantial impact on the absolute values of the interaction strengths, notably that the scale of the gap closing, with as compared with , becomes significantly larger in 2D.

Spiral ordered interface. For a spiral ordered magnetic interface with the LWA is further restricted to , and at this upper limit its sub-gap bands always join the continuum (see Fig. 6). Thus, using the exact Green’s function becomes a necessity at large . The exact sub-gap dispersion of Eq. (14) becomes then rather complicated but can be straightforwardly computed using computer algebra software. Again only 2 bands contribute to sub-gap energies. Fig. 3 shows that with increasing essentially only the value of the LWA result remains accurate. The difference to the exact result becomes particularly striking for at which the LWA would predict the absence of any sub-gap states but the exact solution provides a gap closure with Dirac type bands at the time reversal symmetric momentum (in the gauge transformed basis) for (and at for general ). For larger this gap reopens. In addition, when tuning from to the closed gap at first opens and splits into two band minima, associated with the features at and . The minimum at then quickly rises into the continuum while the minimum at evolves smoothly into the Dirac type band at .

Conclusions. We have used a continuum approach to appropriately model the dense impurity case . This allows us, by removing lattice spacing related issues, to obtain much more tractable solutions as compared to approaches relying on a finite spacing (12); (13); (53). As an example we have shown that the identification of interaction strengths for the various gap closings is straightforward and, to highlight their physical origin, that gap closures at large momentum are strongly affected by non-superconducting impurity scattering processes. This improved understanding extends from ferromagnetic to spiral magnetic chains, and the results are fully consistent with a self-consistent numerical solution of a tight binding model (see App. B). Similar features are also observed in the 4-band model of (53). Though we focus on sub-gap states, in this framework we obtain the full Green’s function and thus access to all electronic properties.

We have furthermore highlighted with the LWA that the seemingly reasonable approximation of Eq. (5) leads to incorrect results. In particular, the LWA would suggest superconductors are never gapped due to the nonanalytic behavior of the gap closure condition in the limit . In contrast the exact solution shows a persistent gap closure only for and allows us to trace it through all spirals up to the Dirac type closure at for .

Gap closures at high momentum can have a topological significance and have been discussed before. For instance, Ref. (52) considers ferromagnetic chains on a superconducting surface with Rashba spin-orbit interaction, and various topological phase transitions were identified from the different gap closings. Through the gauge transformation the effect of a spiral interface is equivalent to a ‘unidirectional’ Rashba spin-orbit interaction (54), affecting only the momenta while the magnetic spiral unwinds to a ferromagnetic interface. It is thus likely that the gap closures for the various and have similar topological significance. This would also be suggested by the phase transitions in hybridized spiral Shiba chains or nanowires with magnetic spirals (10); (11); (42); (12); (13); (14); (15); (44); (16); (43); (45); (53); (17); (38); (39); (46); (18); (35); (36); (37); (40). However, the topological classification in these references are based on 1D models. With the strong influence of the 2D substrate we demonstrate here it is not easily achievable to single out a 1D subsystem that would allow us to use the same classification methods. In particular, since we can fully trace all bands, we see that upon varying all sub-gap bands connect to the continuum, thereby continuously changing between a 1D and 2D character. Topological indicators that have been designed for a distinct spatial dimension should thus not be used without further examination.

To single out purely 1D subbands we have to rely on a finite , and zone folding tells us that must be large enough to clip off the connection to the continuum. The lowest subbands then take the role of the Shiba bands and can be analyzed by the normal topological tools. This leads to the condition used for Shiba chains (12); (13); (15); (16); (17). But our analysis shows that even then there is always a significant contribution to the sub-gap states from the zone folding of the higher momenta which will produce a variety of further sub-gap bands. These may obscure the Shiba bands and even close the gap, but our results show that they have distinct shapes that should be detectable experimentally and thus can be used to strengthen the identification of the low-lying Shiba bands and their potential for topological states.

Acknowledgments. We thank T. Cren, R. Queiroz, and P. Simon for helpful discussions, and A. V. Balatsky and T. Ojanen for discussions during the early stage of this work. CJFC acknowledges studentship funding from EPSRC under Grant No. EP/M506631/1. Open data compliance: This is a theory paper and all results are reproducible from the given formulas.

## Appendix A Integrals required for exact Green’s function

As noted in the main text, the novel features in our approach stem from the 1D Fourier transform which can be performed exactly for quadratic dispersions. This requires solving two integrals

 I1 =∫∞−∞dkyeiyky1k4y+α2k2y+α0, (9) I2 =∫∞−∞dkyeiykyk2yk4y+α2k2y+α0, (10)

where the constants are functions of . Both integrals are convergent with a choice of semi-circular contour closing for () in the upper (lower) half complex plane, thus the problem reduces to identifying which of the four poles lay within these regions.

The four poles are located at

 ky=±1kF√1−(kx+σkmkF)2±2√−~Δ2+~ω2+iη (11)

in the complex plane, for independent signs and , which can be written in the form by repeated use of the identity for . In this way we find that the poles and always lay in the upper half complex plane, and the other two always in the lower half plane. With this identification the exact Green’s function, Eq. (4), immediately follows from the residue theorem.

## Appendix B Sub-gap bands from self-consistent numerics

We compare the analytic results with the results from a fully self-consistent numerical computation. For the numerical simulation we use a single band tight binding model on a 2D square lattice described by the Hamiltonian

 H=∑⟨i,j⟩,σ(−t−μ)c†i,σcj,σ+∑i[Δici,↓ci,↑+h.c.% ], (12)

where label the lattice sites, and are the spin projections. The kinetic energy is described by the hopping between nearest neighbor sites, denoted by in the summation, with hopping integral and chemical potential . To keep the model assumptions at the necessary minimum we will consider an on-site singlet pairing only, but extensions to allow for triplet pairing have shown no significant difference (triplet pairing amplitudes are generated in any case from the magnetic scattering but have only little influence on the Hamiltonian). The gap function is given by , where is the BCS ground state and the interaction potential strength.

The interaction with the magnetic impurities is given by the Hamiltonian

 Hm=Vm∑i=(ix,0)Mi⋅Si, (13)

where is the interaction strength, are unit vectors that are either ferromagnetically aligned or twisted into a spiral, and is the electron spin operator, with the vector of Pauli matrices. We use a 2D lattice of dimensions with periodic boundary conditions in both directions. The lattice sites are labeled by , but due to the translational symmetry along the direction for the ferromagnetic case we can partially diagonalize the system by making the Fourier transformation . For spiral with winding wave number we choose the spin axes such that rotates in the spin- plane which allows us to treat the impurity scattering as a ferromagnetic interface by letting as for the continuum model. The periodic boundary conditions along are always applied to the gauge transformed basis and lead just to the standard quantization of the , connecting smoothly to the infinite system for . All dependence in the Figures is with respect to this shifted momentum basis, in which the system becomes again diagonal in .

The translational symmetry is broken along the direction by the line of impurity scatterers, which we place at the center position . The system size is chosen large enough (e.g. ) to suppress any influence of the periodicity along the direction, and due to the partial diagonalization we can choose large such as .

The self-consistent solutions is obtained by a fixed point iteration. We choose arbitrary initial values for , diagonalize the Hamiltonian and compute new from the resulting eigenstates. This procedure is repeated until convergence, and a relative error of is typically obtained after 10–20 iterations. In such a solution is suppressed in the vicinity of the magnetic scatterers. To probe the influence of this nonuniformity on the sub-gap features we also compare with a non-selfconsistent solution, obtained by imposing a uniform (corresponding to the self-consistent result far from the interface) and diagonalizing the Hamiltonian only once.

In Fig. 4 we show the resultant sub-gap energy bands obtained from the continuum model and the numerics in the appropriate limit where (compare to Fig. 3). The critical values are decided by analogy with the continuum model as the critical wave vector is precisely the point at which the touching point occurs for a ferromagnetic chain (), and we confirm that this occurs at as expected from the tight binding dispersion. We note in particular that the closing point evolves smoothly to by tuning from to in both the numerics and the continuum model. Interestingly this result appears to be independent of doping so long as one remains within the bandwidth. The gap closure displays the expected dependence, , where depends on the system parameters (see Fig. 4). These characteristics appear in both the self-consistent and non-selfconsistent results, but the former requires a slight rescaling of the values at which the various gap closings occur.

## Appendix C Additional solutions to exact sub-gap bands for ferromagnetic line impurity

The exact solution of the ferromagnetic line impurity subgap bands has three bands, each with a particle hole symmetric pair. These generalize the result in Eq. (7) as

 ~ω=±√−(b+ζqγ+δ0/ζqγ)/3, (14)

where , , , with , , , and , and where are the third roots of unity for . The additional solutions are always above the gap and are plotted in Fig. 5. We note that the extra solutions can be seen to make up the underlying electronic band structure in Fig. 5 (a) which, at larger interactions strengths, combine with the solution plotted in Fig. 2 to make bands similar to those expected of Zeeman shifted electronic bands. We have plotted only places where the exact solutions are purely real but note that there are also complex resonances above the gap. These follow the rest of the expected electronic band displayed in the limit in Fig. 5 (a).

## Appendix D Implicit sub-gap bands for spiral impurity

In the case of a ferromagnetic line impurity () the condition can be reduced to a squared, cubic equation with the solution shown in Eq. (7). Due to the additional structure added by the spiral wavevector such a simple exact solution is no longer possible for . The equation to solve becomes

 (~Δ2−~ω2)(S+S−+C4m)−C2m√(E−−S−)(E+−S+) ×[(E−+S−)(E++S+)+~Δ2+~ω2]=0, (15)

where , . The latter equation corresponds to a polynomial of order 16 in , without any evident symmetry that can be used to reduce it to a simpler form as for . Lengthy but explicit solutions are nonetheless obtainable with the help of standard computer algebra systems, and similar to the analysis shown in Fig. 5, it is straightforward to identify the relevant sub-gap bands, leading to the solutions presented in Fig. 3.

## Appendix E Long wavelength approximation sub-gap bands for spiral impurity

In the main text the comparison to the LWA was restricted only to the ferromagnetic line impurity . In Fig. 6 we present exact sub-gap bands versus the LWA bands as is tuned. In Fig. 6 (a)–(d) we reproduce Fig. 2 for comparison. Figs. 6 (e)–(h) display the restricted region of validity for the LWA if and that it falsely suggests a gap closure at nonzero wave vector until it re-opens near the exact closure derived from the model, at . This is to be contrasted with the exact description of superconductivity surviving for small , reaching a minimum rather than a closing point due to finite lifting degeneracy and the gap closing at the aforementioned point.

In Figs. 6 (i)–(l) () we note that the LWA predicts there are never sub-gap states. This bears stark contrast to the plotted exact bands which tune through a Dirac like closure at .

For there is naturally a physical difficulty in forming sub-gap states due to the separation of spin up and down Fermi surfaces. Nonetheless, it is interesting to note that the closure at persists in this limit and is plotted in Fig. 6 (p).

## Appendix F Comparison to 1D model

The set of sub-gap band features is not unique to 2D but appears in a purely 1D superconductor as well, corresponding to a superconducting chain with fixed (induced) superconductivity with embedded magnetic scatterers. This model is described by the Hamiltonian with

 H =∑kx{∑σϵkx+σkmc†kx,σckx,σ+[Δc−kx,↓ckx,↑+h.c.] +ˇVm[c†kx,↑ckx,↓+h.c.]}, (16)

in which we have already unwound the spiral of magnetic scatterers to a ferromagnetic alignment by the same spin-dependent gauge transformation used in 2D. The result is equivalent to a 1D wire with spin-orbit interaction in a uniform magnetic field. This model can be readily diagonalized, and it reproduces all sub-gap features found in 2D. In Table 1 we outline the gap closures present in this model alongside the equivalent closures in the full 2D model discussed above. In particular, the closure in 1D occurs at as one would expect from a global magnetic field. In 2D, the same closure occurs but requires a comparatively larger magnetic field as the magnetism only exists locally. In this manner, we can associate this first gap closure as being due to the magnetism effecting the underlying normal state by analogy as discussed in the main text.

### References

1. J. K. Pachos, Introduction to Topolological Quantum Computation (Cambridge University Press, Cambridge, UK, 2012).
2. C. Nayak, S. H. Simon, A. Stern, M. Freedman, and S. Das Sarma, Rev. Mod. Phys. 80, 1083 (2008).
3. J. Alicea, Rep. Prog. Phys. 75, 076501 (2012)
4. C. W. J. Beenakker, Annu. Rev. Condens. Matter Phys. 4, 113 (2013).
5. D. Aasen, M. Hell, R. V. Mishmash, A. Higginbotham, J. Danon, M. Leijnse, T. S. Jespersen, J. A. Folk, C. M. Marcus, K. Flensberg, and J. Alicea, Phys. Rev. X 6, 031016 (2016).
6. L. Yu, Acta Phys. Sin. 21, 75 (1965).
7. H. Shiba, Prog. Theor. Phys. 40, 435 (1968).
8. A. I. Rusinov, JETP Lett. 9, 85 (1969).
9. A. V. Balatsky, I. Vekhter, and Jian-Xin Zhu, Rev. Mod. Phys. 78, 373 (2006).
10. T. P. Choy, J. M. Edge, A. R. Akhmerov, and C. W. J. Beenakker, Phys. Rev. B 84, 195442 (2011).
11. M. Kjaergaard, K. Wölms, and K. Flensberg, Phys. Rev. B 85, 020503(R) (2012).
12. F. Pientka, L. I. Glazman, and F. von Oppen, Phys. Rev. B 88, 155420 (2013).
13. F. Pientka, L. I. Glazman, and F. von Oppen, Phys. Rev. B 89, 180505(R) (2014).
14. K. Pöyhönen, A. Westström, J. Röntynen, and T. Ojanen, Phys. Rev. B 89, 115109 (2014).
15. J. Röntynen and T. Ojanen, Phys. Rev. B 90, 180503 (2014).
16. A. Heimes, P. Kotetes, and G. Schön, Phys. Rev. B 90, 060507(R) (2014).
17. P. M. R. Brydon, S. Das Sarma, Hoi-Yin Hui, and Jay D. Sau, Phys. Rev. B 91, 064505 (2015).
18. K. Pöyhönen, A Westström, and T. Ojanen, Phys. Rev. B 93, 014517 (2016).
19. S. Nadj-Perge, I. K. Drozdov, J. Li, H. Chen, S. Jeon, J. Seo, A. H. MacDonald, B. A. Bernevig, and A. Yazdani. Science 346, 6209 (2014).
20. M. Ruby, F. Pientka, Y. Peng, F. von Oppen, B. W. Heinrich, and K. J. Franke, Phys Rev. Lett. 115, 197204 (2015).
21. R. Pawlak, M. Kisiel, J. Klinovaja, T. Meier, S. Kawai, T. Glatzel, D. Loss, and R. Meyer, npj Quant. Inf. 2, 16035 (2016).
22. M. Sato and S. Fujimoto, Phys. Rev. B 79, 094504 (2009).
23. J. D. Sau, R. M. Lutchyn, S. Tewari, and S. Das Sarma, Phys. Rev. Lett 104, 040502 (2010).
24. Y. Oreg, G. Refael, and F. von Oppen, Phys. Rev. Lett. 105, 177002 (2010).
25. R. M. Lutchyn, J. D. Sau, and S. Das Sarma, Phys. Rev. Lett. 105, 077001 (2010).
26. V. Mourik, K. Zuo, S. M. Frolov, S. R. Plissard, E. P. A. M. Bakkers, and L. P. Kouwenhoven, Science 336, 1003 (2012).
27. M. T. Deng, C. L. Yu, G. Y. Huang. M. Larsson, P. Caroff, and H. Q. Xu, Nano Lett. 12, 6414 (2012).
28. A. Das, Y. Ronen, Y. Most, Y. Oreg, M. Heiblum, and H. Shtrikman, Nano Lett. 8, 887 (2012).
29. L. P. Rokhinson, X. Liu, and J. K. Furdyna, Nat. Phys. 8, 795 (2012).
30. L. Fu and C. L. Kane, Phys. Rev. Lett 100, 096407 (2008)
31. B. van Heck, F. Hassler, A. R. Akhmerov, and C. W. J. Beenakker, Phys. Rev B 84, 180502 (2011).
32. B. van Heck, A. R. Akhmerov, F. Hassler, M. Burrello, and C. W. J. Beenakker, New J. Phys. 14, 035019 (2012).
33. F. Hassler and D. Schuricht, New J. Phys. 14, 125018 (2012).
34. G. C. Ménard, S. Guissart, C. Brun, M. Trif, F. Debontridder, R. T. Leriche, D. Demaille, D. Roditchev, P. Simon, and T. Cren, arXiv:1607.06353.
35. B. Braunecker and P. Simon, Phys. Rev. Lett. 111, 147202 (2013).
36. J. Klinovaja, P. Stano, A. Yazdani, and D. Loss, Phys. Rev. Lett 111 186805 (2013).
37. M. M. Vazifeh and M. Franz, Phys. Rev. Lett. 111 206802 (2013).
38. M. Schecter, M. S. Rudner, and K. Flensberg, Phys. Rev. Lett. 114, 247205 (2015).
39. W. Hu, R. T. Scalettar, and R. R. P. Singh, Phys. Rev. B 92, 115133 (2015).
40. B. Braunecker and P. Simon, Phys. Rev. B 92, 241410(R) (2015).
41. I. Martin and A. F. Morpurgo, Phys. Rev. B 85, 144505 (2012).
42. S. Nadj-Perge, I. K. Drozdov, B. A. Bernevig, and A. Yazdani, Phys. Rev. B 88, 020407 (2013).
43. N. Y. Yao, L. I. Glazman, E. A. Demler, M. D. Lukin, and J. D. Sau, Phys. Rev. Lett. 113, 087202 (2014).
44. Y. Kim, M. Cheng, B. Bauer, R. M. Lutchyn, and S. Das Sarma, Phys. Rev. B 90, 060401(R) (2014).
45. I. Reis, D. J. J. Marchand, and M. Franz. Phys. Rev. B 90, 085124 (2014).
46. M. H. Christensen, M. Schecter, K. Flensberg, B. M. Andersen, and J. Paaske, Phys. Rev. B 94, 144509 (2016).
47. E. Lieb, T. Schultz, and D. Mattis, Ann. Phys. 16, 407 (1961).
48. A. Y. Kitaev, Phys. Usp. 44, 131 (2001).
49. G. C. Ménard, S. Guissart, C. Brun, S. Pons, V. S. Stolyarov, F. Debontridder, M. V. Leclerc, E. Janod, L. Cario, D. Roditchev, P. Simon, and T. Cren, Nat. Phys. 11, 1013 (2015).
50. B. Braunecker, P. Simon, and D. Loss, Phys. Rev. Lett. 102, 116403 (2009).
51. B. Braunecker, P. Simon, and D. Loss, Phys. Rev. B 80, 165119 (2009).
52. A. Heimes, D. Mendler, and P. Kotetes, New J. Phys. 17, 023051 (2015).
53. A. Westström, K. Pöyhönen, and T. Ojanen, Phys. Rev. B 91, 064502 (2015).
54. B. Braunecker, G. I. Japaridze, J. Klinovaja, and D. Loss, Phys. Rev. B 82, 045127 (2010).
55. Q. H. Wang and Z. D. Wang, Phys. Rev. B 69, 092502 (2004).
56. B. Braunecker, P. A. Lee, and Z. Wang, Phys. Rev. Lett. 95, 017004 (2005).
You are adding the first comment!
How to quickly get a good reply:
• Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
• Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
• Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minumum 40 characters