key: cord-0982497-w3sgpody
authors: Alsharari, Abdulrhman M.; Ulloa, Sergio E.
title: Inducing chiral superconductivity on honeycomb lattice systems
date: 2021-07-07
journal: Journal of physics. Condensed matter : an Institute of Physics journal
DOI: 10.1088/1361-648x/ac5a03
sha: e8405d6eb024993c6acfbfaa6264757fe1431489
doc_id: 982497
cord_uid: w3sgpody

Superconductivity in graphene-based systems has recently attracted much attention, as either intrinsic behavior or induced by proximity to a superconductor may lead to interesting topological phases and symmetries of the pairing function. A prominent system considers the pairing to have chiral symmetry. The question arises as to the effect of possible spin-orbit coupling on the resulting superconducting quasiparticle spectrum. Utilizing a Bogolyubov-de Gennes (BdG) Hamiltonian, we explore the interplay of different interaction terms in the system, and their role in generating complex Berry curvatures in the quasiparticle spectrum, as well as non-trivial topological behavior. We demonstrate that the topology of the BdG Hamiltonian in these systems may result in the appearance of edge states along the zigzag edges of nanoribbons in the appropriate regime.

The intriguing properties of topological matter, especially superconductivity, have greatly motivated the ongoing search for such behavior in different material classes [1, 2] . Among these searches, interest in two-dimensional (2D) materials has been intense [3] . Graphene-based 2D systems promise unusual properties and applications in a variety of spintronic designs. Achieving this behavior involves non-trivial properties, such as the quantum anomalous Hall effect, topological phases with complex Berry curvature and associated Chern numbers, helical states that may support Majorana fermion excitations, as well as possible chiral currents along the edges of 2D ribbons [4] [5] [6] [7] .

The quantum anomalous Hall (QAH) regime has received a great deal of attention at the level of singleparticle behavior. An intrinsic magnetization in the system resembles to some extent the effect of a magnetic field on the electronic dynamics with peculiar results [8] [9] [10] . In particular, the magnetization breaks time reversal symmetry (TRS) and lifts the spin degeneracy. A quantized charge Hall conductance characterizes this phase, displaying an insulating state in the bulk and edge states at the boundary of a nanoribbon system.

The addition of superconducting correlations on a QAH system has been predicted to give rise to intriguing new properties and possible Majorana edge states which may be used as robust platforms in quantum computation [11] . The induced pairing correlations have naturally a direct impact on the appearance and specifics of these properties [2] . The symmetries of the resulting pairing channel are closely related to underlying lattice symme- * aalsharari@ut.edu.sa tries as well as to the mechanism mediating the interaction. This is turn produces different quasiparticle spectra with markedly different features.

In contrast, the proximitized honeycomb lattice of graphene may result in broken lattice inversion symmetry and strong spin-orbit coupling that induces an effective Ising (or Zeeman-like) spin field that preserves TRS. Such an effect results in the spin projection perpendicular to the plane being a good quantum number but with opposite value in each valley to preserve TRS [12] [13] [14] [15] [16] . Such a situation could be induced in graphene deposited on a transition metal dichalcogenide substrate, for example [14] [15] [16] . Considering the effect of superconducting correlations on such Zeeman-like system is expected to provide non-trivial Chern numbers and complex Berry curvatures that should affect the resulting symmetries of the quasiparticle excitations [10] .

The classification of superconducting pairing symmetries residing on the honeycomb lattice has been discussed previously [2, 4, 7, 17] . Although deciding which symmetry dominates remains somewhat controversial, an intriguing possibility of superconductivity induced onto proximitized graphene is the appearance of chiral pairing symmetries and their associated broken TRS. For instance, analysis of various symmetry channels in a Hubbard model on a honeycomb lattice shows that a singlet 'd+id'-wave is the dominant pairing symmetry [18] . Highly doped monolayer graphene has been proposed to exhibit these features, with unusual different properties [19] [20] [21] [22] [23] .

The study of either a QAH or a Zeeman-like sytem in the presence of induced chiral d-wave superconductivity has received little attention, and we focus on this issue here. Our primary interest is to explore the relationship between the chiral superconducting pairing functions and topological properties of the corresponding Bogoliubovde Gennes (BdG) Hamiltonian that couples electron and hole subsystems through the superconducting pairing. We explore this effect with a focus on inducing topologically nontrivial character in the proximitized system. We study the quasiparticle (QP) band spectra of the superconducting states, and demonstrate that the topology of the BdG Hamiltonian is different from the constituent subsystems. To this end, we study the topological properties using both the Berry curvature of the QP dispersion and their relation to the appearance of edge states in finite nanoribbon systems. For suitable chemical potential and superconducting pairing strength, we find the appearance of robust midgap states at zigzag edges, well protected by large excitation gaps and momentum transfer.

In Sec. II below we describe the low-energy electronic Hamiltonian for the QAH and Zeeman-like regimes of the honeycomb lattice system, and describe the main features and differences of the resulting single-particle band structure. We then discuss the effect of spatial range on chiral superconducting pairing in Sec. II B, by considering coupling for nearest and next-nearest-neighbors in the lattice, and what symmetry this entails. Section III explores the main features of the QP spectra for the different regimes and superconducting symmetries, including the appearance of topologically protected edge states on zigzag-terminated nanoribbons. Section III D is devoted to a discussion of realistic model parameters and the possibility of realizing the discussed topological phases experimentally. Finally, we discuss the importance and the prospects of these findings in the concluding section.

A. Electronic states of proximitized honeycomb lattice

The Bogoliubov-de Gennes (BdG) Hamiltonian for superconducting pairing on the graphene honeycomb lattice is described in momentum space as [24] 

where s y is a Pauli matrix acting on spin space, η is the chemical potential of the system, and H e is the effective electronic Hamiltonian for proximitized graphene (without superconductivity), which in the basis Φ = (A ↑,B ↑,A ↓,B ↓) T can be written as [15, 16] H e (k) =

Here, k is the 2D momentum vector, and

where a 1 , a 2 , and a 3 are the real-space vectors of graphene (defined explicitly in the next section). The A constant in Eq. (2) represents the effect of a uniform magnetic or exchange field that induces the system to be in a QAH regime. This term is in clear competition with the effective Zeeman field Z(k) induced by the spin-orbit coupling constant λ that embodies the underlying lattice structure and preserves TRS. The R constant represents a possible Rashba spin-orbit coupling, allowed as well by the broken inversion symmetry. For convenience, we set the nearest neighbor hopping as the unit of energy, t = 1, throughout the paper. We now illustrate the differences between the QAH and Ising spin-orbit honeycomb systems before the onset of superconductivity. The electronic Hamiltonian in the QAH regime, which we will denote as an A-System, uses A = 0.05 and λ = 0. Similarly, we will look at a Z-System, where λ = 0.05 and A = 0. In both cases we assume the presence of a Rashba field with R = 0.05, which has as main effect to mix different spins, splitting their degeneracy (increasing the Rashba coupling does not produce a transition to a topological phase [25] ). The parameters chosen capture the essential characteristics of each regime; systems with other coupling constant values exhibit the same qualitative behavior presented here. We note that other terms in the Hamiltonian may exist in a proximitized graphene system, such as that produced by a difference between A and B sublattices (the 'staggered' field) or an intrinsic spin-orbit coupling term. However, we have verified that as long as such parameters remain small, as it is typical in these structures, their inclusion in the model has only a small quantitative effect, with no change in the qualitative behavior. Panels (a) and (b) of Fig. 1 show the electronic band structure of both the Z-and A-Systems, with their corresponding spin S z -projection, along the K-M-K' path in the Brillouin zone. Notice that K and K' states have the same spin ordering in the A-System whereas the spin is reversed in the Z-System. This is understood since QAH breaks TRS, whereas this symmetry is preserved by the Ising spin-orbit coupling in the Z-System. As we will see below, the presence of chiral superconducting order parameters makes the BdG Hamiltonian of the Z-System odd under time-reversal. A chiral gap function breaks TRS and hence allows splitting of bands that were once degenerate.

We want to consider a chiral pairing amplitude Γ(k) that breaks time-reversal symmetry in Eq. 1. Moreover, we want to explore the differences on superconducting QP spectra arising from either nearest neighbor (NN) and next nearest neighbor (NNN) coupling, as the spatial content results in different momentum-space pairing symmetries. The pairing intensities are accompanied with different variables/factors for different directions, as shown in Fig. 1 (c) and (d). As a result, the non-zero elements of Γ(k) in Eq. (1) for d x 2 −y 2 +id xy superconducting order can be expressed as

where k = (k x , k y ), and the graphene lattice constant a is set to 1. We have defined the honeycomb vectors connecting nearest neighbors (sublattice A to sublattice B) as

The corresponding lattice vectors connecting sublattice A (or B) to its NNN are given by

Notice that the symmetries of the singlet pairing function dictate that

Figure 1(e) and (f) show the k-dependent map of the superconducting pairing amplitude for NNN and NN couplings, respectively. Notice that the d-wave chiral amplitude has nodal points or lines when the pairing amplitude transitions between negative and positive values in momentum space. In both cases, the Γ(k) function has at least C 3 symmetry with respect to the Γ point [26] . However, the symmetry at the K and K' valleys is very different in both pairings. This would have significant impact on the QP spectra and associated Berry curvature, as will be seen below.

We study the role of the chemical potential in the QP spectrum as the superconducting pairing strength changes. We focus on energies close to the Dirac points and keep both η (chemical potential) and γ (superconducting pairing strength) relatively small, η, γ < t/4. This focuses on doping in the vicinity of the K valleys. Notice that heavy doping requires a more complicated structure, with consideration of different regions in momentum space. For large doping, the topological properties are dependent on the Γ and/or M points near van Hove singularities, as well as possible band renormalization effects [6, 7, 23] .

As illustrated in Fig. 1(a) and (b) , both the A-and Z-systems exhibit an insulating gap around the charge neutrality point that is proportional to the relevant structure parameter (A, λ and R). To understand the effect of the chemical potential and superconducting pairing function on the different phases of the BdG Hamiltonian, Eq. (1), we show a map of the excitation gap (minimal electron-hole excitation, regardless of its momentum position) near the chemical potential (zero excitation energy) level. Notice that the gap can close at one of the symmetry points (K or K') of the Brillouin zone or away from them, and either way appear gapless in the map. Figure 2 shows the gap behavior in different regimes as function both of chemical potential η and superconducting strength γ. Although the latter is set by the structure and pairing mechanism involved in a given material system, η could be modified somewhat by direct doping or even applied gate voltages. As the gap diagram is symmetric around zero with respect to both η and γ, we plot only the positive quadrant for each case. For γ = 0, the gapless states reflect the single-particle spectral gap seen in Fig. 1(a) and (b). As either η or γ increase, however, the gap behavior changes substantially for the different regimes. The systems studied show distinct, separate regions with contrasting QP behavior.

In the NNN A-System, an excitation gap persists even for large γ, as long as the chemical potential is close to zero, Fig. 2(a) . Increasing η leads to gapless phases, as the pairing function nodal structure dominates. Larger η and γ, however, allows other states to participate in the pairing correlations, which results in extended phases with gap size that is proportional to η and γ, both exceeding the value of A (= 0.05 here).

To further distinguish/identify the characteristic features of these phases, we plot the QP band structure at selected parameters in these regions in Fig. 2 . We focus on one valley for clarity and show band structure and density of states (DOS) of selected points. Each panel to the right of panel Fig. 2(a) corresponds to different (η, γ) values as indicated by the blue/red dots and associated frames. The rightmost panel in each case shows the total density of states (DOS) for each spectrum. The dispersion is shown along the Γ-K-M path and curves are colored according to their S z spin projection. The excitation spectra in all systems here is symmetric around zero, as dictated by particle-hole symmetry. As the A-System is constructed to break TRS, this is carried over after the onset of superconducting correlations, and we find that the spin of electrons and holes are the same at each K valley, although reversed in opposite valleys about the M point (not shown). Figure 2 (b) shows the gap map for the NNN Z-System; an excitation gap in the superconducting phase appears immediately as either chemical potential or γ increases, with a non-monotonic dependence on η. However, beyond η 0.05 = λ, the superconducting phase has a gap that increases with either η or γ. The non-monotonic gap dependence, especially the appearance of a gap closing at η λ, suggests that the system goes through a phase transition. We will analyze this statement further below, as it relates to changes in Berry curvature and Chern numbers. As before, the right frames in Fig. 2(b) show typical QP dispersion and DOS for two cases, indicated by the red/blue dots. Notice how the gap at (η, γ) = (0.06, 0.2), blue dot, appears right at the K point, while for (0.15, 0.15), red dot, the gap has shifted away, suggesting band inversion in the QP states. Notice that since the Z-System preserves TRS even as the pairing sets in, the spin states of the bands around the zero energy are reversed in the K and K' valleys.

In the case of NN coupling shown in Fig. 2 (c) and (d) the scenario is different. For NN-coupling in the A-System, Fig. 2(c) shows that the normal single-particle gap closes for small values of η and γ, while a gapless regime persists for η A and all values of pairing strength γ. For larger η, a superconducting gap emerges and grows with both η and γ values, similar to the Z-System. This indicates that in these regimes the bands near the Fermi level couple strongly through the pairing potential, opening a sizable superconducting gap. The opening large gap as the chemical potential increases can be seen as arising from the overlap of singleparticle states with different symmetry. These bands couple through the pairing function and open a broad QP excitation gap.

The Z-System in 2(d) exhibits a sizable gap that persists for all non-zero superconducting coupling and small chemical potential, η < λ = 0.05; the gap is associated with the single-particle gap in Fig. 1(b) . In other words, this system exhibits no phase transition for small η as the pairing strength increases. Notice that even though there are gapless lines, they do not fully enclose/surround a region. One could then continuously change these parameters without having to close the excitation gap. As η increases beyond λ, the system closes the single-particle gap and reopens in a superconducting phase with a gap that increases monotonically as both η and γ increase.

The QP band structure behaves differently as well. The bands around zero excitation energy are nearly fully spin-polarized around the K (or K') valley, whereas bands at higher energies mix. This is a direct consequence of the NN pairing that couples A and B sublattices. The total density of states acquires sharp peaks around the Fermi level, as well as for energies γ, indicating the emergence of low-dispersive bands, as shown on the right panel of 2(d).

Comparison of the QP dispersions in both regimes further illustrates differences, such as the spin content of the states. While the A-System shows strongly mixed states, the QP bands in the Z-System are fully polarized, except around zero energy, where the states mix due to the Ising and Rashba spin-orbit coupling terms. As we will see be- low, the superconducting gap that appears on opposite sides of gap closing regimes have different characteristics.

Since the pairing function Γ(k) has nodal points in the first BZ, Fig. 1 , shifting η to regions in the spectrum with overlapping particle and hole bands leads to gapless phases in Fig. 2 . It is also understandable from the symmetry of particle and hole bands that large regions are gapless, even for non-zero coupling parameter γ. Lastly, notice that the chiral NNN d-wave coupling appears to have small effect on the bands of the A-System at small η regimes, i.e., the QP band structure does not change much as a function of γ.

Notice the expected enhancement of the DOS near the superconducting QP gap, generally larger in the Z-System than in the A-System. [The latter is not evident in Fig. 2 , as the DOS are enlarged for clarity in the A-Systems to emphasize qualitative behavior.] The larger DOS is accompanied by low band velocity (flatter QP dispersion).

We analyze further the topological properties of these systems using the Berry curvature of the QP dispersions. Differences in Berry curvature would help us distinguish their physical features, as it has implications on the linear response and overall spintronic and valleytronic behavior [27] . One can analyze the Berry curvature of the hole bands per valley, Ω n (k), using [27] Ω n (k) = − n =n 2 Im Ψ n k |v x |Ψ nk Ψ nk |v y |Ψ n k ( n − n ) 2 ,

where v x (v y ) is the velocity operator along the x(y) direction, and n, n are band numbers [15, 27, 28] . The resulting Berry curvature for each of the hole QP bands for the systems in Fig. 2 , as well as the total hole curvature near the valley at K = (−4π/3a, 0) are displayed in Fig.  3 . The total Berry curvature Ω tot shows clear changes in structure, associated with each gap closing event and change in the topological character of the bands seen in Fig. 2 . As the BdG Hamiltonian of these systems has particle-hole symmetry by construction, the curvature plots for the other valley are simply related as Ω tot (Kvalley) = Ω tot (K'-valley). The total hole Berry curvature in the NNN phases is circularly symmetric, reflecting the symmetry of the corresponding superconducting pairing functions near the K valleys. In contrast, the three-fold symmetry of the NN pairing function is reflected in the total Berry curvature for those structures, both in the A-or Z-System, as seen in Fig. 3(e) and (f). As stated earlier, at this level of doping the physics of the honeycomb lattice system is controlled by the two nonequivalent valleys, with identical Berry curvature content. One can calculate the Chern number per valley and double it, or equivalently define the total Chern number by integrating over the entire Brillouin zone,

We notice that in systems of low-doping graphene, the BdG Chern number changes by even numbers across topological phase transitions. This typically indicates the appearance (or suppression) of pairs of propagating edge states whenever a gap closing transition occurs. [29] In the systems we study, despite the gap-closing transitions in the NNN A-System, we find that C BdG = 4 remains unchanged for all gapped phases. This is the result of the single-particle QAH regime having a Chern number of 2, which doubles with the inclusion of the hole sector. After the gap closing topological transition in Fig. 2 one still gets C BdG = 4 but with a different nature, as we explain below in the corresponding edge state analysis.

The NNN Z-System shows C BdG = 0 for small γ values, which then changes to C BdG = 2 at higher η over a small window, and finally goes to C BdG = 4 after gap opening at higher η. Similarly, both the NN A-and Z-Systems show a BdG Chern number of 4 after the gap closing transition.

One finds in the Z-System that the gapped region at low η yields zero Chern number-in contrast to the A-System-as expected from the single-particle picture that respects TRS; the chiral superconducting parameter seems to produce an effective weak electron-hole coupling.

The Chern number in single-particle systems represents the number of chiral fermionic edge states [30] . Similarly, the BdG Chern number indicates the number of chiral edge states in a finite size superconducting system [10] . To better identify the topology of these systems, we have studied the corresponding edge states in a finite ribbon with edges along the zigzag direction.

It is known that a QAH system with onsite superconducting singlet pairing yields a topological phase with edge states in a finite size system that cross the excitation gap [6] . An interesting open question is whether the topological phase can be destroyed by the introduction of chiral pairing interactions. As we will see, these systems exhibit edge states in different regimes.

In the NNN pairing for the A-System, the region of the phase diagram in Fig. 2(a) near η 0 has small gaps and show no well-developed midgap edge states. This is shown in Fig. 4(a) , where the QP dispersions for a zigzag nanoribbon with (η, γ) = (0, 0.2) shows a small gap and edge states (as indicated by color, blue = rightmost, red = leftmost) that promptly hybridize with the bulk states. Well-defined edge states appear after the first gap-closing transition in the system, illustrated in Fig. 4(b) for (η, γ) = (0.2, 0.2). The edge states connect particle and hole bands as two right-and left-mover modes localized around each of the K and K' points. The existence of edge modes in the NNN A-System does not depend on the specific coupling parameter set (i.e., γ, η, R, and A), once the system is in the upper right sector of Fig. 2(a) , which highlights their topological origin. The modes are robust and propagate in opposite direction on alternate edges of the nanoribbon, overlapping spatially without coupling. We emphasize that this results in four (two per valley) co-propagating modes on each side of the ribbon.

The states change drastically for the NNN Z-System. Although edge-less for most of the phases shown here, this regime shows a complex structure of edge states for large γ values, as seen in Fig. 4(c) and (d). The small gap in (c) is due mostly to the presence of Rashba coupling. It is then expected that this structure would not be as robust against perturbation, and that its edge states (either 2 or 4) would be easily perturbed/gapped in the presence of impurities or other perturbations.

A more complicated edge state configuration is found in cases of NN pairing. In both A-and Z-Systems, edge states are present, dispersing across the gap and connecting the two non-equivalent valleys, before eventually merging with the bulk bands. The edge states are well protected by the bulk excitation gap, which makes them robust against perturbations. The phase appears to be topological for large enough values of η and γ that drive the system through the gap-closing phase transition.

There are four edge states crossing the gap, with wavefunctions localized at the left (right) boundaries of the nanoribbon, as shown by the red (blue) colors in Fig.  4 (e) and (f), i.e., there are four modes propagating on the same edge and the same direction. These states are protected by the large momentum and space separation, as well as sizable gap.

Experimental realization of chiral topological superconductors is actively pursued. New opportunities to study and observe topological superconductivity have been identified in different material systems, especially in 2D crystals. Localized Majorana states exist at the boundaries of topological superconductors, and can be created by placing chains or islands of magnetic atoms on an s-wave superconductor [31, 32] . Realization of superconductivity by proximity using hybrid systems, combining magnetism and strong spin orbit coupling, creates novel quasi-particle excitations that open new possibilities for the realization of topological superconductor phases in experiments [33] [34] [35] [36] . The ability to control and tune different phases can be experimentally achieved, allowing their characterization through measurements such as magneto-resistance and related transport probes.

Superconducting states at low electron doping are observed below 1 K in WTe 2 monolayers [37] . Similarly, superconductivity in high-quality 2M-WS 2 single crystals has been reported by several groups experimentally [38, 39] . This unsual material is predicted to exhibit topological nodeless superconductivity, with probable Majorana states bound to vortex cores. Similarly promising chiral-triplet superconductivity is predicted in heavyfermion dichalcogenide UTe 2 [40] . Possible topological superconductivity is predicted theoretically in other 2D transition metal dichalcogenides, as well as other correlated 2D materials [41] .

Hybrid systems such as WSe 2 /bilayer-graphene/ WSe 2 [42] have been shown to exhibit helical edge states arising from the strong spin-orbit orbit under proper gate and displacement field conditions. As bilayer graphene can be made superconducting, this arrangement may result in unusual phases with very interesting properties, similar to those we discuss here. Our estimates suggest reaching doping levels η 0.1t 0.3 eV, which might be possible with electrolyte-gating techniques.

These works and others continue to open novel possibilities for topological superconductivity in van der Waals materials with spin-orbit coupling, representing promising platforms to realize and probe non-trivial phases as the ones we have discussed.

Starting from graphene-like honeycomb lattice models, we study systems with singlet chiral superconductivity to explore the interplay between different topological phases and superconducting range. The NNN correlations result in pairing channels that are rotationally symmetric in momentum space near each valley, whereas for NN correlations the symmetry is three-fold. The chiral superconducting pairing function studied here is the result of the complex orbital combination of d xy and d x 2 −y 2 . This combination breaks time reversal symmetry and is shown to be the preferred channel in some graphene models, including those with strong doping (or intercalation).

Tuning the chemical potential produces possible phase transitions that are accompanied by gapless regimes. The topological aspects of these phases are identified and studied by looking into the formation of edge states in finite size systems, as well as their connection to the Berry curvature of the quasi-particle spectrum and associated Chern numbers in the 2D bulk. Across topological transition points, the system changes between trivial and non-trivial phases, with different features of edge states in different parameter ranges. We trust the present work would lead to further interest and work in these systems. The general nature of the model suggests that it can be applicable, with suitable variations, to a wide class of honeycomb lattice systems that involve different proximity effects.

For discussion and references on this topic see M. M. Scherer

The topological characteristics of models that break TRS while preserving particle-hole symmetry, such as the ones we have here

SEU appreciates the hospitality at CNG/DTU and NBI/KU, as well as the support of Nordea Fonden and the Otto Moensted Visiting Professorship Program.