Majorana Zero Modes in Graphene

A clear demonstration of topological superconductivity (TS) and Majorana zero modes remains one of the major pending goal in the field of topological materials. One common strategy to generate TS is through the coupling of an s-wave superconductor to a helical half-metallic system. Numerous proposals for the latter have been put forward in the literature, most of them based on semiconductors or topological insulators with strong spin-orbit coupling. Here we demonstrate an alternative approach for the creation of TS in graphene/superconductor junctions without the need of spin-orbit coupling. Our prediction stems from the helicity of graphene's zero Landau level edge states in the presence of interactions, and on the possibility, experimentally demonstrated, to tune their magnetic properties with in-plane magnetic fields. We show how canted antiferromagnetic ordering in the graphene bulk close to neutrality induces TS along the junction, and gives rise to isolated, topologically protected Majorana bound states at either end. We also discuss possible strategies to detect their presence in graphene Josephson junctions through Fraunhofer pattern anomalies and Andreev spectroscopy. The latter in particular exhibits strong unambiguous signatures of the presence of the Majorana states in the form of universal zero bias anomalies. Remarkable progress has recently been reported in the fabrication of the proposed type of junctions, which offers a promising outlook for Majorana physics in graphene systems.


I. INTRODUCTION
The realisation of topological superconductivity (TS), a novel electronic phase characterised by Majorana excitations, has become a major goal in modern condensed matter research.
We here report on a new approach to obtain TS and Majorana states in graphene/superconductor junctions. Key to our proposal is the interaction-induced magnetic ordering of graphene's zero Landau level (ZLL). Coupling this unique state to a conventional superconductor gives rise to novel edge states whose properties depend on the type of magnetic order. In particular, the canted antiferromagnetic phase is a natural host for Majorana bound states. Our proposal combines effects that were recently demonstrated experimentally (tunable spin ordering of the ZLL [18,19] and ballistic [20,21] graphene/superconductor junctions of high-transparency [21] operating in the Quantum Hall regime [20]), and is thus ready to be tested.
While intrinsic TS is rare, it can be synthesised effectively through the coupling of a conventional s-wave superconductor (SC) and tailored electronic gases with spin-momentum locking. Using this recipe, it has been predicted that Majorana excitations should emerge when one induces superconductivity onto topological insulators [14] or semiconductors with strong spin-orbit coupling [15]. Particularly attractive are implementations of one-dimensional TS using either semiconducting nanowires [16,17] or edge states in two-dimensional Quantum Spin Hall (QSH) insulators, since the main ingredients are already in hand. These ideas have spurred a great deal of experimental activity [1][2][3][4][5][6][7][8][9][10][11][12][13]. Despite this progress, however, an unambiguous demonstration of TS is, arguably, still missing. Important limitations of these systems include disorder, bulk leakage, or imperfect proximity effect (the so-called soft gap problem). Thus, it is worthwhile to explore alternative materials.
One particularly interesting option is graphene [22], which exhibits very large mobilities even in ambient conditions, and where a ballistic proximity effect has been recently demonstrated [21]. Graphene was the first material where a topological insulating phase was proposed [23] in the presence of a finite intrinsic spin-orbit coupling. Kane and Mele showed that graphene then becomes gapped around neutrality, and a single helical edge mode with spin locked to propagation direction develops at each edge. In such QSH regime, gapping the edge states through proximity to a conventional superconductor gives rise to a one-dimensional TS along the interface [14]. Graphene's negligible spin-orbit coupling, however, has proved to be a fundamental roadblock in this programme.
In this work we present a simple mechanism to realise the above situation in graphene without recourse to spin-orbit coupling. We consider a graphene ribbon in the Quantum Hall (QH) regime, in which, unlike in the QSH case, time-reversal symmetry is broken by a strong magnetic flux (Appendix A). In contrast to con- ventional two-dimensional electron gases, graphene develops a zero-Landau level at the Dirac point, which has been shown in pristine samples to become split due to electronic interactions [24][25][26][27][28][29]. Experimental evidence [19] points towards spontaneous antiferromagnetic ordering [30,31], although other broken symmetries have been discussed [18]. In this work we consider all possible magnetic orders. Fig. 1 summarises the different possibilities, ranging from ferromagnetic (F) to antiferromagnetic (AF) ordering [32], including canted AF which may be controlled by an external in plane Zeeman field as argued in Ref. 19. The different orders are parameterised by the angle θ between the spin orientation of the ZLL in the two graphene sublattices, so that θ = 0 for F and θ = π for AF. While θ is considered in our model as an externally tuneable parameter as in Ref. 19, we have checked that a mean-field calculation in a honeycomb Hubbard model under an in-plane Zeeman field (Appendix A) yields the same bulk and edge phenomenology presented in this work [31]. Corrections beyond mean field and the Hubbard model have been explored theoretically in the past, and have predicted the formation of a Luttinger liquid domain on an infinite vacuum edge [33]. The corresponding excitation density resembles the non-interacting edge for F order, rather than the AF case. The interacting problem at a highly transparent superconducting contact remains an open problem. We conjecture that, given their robust topological origin, the Majorana phenomena de-scribed here at a mean field level would survive in the Luttinger regime, at least within a limited range of parameters, and with power-law corrections to the transport results. These issues, however, remain beyond the scope of this work.

II. TOPOLOGICAL SUPERCONDUCTIVITY IN QUANTUM HALL GRAPHENE
An infinite ferromagnetically ordered (θ = 0) ribbon in vacuum has a QH mean-field bandstructure [34] as shown in Fig. 1b. (Details on the modelling are given in Appendix A). The ZLL is spin-split into the two spin sectors in the direction of the ferromagnetic order, denoted by |↑ , |↓ at energies ±∆ ZLL /2 respect to the Dirac point. For energies within this gap, a single pair of gapless counterpropagating spin-polarised edge states develop, shown in red for states at the upper vacuum edge of the sample. Note the peculiar situation created in this energy window: edge states are not chiral like in the conventional QH regime, but may rather propagate in both directions with opposite spins, like in the QSH regime. Also, valley degeneracy is lifted at any given edge, which hosts a single state per propagating direction. Upon contacting one edge to a conventional superconductor of gap ∆ SC (as in the sketch of Fig. 1a), while keeping the chemical potential (dotted line) within the ZLL gap, the edge states along the interface develop an induced gap ∆ * SC (see the corresponding Nambu bandstructure of Fig. 1c -red lines, once more, indicate states localised at the upper edge of the graphene ribbon, including now the dense quasiparticle spectrum of the superconductor). The gap ∆ * SC is an important scale in this problem, since it turns out to be a topologically non-trivial gap. This is confirmed by computing the bandstructure's Z 2 topological invariant, Eq. (B2), relevant for quasi-one dimensional D-class systems [35]. Unlike in a conventional QSH system, where time-reversal symmetry is required, the vacuum edge states do not immediately develop a topological gap when contacted to the supercoductor, but requires a reasonably good contact instead. On the other hand, while the conventional QSH metal becomes destroyed by any time-reversal-breaking perturbation (such as inelastic scattering or magnetic impurities), which in turn spoil the non-trivial superconducting gap, this is not the case of the the present implementation, which has a broken time-reversal symmetry from the start. Moreover, we emphasize once more that no spin-orbit coupling at all is necessary for ∆ * SC to develop. The immediate consequence of a non-trivial gap topology is the appearance of zero energy Majorana bound states (MBSs) at an interface with a trivial insulator, following the bulk-boundary correspondence principle. In this case, however, both ends of the topologically gapped superconducting interface are coupled to a gapless vacuum edge, so that the zero modes become delocalised into the continuum away from the interface (see the zero energy local density of states [LDOS] in blue in Fig. 1a). This situation is similar to the fate of MBSs at the ends of a topological proximised semiconductor nanowires when strongly coupled to a metallic environment. [36] The electronic structure associated to an antiferromagnetic ribbon (θ = π, Fig. 1d) is the opposite. The states along a vacuum edge are now (trivially) gapped [31] like the ZLL itself (Fig. 1e). Surprisingly, when contacting the edge to a conventional superconductor two pairs of gapless helical edge modes emerge with spinmomentum locking around conjugate momenta K and K * (Fig. 1f). These unexpected states, spatially spread along the interface (see the blue LDOS in Fig. 1d), are decoupled electron-hole (e-h) superpositions with orthogonal and well defined spin orientation along the AF axis, |K ( * ) ↑ = a|φ e↑ + b|φ h↓ and |K ( * ) ↓ = a |φ e↓ + b |φ h↑ . A full discussion of these states is presented in Appendix D. The two helical edge modes remain gapless as long as no AF canting is present in graphene (θ = π) and intervalley scattering is zero at the interface. We next consider deviations from these two assumptions.
Canting of the AF order may be induced by means of a large enough in-plane Zeeman field, and is thus to some extent externally tunable. This idea was employed in Ref. 19 to tune a graphene QH bar in vacuum between the AF and F regimes, leading to a insulator-to-helical metal transition in edge transport (evolution from Figs. 1e to 1b). A typical canted AF bandstructure (θ = π/2) is shown in Fig. 1h. The vacuum edge states exhibit a topologically trivial and θ-dependent gap, smaller than the bulk ∆ ZLL . Along a superconductor interface, the canted AF helical states are also gapped, Fig. 1i. Like in the ferromagnetic case, this gap is topologically nontrivial. This situation allows for the emergence of true localised zero-energy MBSs at the ends of the superconductor interface, where the edge gap changes topology, see Fig. 1g. The MBSs are topologically protected, and are not destroyed by any small perturbations, or even by modifying the crystal structure of the superconductor (Appendix E).
To understand the full phase diagram of the graphene/superconductor interface quantitatively, it is useful to employ a simplified description in terms of an effective low-energy Hamiltonian for the edge states (see Appendix B). The model is valid for a chemical potential tuned to the ZLL gap, and has the advantage of allowing us to incorporate the effects of atomic disorder along the junction (encoded in an intervalley coupling w, where 'valley' here refers to the conjugate K and K * momenta) and arbitrary spin canting (encoded in an intravalley splitting b θ = ∆ * SC cos[θ/2]). It correctly describes the three possible phases for low-energy interface modes: gapless, trivially gapped and non-trivially gapped. The corresponding phase diagram is shown in Fig. 1k. Typical edge-mode dispersions within each phase (with K and K * points folded onto the Γ point) are shown in Fig. 1j, and are characterised by their energies µ 1,2 < ∆ * SC at k x = 0 for θ = π [AF], w = 0, and their corresponding velocities v 1,2 > 0 (left panel).
The gapless interface regime (light blue in Fig. 1k) is achieved for |w| < w ins ≡ 1 The gapped regimes are characterised by the Z 2 topological invariant of the system, which reads ν = sign w 2 + µ 1 µ 2 − b 2 θ (see Appendix B). A trivially gapped phase ν = +1 is reached for strong intervalley coupling w at the interface, while for intervalley scattering below a threshold w < b 2 θ − µ 1 µ 2 , the interface is one-dimensional TS with invariant ν = −1. The nontrivially gapped regime is most robust against disorder for F order, for which the threshold w reaches its maximum w 0 = (∆ * SC ) 2 − µ 1 µ 2 . Note that to achieve a nontrivial TS interface, the intervalley coupling w should therefore never exceed the induced gap ∆ * SC (this is always the case for sufficiently transparent junctions).

III. EXPERIMENTAL SIGNATURES OF GRAPHENE MAJORANAS
We finally consider measurable signatures of the MBSs in the system. A powerful probe whose feasibility has been recently demonstrated experimentally [20,21] involves interferometry of critical currents in Josephson junctions, and Fraunhofer pattern anomalies in particular [9,12,37]. Lee et al. predicted [38] that topological superconductivity in a short and wide Josephson junction

FIG. 2. Normalised critical current
Ic/I 0 c as a function of magnetic flux through a Josephson junction. Curves from bottom to top (shifted for better visibility) correspond to the non-interacting case (black), gapless AF phase (purple), a trivially gapped AF phase (light gray) and a ferromagnetic phase with a topologically non-trivial gap (orange). Only the latter shows non-decaying non-zero minima in Ic, a consequence of an Andreev spectrum (inset a) with an odd number (one) of edge-resolved zero energy crossings (bottom-right panel). States coloured in red and blue are located at the top and bottom edges, respectively (inset b), while states in gray are spread across the width of the junction (inset c).
could be directly detected in its Fraunhofer pattern, in the form of non-vanishing minima of the critical current for arbitrary magnetic flux through the junction. Such Fraunhofer anomaly was recently observed in a threedimensional topological insulator, and was interpreted as possible evidence of MBSs [12]. Fig. 2 shows the Fraunhofer pattern in a short and wide (10nm × 3µm) graphene Josephson junction in various regimes. The black (bottom) curve corresponds to the non-interacting case (no magnetic ordering, ∆ ZLL = 0), which exhibits the conventional I c (Φ) = I 0 c | sin(πΦ/Φ 0 )|/(πΦ/Φ 0 ) critical current that decays as the inverse magnetic flux Φ through the junction I c ∼ 1/Φ and vanishes at multiples of the flux quantum Φ 0 = h/2e. A similar behaviour is obtained in the gapless regime θ > θ ins , w < w ins (purple curve, with θ = π [AF] and w = 0). In both cases, the junction is host to a narrow quasi-continuum of Andreev bound states around the Fermi energy for any value of the superconducting phase difference φ. The corresponding spectra are shown in the top row of inset a (for Φ = 15.5Φ 0 ). States in red and blue are localised at the top and bottom edges of the junction (inset b), respectively, with gray denoting states spread across the junction (inset c). The case with trivially gapped interfaces (strong intervalley scattering) also exhibits a generic Fraunhofer pattern with vanishing minima (light gray curve, θ = 0 [F] and w > w 0 ). The minima, however, occur at fluxes that are shifted away from integer Φ/Φ 0 at high Φ, while the maxima do not decay like the conventional I c ∼ 1/Φ pattern, which is connected to non-uniform currents across the junction [37]. The Andreev spectrum of the trivially gapped phase is qualitatively different from the gapless spectra, and generally shows a distinct gap devoid of any edge states (inset a, bottom left). For certain values of parameters, it may exhibit zero-energy crossings inside the gap, but in such cases these crossings are accidental (not topologically protected) and there is always an even number of them at a given edge.
The Josephson junction with a non-trivial gap along the contacts is distinctly different from all previous cases. This phase develops two MBSs at each edge (top and bottom) which hybridise to carry a finite supercurrent that never vanishes as the flux increases. The corresponding Fraunhofer pattern thus exhibits a finite background with a superimposed non-decaying oscillation (orange curve in Fig. 2), as described in Ref. 38. The finite minima are roughly one half of the maxima at large flux, and occur away from integer Φ/Φ 0 , at values very close to the zeroes of the trivially gapped phase. (Note, however, that such distinctive pattern develops only in junctions shorter than the Majorana localization length and pierced by a large number of flux quanta, so that the Majoranas are well developed and opposite edges are decoupled). These Fraunhofer anomalies, although not completely unambiguous, thus constitute a measurable hint of the presence of MBSs at the end of a graphene/superconductor interface. Unfortunately, fabricating a very short junction is challenging, in particular due to charge-transfer effects from the superconductors which were neglected here, and which will dope graphene away from neutrality within a few nanometers of the contacts. It is thus important to explore other less stringent experimental schemes that are at the same time not ambiguous. The key is to probe the Andreev spectrum directly for signatures of Majoranas and non-trivial topology.
The presence of the two hybridised Majoranas per vacuum edge in the non-trivial phase manifests in the Andreev spectrum as a single topologically protected zero energy crossing at each edge as φ is increased by 2π (one red and one blue crossing, see bottom-right inset a in Fig. 2). An odd number of such zero energy crossings has been shown [39] to be an direct manifestation of nontrivial topological order ν = −1, and is the underlying reason for the anomalous Fraunhofer pattern of the junc- tion. A completely non-ambiguous demonstration of the presence of MBSs is also thus possible in principle, by directly counting edge-resolved zero-energy crossings using Andreev spectroscopy [40] in a phase-controlled Josephson junction. This may be achieved by measuring differential conductance dI/dV through a normal point contact attached to one edge of the junction, as sketched in Fig. 3(a). We assume a single spinful channel is open. Each Andreev level of energy in the junction is detected as a dI/dV resonance through the probe at bias V = /e, with a resonance width that measures the state's probability density at the point contact. Figs. 3(b,d) show a simulation (see Appendix C for details) of such a dI/dV as a function of φ (at θ = 0 and θ = 0.1π respectively) for a square 500 nm× 500 nm Josephson junction. Note that, unlike in the Fraunhofer simulation, this is not a short junction, since that is no longer a desirable or realistic requirement in the context of Andreev spectroscopy. The number of edge-resolved zero energy crossings as φ is swept from 0 to 2π is one, as corresponds to Majoranahosting SC contacts (compare this dI/dV (φ) to the blue lines in bottom-right inset Fig. 2a [short junction]). Incidentally, the dI/dV (φ) profiles at zero (non-zero) θ, panel b (d), follow the characteristic Andreev level spectra in topological Josephson junctions through long semiconducting nanowires at perfect (non-perfect) transparency [41][42][43][44][45]. The conductance at the crossings is pinned to a universal value 4e 2 /h (white spot), see Appendix C.
While a φ-controlled junction typically requires a SQUID-like geometry and may be experimentally challenging, the existence of MBSs in the junction may be detected even more simply by varying the out of plane magnetic field B z , and hence the total flux Φ through graphene. An increase of Φ by Φ 0 is equivalent to increasing φ by 2π for large Φ/Φ 0 . This is shown in Fig.  3e (compare to Fig. 3d). Moreover, at fixed B z and φ = 0, the MBSs also show up as a zero-bias dI/dV peak (Fig. 3c) as θ is tuned -by the in-plane Zeeman field B -within a range θ L < θ < θ W ≤ θ ins (dashed lines). For the realistic parameters used in Fig. 3 (B z = 1T, L = W = 500nm), and using B 0 = 10T as the typical in-plane field B for complete Ferro polarization, θ L ≈ 0.17π and θ W ≈ 0.5π are reached for B ≈ 9.6T and B ≈ 7.1T, respectively (see Appendix C). Inside this window (point B in Figs. 3a,c), the four Majoranas are concentrated at the corners of the junction and do not overlap, as they decay within a distance smaller than both the width W and length L of the junction, and hence appear as a sharp zero-bias resonance in the dI/dV with exponentially small splitting. The resonance has a universal magnitude of 2e 2 /h at low temperatures, see Appendix C. This is analogous to the zero-bias anomalies reported in pioneering Majorana experiments on semiconducting nanowires [1]. In contrast, for θ < θ L (point A, bulk approaching ferro-ordering), Majorana pairs overlap along the vacuum edge and produce a split resonance. Likewise for θ > θ W (e.g. point C, SC contact close or inside the gapless regime), Majoranas overlap along the SC contact, and develop a (roughly φ-independent) splitting, making the θ > θ L junction strictly trivial. Note also the strong suppression of the width in the dI/dV resonances in this case, due to the exponentially small wavefunction amplitude at the point contact in this geometry (dotted lines overlaid for θ > 0.6π to improve visibility). Essentially the same dI/dV (θ) phenomenology is obtained in a setup with a single superconducting contact (which hosts two Majoranas instead of four), see Appendix C.

IV. DISCUSSION
Our results show that the spontaneous magnetic ordering of the ZLL in graphene enables the creation of topological superconductivity and Majorana states at an interface with a conventional superconductor, even in the absence of spin-orbit coupling in the system. The key is to tune the Fermi energy in the contact into the ZLL gap, and to achieve a good proximity effect therein. The recently characterised samples of Ref. 20 are good candidates to realise our proposal. Impressive progress in controlling graphene filling into proximity gaps has also been reported [46]. We furthermore showed that nonvanishing and non-decaying supercurrent minima in the Fraunhofer pattern across a depleted graphene Josephson junction constitute a characteristic signal of topological order and the presence of Majorana bound states in the junction [38]. Fraunhofer patterns of extraordinary quality have been recently reported in high-transparency ballistic graphene Josephson junctions [21]. We predict even stronger observable signatures of non-trivial topology in Andreev transport spectroscopy, both in the form of an odd number of 4e 2 /h edge-resolved zero-bias crossings versus junction phase difference φ or out of plane magnetic field B, and extended 2e 2 /h zero bias anomalies versus canting angle θ. These experimental probes and the required device parameters are within reach in top laboratories today. We thus expect that the possibility of tuning graphene/superconducting interfaces into a topological phase hosting Majorana bound states could be tested soon. In this section we present the system models employed in this work. A non-interacting graphene flake may be modelled by a nearest-neighbour tight-binding Hamiltonian in an honeycomb lattice, with lattice constant a 0 where n is = c † is c is , i is the site index, s is spin, and φ ij = − e rj ri d r · A( r) is the Peierls phase due to the magnetic flux B z =ẑ · ( ∇ × A). We consider a perturbation i B · S i arising from an external Zeeman field B, where S i = ss c † is σ ss c is is the spin at site i. Intrinsic electron-electron interactions is furthermore included in the local Hartree-Fock approximation [31], which then take the form of a self-consistent Zeeman-like field B U ( r i ) that is different in the two honeycomb sublattices. The total Zeeman-like perturbation H Z thus reads While B is uniform, favoring ferromagnetic (F) ordering, the self-consistent B U ( r i ) is generally opposite for nearest neighbours, favouring antiferromagnetic ordering of the bulk. The combination of the two leads to a B-tuneable, spin-ordered ZLL that can be tuned from AF to F, as discussed in the main text. Further details on the mean field numerics and results can be found in Ref. 31. The hybrid graphene/superconductor (SC) system is described by the Hamiltonian The Fermi energy in H 0 above is µ N for r i / ∈ SC in graphene and µ S for r i ∈ SC at the superconductor (also a honeycomb lattice here, although this is not essential, see Appendix E). Similarly B, B U and B z are zero in the superconductor. However, gauge invariance demands that for a finite magnetic flux in graphene, ∆ SC ( r) = ∆ SC exp − 2e r d r · A( r) , where ∆ SC is the pairing for zero flux at a given superconductor.
The critical currents for the Fraunhofer patterns have been calculated in a wide and short graphene Josephson junction, described by a discretized H using an up-scaled a 0 for numerical efficiency, following the ideas of Ref. 47. The critical current for each magnetic flux is calculated as I c = 2e/ × max φ (dF/dφ), where φ is the superconducting phase difference and F (φ) is the free energy of the junction. [48] Exact diagonalization of the Hamiltonian is used to evaluate F (φ) at zero temperature.
The method to compute the differential conductance for transport spectroscopy is explained in Appendix C.

Appendix B: Low-energy description of a graphene/superconductor junction
In this Appendix we derive a simplified effective Hamiltonian H eff for the two helical edge modes below the superconducting gap ∆ SC that arise in a generic Quantum Hall graphene ribbon and a superconductor. We also characterise its topology by deriving expressions for the relevant topological invariant as a function of model parameters.
We assume the chemical potential lies within the gap ∆ ZLL induced by interactions in the zero Landau Level (ZLL). The gap is associated to a bulk spin ordering described by a canting angle θ between the two sublattice. The effective model incorporates an arbitrary value for θ in graphene and also intervalley coupling due to atomic disorder along the interface. Formally, H eff is a projection of the microscopic Hamiltonian on the basis {|K ↑ , |K ↓ , |K * ↑ , |K * ↓ } of the four AF helical edge states (see Fig. 1f in the main text, and their analytical description in the preceding section) along the junction. We furthermore consider a linearisation of their dispersion in the AF case around 'valleys' K and K * points. These two valleys are folded onto the Γ point by appropriately expanding the ribbon unit cell. Such folding allows us to include intervalley scattering into H eff in a simple way. H eff then takes the form ( = 1) (B1) Here, the intravalley coupling b θ = ∆ * SC cos(θ/2) implements AF canting, and couples opposite spins within the same valley. The intervalley coupling w is spin-independent, and corresponds to the harmonic of wavenumber ∆K = ( K * − K) ·x of any disorder term W close to the interface, w = φ ↑K * |W (∆K)|φ ↑K . Both b θ and w can be chosen real without loss of generality. v 1,2 > 0 are the velocities of the counter-propagating helical states, and µ 1,2 < ∆ * SC are their energy, relative to the Fermi energy, at the Γ point (see Fig. 1j in the main text). The overall structure apparent in H eff is fully determined by the particle-hole symmetry of the underlying Nambu description. Note that H eff only retains terms linear in momentum k x along the interface.
The simplicity of the linearised low energy model above allows us to compute analytical expressions for the topological invariant and the boundaries that separate different phases of the junction (gapless, trivially gapped, non-trivially gapped), as summarised in the main text. We now briefly sketch their derivation.

Topological invariant of edge Hamiltonian
In symmetry class D (superconductors without time reversal symmetry), the one-dimensional topological invariant is Z 2 . It is conventionally defined, for periodic systems, as [49] Hamiltonian for momentum k x , a 0 is the lattice constant of the ribbon lattice, τ x is the first Pauli matrix in the electron-hole sector, and Pf is the Pfaffian. It can be shown in general that Hτ x is antisymmetric at the high-symmetry points k x a 0 = 0, π. The invariant ν is thus fully determined by the structure of H at these two points. The effective Hamiltonian H eff , Eq. (B1), only gives a faithful representation of the full microscopic Hamiltonian H around one of them, the folded Γ point k x = 0 (Fig. 1j, main text). To extract analytic results for ν for the full H using H eff , we must ensure that at the π-point s π does not change when sweeping the parameter space (since this sector of states is not described by H eff ). This is indeed the case in our system for the chosen basis, for which s π = 1. The changes in topology stem from the reconnections of the low-energy edge states that are concentrated around Γ, and are well described by H eff . Higher excited states not included in H eff never cross zero energy, and therefore cannot affect the sign of the Pfaffian of H(0)τ x . One can thus write We have numerically verified the above result by evaluating the Z 2 invariant exactly from the microscopic Hamiltonian H. We finally sketch the derivation of the insulating thresholds w ins and θ ins . These are extracted by computing the solutions for the wavenumber k x of H eff modes at zero energy. Since H eff (k x ) is linear in k x , said k x solutions at = 0 can be obtained as eigenvalues of a matrix −(∂ kx H eff (0)) −1 H eff (0). These can be worked out analytically, and turn out to be all complex (i.e. the interface becomes gapped, see e.g. Fig. 7f) if w > w ins or θ < θ ins , with the expressions given in the main text, Appendix C: Transport spectroscopy 1. Computing differential conductance The computation of the transport spectroscopy results shown in Fig. 3 follows standard techniques of quantum transport. All interactions are incorporated into the mean field solution for B U , Eq. (A2). In this case, the differential conductance dI/dV through a normal contact with a single spinful channel at a certain bias V is given by the BTK formula dI/dV = (2 − R ee + R he )e 2 /h, where the total electron-electron (R ee ) and electron-hole (R he ) reflection probabilities of free electrons incident on the contact are evaluated at energy = eV . The full R matrix, including both electron and hole sectors, may be obtained in terms of Caroli's formula R = Tr (G r ΓG a Γ), where G r/a are the (dressed) retarded/advanced Green functions in the sample, and Γ = −i(Σ−Σ † ) is (twice) the decay rate matrix into the normal probe. All these operators are defined in the Nambu basis, just like e.g. Eq. (D5). The Green's function matrices where solved by inverting the Dyson equation ( ± i0 + − H − Σ)G r/a ( ) = I using efficient linear algebra routines, where H is the Nambu Hamiltonian of graphene, including the superconductors. The self-energy Σ = V † g r V , defined in terms of the hopping matrix V between graphene and the semiinfinite normal lead, and the retarded Green function g r of the latter. g r is obtained by solving the corresponding (self-consistent) Dyson equation ( +i0 + −h−v † g r v)g r = I, where h and v are the on-site and hopping matrices acting on the lead's constituent unit cells.

Estimates for canting angles θL,W
The canting angle θ L is defined as the θ such that the corresponding decay length of edge states along a vacuum edge equals the length L of the Josephson junction, see Fig. 3a. Likewise, θ W is defined as the θ such that the corresponding decay length of edge states along a superconducting edge equals the width W (hence θ W ≤ θ ins ).
The evaluation of θ W can be made by extracting the decay length 1/Imk x of gap states at zero energy from the effective Hamiltonian of Eq. (B1) without disorder (w = 0), and equating that to W . The result comes out simply as where, recall, = 1 and cos(θ ins /2) = 1 2 Note that in the limit W → ∞, θ W = θ ins , as expected. An analogous calculation can be done for the vacuum edge along the y direction, whose effective (normal) Hamiltonian can be written in analogy to Eq. (B1) as This leads to the estimate Note that as L → ∞, θ L reaches a minimum value that corresponds to the threshold where the vacuum edge becomes gapped (non-zero for µ N = 0). Apart from a small oscillatory splitting that decays exponentially with TS length, the dI/dV in the latter case shows a universal 2e 2 /h zero-bias anomaly (orange).
Mean field results within the Hubbard model (see Appendix A and Ref. 31) yield a dependence of canting angle θ with in-plane magnetic field B of the form sin(θ/2) ≈ 1 − (B /B 0 ) 2 , where B 0 is the in-plane field that achieves complete Ferromagnetic polarization.

Differential conductance with a single superconducting contact
The formation and detection of graphene-based Majoranas only requires a single superconducting contact. In that sense, the geometry discussed in Fig. 3 of the main text, while relevant in the context of the Josephson effect, Fraunhofer patterns and parity crossings, is not minimal. A simpler geometry with a single superconducting contact allows for the detection of two Majoranas (instead of four) in the form of a zero-bias anomaly, analogous to that of Fig. 3c. In Fig. 4 we present the dI/dV in such a geometry. Panel b shows the formation of a zero bias 2e 2 /h anomaly within θ W < θ < θ L * (with a L * = 2L + W that now corresponds to the total length of the vacuum edge). Note that, while this setup (panel a) does not allow for a φ-controlled modulation, its θ dependence is qualitative the same as for a Josephson junction. Interestingly, moreover, the analogous to the φ modulation induced by changing the flux B z does operate in this setup just like in a Josephson junction. The reason is that the flux changes the relative phase of the two Majoranas in the SC contact, making them cross at zero energy each time the flux Φ is increased by Φ 0 . This is shown in Fig. 4c.

Relation of the graphene dI/dV to junctions of topological nanowires
The transport spectroscopy results for Majoranas in graphene presented in the main text exhibit three different regimes, labeled as A, B and C in Figs. 3 and 4. These arise due to the interplay, as a function of canting angle θ between delocalization of Majorana bound states along either a vacuum edge and the superconducting contact. From the point of view of the normal point contact, the former case (A) is analogous to a onedimensional TS/normal/TS Josephson junction, where the normal probe is an extra lead coupled to the normal section for spectrocopy. In contrast, in the case for which the Majoranas remain bound to the corners of the sample and do not delocalize (B), the probe is tunnel coupled to the closest Majorana bound state (top-left corner of the sample), and plays the role of a tunnel normal/TS junction like in the zero-bias anomaly experiment of Ref. [1]. Finally, case C is like case B albeit for a trivial normal/S junction without Majorana bound states. These mappings to well-understood systems are useful to understand the universal values of the dI/dV obtained in each case.
For case analogous to B, a one-dimensional normal/TS junction, it is well known [50][51][52][53][54][55][56] that the dI/dV for a long enough TS yields a universal 2e 2 /h zero-bias conductance resonance, a telltale signal of the presence of a Majorana bound state at the contact. If the TS has a finite length, the overlap of the Majorana with its sibling at the opposite end of the TS section gives rise to a splitting of the zero-bias resonance by an energy exponentially small in the TS length divided by the spin-orbit length (or, more precisely, the Majorana localization length [57]). A simple model for such a N/TS junctions was devised by Oreg et al. and by Lutchyn et al. in Refs. 16 and 17. The model is based on semiconducting wires that exhibit a TS phase when proximized to an s-wave SC while under a longitudinal Zeeman V Z exceeding a critical value V c Z = µ N + ∆ 2 . The typical dI/dV in such a junction as a function of V Z and bias V is shown in Fig. 5c. Note that, indeed, for V Z > V c Z , a zero bias anomaly of magnitude 2e 2 /h develops, with a small splitting (oscillating in V Z ) [58][59][60][61].
In the case analogous to A, we have a TS/normal/TS, where the normal portion represents the graphene vacuum edge over which the Majoranas delocalize, and probed with an additional point contact. The Oreg-Lutchyn model in the topological phase yield a transport spectroscopy map, shown in Fig. 5b, with a single zero energy crossing as the junction phase difference φ increases by 2π, just like in Fig. 3(b,d), by virtue of the non-trivial junction topology. Moreover, the zerobias dI/dV at the crossing is pinned to 4e 2 /h (white), again like in the graphene case Fig. 3(b,d). This universal value can be understood intuitively as the addition in parallel of two 2e 2 /h normal/TS zero bias anomalies, one per Majorana (both are coupled to the probe in this geometry), when the two become decoupled at the appropriate φ.
Finally, note that the trivial S/normal/S Josephson junction, corresponding to case C, does not yield a universal zero-bias anomaly at any phase φ, see Fig. 5a. Instead, the finite energy Andreev levels yield a dI/dV at finite bias that approaches (non-universal) 4e 2 /h. Unlike for the topological case A above, however, perturbations to the system may introduce additional normalreflection component to said Andreev levels that suppress this value.
Appendix D: Helical edge states in an AF graphene/superconducting interface

Interface states without Landau levels
The interface states between a superconductor and an antiferromagnetic honeycomb lattice are not related to the Landau level structure. In the particular case of graphene, the magnetic field is the key ingredient to develop magnetic order (due to the large kinetic energy of electrons), which is developed when the kinetic energy is quenched by the magnetic field. However, in a general honeycomb lattice, provided the interactions are large enough at B = 0, electrons might be able to develop AF order, and thus create interface states between a SC even at B = 0. This behaviour might be relevant for other honeycomb systems with smaller hopping strengths than graphene, as silicene, germanene, stanene or honeycomb oxides.
We show in Fig. 6, that such interface sustains the same kind of states as in the case of the antiferromagnetic quantum Hall state. In the case of zigzag interfaces (6c,6d), each valley supports its own set of interface states, whereas for an armchair interface (6b), the two valley are folded and interface states between different valleys can couple. In the former case, if the interface is abrupt enough, a small gap opens up due to intervalley mixing.
If a canting in the magnetic moments is introduced, the zizgag interphase remains gapless. The same happens when only a orbital magnetic field is introduced. Only when both perturbations are present simultaneously, the system is able to enter into the topological superconducting state. Thus, an off-plane magnetic field is mandatory to observe the Majorana bound states. The previous phenomenology, suggests that in order to develop a topological gap, both the spin rotation symmetry and the spatial gauge symmetries have to be broken.

Helical edge states from wavematching
The interface states between a honeycomb antiferromagnet and a superconductor are not intrinsically related to the Landau level spectrum. Although in graphene, the antiferromagnetic state is only expected to arise when the system enters in the quantum Hall regime, a general antiferromagnetic honeycomb lattice might also sustain interface states when attached to a superconductor without a magnetic flux.
In this section, we will show how that interface states naturally arise by an analytic argument in the absence of magnetic field. In particular, a simple wavematching between a E = 0 energy state shows that the boundary between an antiferromagnet and an swave superconductor is able to sustain such state.
To proceed, we will look for bounded solutions such that H|φ = 0 with In the following we will focus in one of the four decoupled sectors, in particular the |e, ↑ K with |h, ↓ K sector

For the antiferromagnet the Hamiltonian reads
On the other hand, for the superconductor the Hamiltonian reads with E = 0 solutions Imposing continuity at the interface x = 0, the full E = 0 solution reads so that a normalizable E = 0 exist for an interface between a trivial antiferromagnet and a trivial Dirac superconductor.

Helical modes from explicit integration
In previous section, we build the E = 0 by wavematching across a sharp interface. However, it is is possible to give a general solution for the interface state between the antiferromagnet and the superconductor. Without loss of generality, in the following we will assume m > 0 and ∆ SC > 0. The Hamiltonian for an arbitrary antiferromagnet and pairing profile for k y = 0 reads with γ 1 , γ 2 , γ 3 defined by Defining γ 4 = −iγ 2 γ 1 and γ 5 = iγ 2 γ 3 , the zero energy equation reads The spinor wavefunction verifies γ 4 φ 0 = γ 5 φ 0 = φ 0 , which allows to build the E = 0 solution which is normalizable provided that m(+∞) > ∆ SC (+∞) and ∆ SC (−∞) > m(−∞), which is the condition of domain wall between superconductor and antiferromagnet. For the case of step profiles, the solution obtained by wavematching is recovered.

Influence of µSC in the critical Zeeman
In the idealised situation in which the superconductor is described as a single-orbital honeycomb lattice at half filling, a arbitrary small Zeeman field is capable of opening the interface topological gap. However, charge transfer processes are expected to shift the chemical potential of the proximized graphene (SC region in Fig. 1). In this situation, band bending of the interfacial states leads to a one dimensional gapless state (see Fig. 7c). In order to reach the interfacial topological superconducting state, the bended bands (Fig. 7c) have to be moved up in energy. This can achieved by increasing the in-plane field, so that θ < θ ins . The critical in-plane field B x as a function of doping, which separates the gapless and topologically gapped states is shown in Figs. 7e,f. For small µ SC , the critical field increases linearly, leading to small critical fields at small doping in the SC, whereas for large doping the critical field saturates.

Emergence of AF helical states from a topological point of view
An analysis of the emergent AF helical edge modes in terms of topology can be made, but it is less rigorous mathematically than the topological superconductor order in the canted AF phase. As shown in Fig. 1k, there is a finite volume in parameter space for which the SC contacts are helical metals. The relevant symmetry class within the ten-fold way [62] for the infinite system is the 2D class D, and its topological invariant is Z, which corresponds to the number of chiral edge states at a surface [63]. Within this language the helical SC contact has a trivial (zero) Z invariant, since the number of right minus left propagating modes is zero. This is actually the reason why crossing the θ = θ ins destroys the helical edge states without an intervening bulk-gap inversion. Therefore, the reason for the existence of helical states for θ > θ ins cannot be found in the standard homotopy classification. It is rather an instance of non-trivial valley Chern number.
If one computes the Z invariant in 2D of our system, both in the graphene side and on the superconductor side, one needs to integrate the Berry curvature of the Nambu bands. For the superconductor one obtains negligible Berry curvature for all momenta. However on the graphene side (and choosing the magnetic unit cell to compute the bands), one finds that while the integrated curvature is zero (hence Z is zero), it is the sum of two integer and opposite contributions from different valleys. For a specific spin sector (e.g. |e ↑ , |h ↓ ), one valley has partial integral 1 and the other -1. Therefore, assuming valley symmetry is preserved at the graphene/SC interface (that is w = 0 in the phase diagram of Fig. 1k, i.e. the contact is transparent), one can invoke a bulk-boundary correspondence principle for each of the two valleys and spin sectors independently, which yields two pairs of counterpropagating (helical) states (one per valley and spin). If the contact is not transparent (w > 0) but intervalley scattering is below the threshold w < w ins , the helical states will split, but the contact will still be metallic (since the splitting is smaller than the energy where they cross). Hence the AF helical states are an instance of weak topology from the two valleys, just like e.g. the helical modes in works like Refs. [64,65]. Note, however, that the mathematical standing of these arguments are less sound than the conventional ones from full-Brillouin-zone invariants, since to our knowledge there is no rigorous theorem that guarantees the existence of surface states from partially integrated (valley) Chern numbers.
Appendix E: Square-lattice superconductor In Appendix D we have considered the superconductor arising from electrons hopping in a honeycomb lattice and subjected to pairing potential. Nevertheless, the fact that the interface states persist even upon doping of the superconductor, suggests that their existence goes far beyond what our analytic argument might suggest. Actually, we here show that a honeycomb superconducting lattice is not mandatory, so that even an interface between a square superconductor will give rise to localized Majorana states.
To illustrate this, we show in Fig. 8 the bandstructure of an interface between the canted Quantum Hall antiferromagnet, and the local density of states for a finite system. In the case of fully collinear antiferromagnetism, the interface sustains a gapless channel. When the moments are canted, a topological gap in the inter-facial bands opens up and localized Majorana modes show up.