Harnessing ultraconfined graphene plasmons to probe the electrodynamics of superconductors

Significance Superconductivity and plasmonics constitute two extremely vibrant research topics, although with often nonoverlapping research communities. Here, we bridge these two active research fields by showing that graphene plasmons’ unprecedented light localization into nanometric scales can be exploited to probe the electrodynamics (including collective excitations) of superconductors. Our findings are important both from a fundamental standpoint, representing a paradigm shift (i.e., probing of Higgs modes by light fields), and also for future explorations interfacing nanophotonics with strongly correlated matter, which holds prospects for fostering additional concepts in emerging quantum technologies.


I. INTRODUCTION
The superconducting state is characterized by a spontaneously-broken continuous symmetry [1].As a consequence of the Nambu-Goldstone theorem, superconductors are expected to display two kinds of elementary excitations: the so-called Nambu-Goldstone (NG) and Higgs modes [2][3][4].The NG mode is associated with fluctuations of the phase of the order parameter, whereas the Higgs mode is related to amplitude fluctuations of the same.In superconductors and electrically charged plasmas, the NG (phase) mode couples to the electromagnetic field and its spectrum effectively acquires a gap (mass) due to the long-range Coulomb interaction (Anderson-Higgs mechanism) [2]; this gap corresponds to the system's plasma frequency [1,5,6].On the other hand, the Higgs (amplitude) mode is always gapped, and in superconductors its minimum energy is equal to twice the superconducting gap [7].Curiously, one often encounters in the literature statements that the Higgs mode does not couple to electromagnetic fields in linear response, making it difficult to observe in optical experiments [2,8].Experimental detection has only been achieved through higher-order response [9], e.g., by pumping the superconductor with intense terahertz (THz) fields and measuring the resulting oscillations in the superfluid density [10][11][12][13][14].
Naturally, the light-Higgs coupling is subjected to conservation laws, whereby translational invariance manifests in the conservation of wave vectors.Since far-field photons carry little momentum, wave vector conservation cannot be satisfied and the coupling is suppressed.However, little attention has been given to the fact that, strictly speaking, the linearresponse coupling of the electromagnetic field to the Higgs mode only effectively vanishes in the q → 0 limit [8,15].Contrasting this, at finite wave vectors-i.e., in the nonlocal regime-, the linear optical conductivity of the superconductor has a contribution associated with the coupling to the Higgs mode [8,15,16].Hence, electromagnetic near-fields provided by, for instance, plasmons, emitters, or small scatterers, can couple to such amplitude fluctuations and therefore constitute a feasible, promising avenue toward experimental observations of the Higgs mode in superconductors.In this context, ultra-confined graphene plasmons [17,18] constitute a new paradigm for probing quantum nonlocal phenomena in nearby metals [18][19][20][21][22][23], while their potential as tools for studying the intriguing electrodynamics of strongly-correlated matter [24][25][26] remains largely virgin territory.
Here, we exploit the unprecedented field confinement yielded by graphene plasmons (GPs) for investigating the nearfield electromagnetic response of a heterostructure composed of a graphene sheet separated from a superconductor by a thin dielectric slab (see Fig. 1).Both the superconductor and the graphene sheet are characterized by their optical conductivity tensors [16,17].The optical conductivity tensor of the superconductor is intrinsically nonlocal [16], whereas for graphene it is possible to employ a local approximation at wave vectors much smaller than graphene's Fermi wave vector [17,20,21].We show that the coupling between the Higgs mode in the superconductor and plasmons in the graphene manifest itself through the existence of an anticrossing-like feature in the near-field reflection coefficient.Furthermore, the energy and wave vector associated with this feature can be continuously tuned using multiple knobs, e.g., by changing (i) the temperature of the superconductor, (ii) the Fermi level of the graphene sheet, or (iii) the graphene-superconductor separation.Finally, we suggest an alternative observation of the GPs-Higgs coupling through the measurement of the Purcell enhancement [18,27,28] near the heterostructure.To that end, we calculate the electromagnetic local density of states (LDOS) above the graphene-dielectric-superconductor heterostructure; our results show that, in the absence of graphene, the coupling between the superconductor's surface polariton and its Higgs mode leads to an enhancement of the LDOS near the frequency of the latter.The presence of graphene changes qualitatively the behavior of the decay rate around the frequency of the Higgs mode, depending strongly on the emitter-graphene distance.

II. COUPLING OF THE HIGGS MODE OF A SUPERCONDUCTOR WITH GRAPHENE PLASMONS
A. Theoretical background

Electrodynamics of BCS-like superconductors
The electrodynamics of superconductors and other stronglycorrelated matter constitutes a fertile research area [24,25].In the following, we assume that the superconducting material is well-described by the Bardeen-Cooper-Schrieffer (BCS) theory of superconductivity [16,29,30].Chiefly, the microscopically-derived linear optical conductivity tensor of a superconductor requires a nonlocal framework due to the finiteness of the Cooper-pair wave function.For homogeneous superconducting media, the longitudinal and transverse components of the nonlocal optical conductivity tensor-while treating nonlocality to leading-order [31]-can be expressed as [16,32,33] respectively, where σ D (ω) = ine 2 m(ω+iγ) is the Drude-like conductivity, and the dimensionless coefficient ᾱ(ω, T ) amounts In the previous expression, 2m is the singleparticle energy of an electron with wave vector k, In possession of the response functions epitomized by Eqs.(1), we employ the semiclassical infinite barrier (SCIB) formalism [18,34] to describe electromagnetic phenomena at a planar dielectric-superconductor interface [32,33,35].Within this framework, the corresponding reflection coefficient for p-polarized waves is given by (see Supplemental Material [36]) [18,34] with and Ξ has the form L (q, ω) where q 2 = q 2 + q 2 ⊥ , and L,T = ∞ + iσ L,T /(ω 0 ) are the components of the superconductor's nonlocal dielectric tensor (we assume ∞ = 1 hereafter).

Electrodynamics in graphene-dielectric-superconductor heterostructures
With knowledge of the reflection coefficient for the dielectric-superconductor interface (2), the overall reflection coefficient, i.e., associated the dielectric-graphenedielectric-superconductor heterostructure, follows from imposing Maxwell's boundary conditions [37] at all the interfaces that make up the layered system.At the interface defined by the two-dimensional graphene sheet, the presence of graphene enters via a surface current with a corresponding surface conductivity [17].
Signatures of system's collective excitations can then be found by analyzing the poles of corresponding reflection coefficient, which are identifiable as features in the imaginary part of the reflection coefficient, Im r p .

B. Signatures of the Higgs mode probed by graphene plasmons
Like ordinary conductors [38], superconductors can also sustain surface plasmon polaritons (SPPs) [39,40].In turn, these collective excitations can couple to the superconductor's Higgs mode [32,33].Typically such interaction is extremely weak due to the large mismatch between superconductor's plasma frequency, ω p , and that of its Higgs mode, ω H = 2∆ 0 / ; for instance, ω H /ω p ∼ 10 −2 in high-T c superconductors, with ω p and ω H falling, respectively, in the visible and THz spectral ranges.As a result, at frequencies around ω H the SPP resembles light in free-space and thus the SPP-Higgs coupling is essentially as weak as when using far-field optics (Fig. 2a).
On the other hand, graphene plasmons not only span the THz regime but also attain sizable plasmon wave vectors at such frequencies [17,18].Moreover, when the graphene sheet is near a metal-or a superconductor for that matter-graphene's plasmons become screened and acquire a nearly linear (acoustic) dispersion, pushing their spectrum further toward lower frequencies (i.e., a few THz) and larger wave vectors [18][19][20][21][22]. Therefore, these properties of acoustic-like GPs can be harnessed by placing a graphene monolayer near a superconducting surface, thereby allowing the interaction of graphene's plasmons with the Higgs mode of the underlying superconductor (Fig. 2b).In this case the plasmon-Higgs interaction is substantially enhanced, a fact that is reflected in the observation of a clear anticrossing in the GP's dispersion near ω H , which, crucially, is orders of magnitude larger than that observed in the absence of graphene (cf.Fig. 2a-b).
Furthermore, the use of graphene plasmons for probing the superconductor's Higgs mode comes with the added ben-efit of control over the plasmon-Higgs coupling by tuning graphene's Fermi energy electrostatically [17,18,[41][42][43].This is explicitly shown in Fig. 3a, for a vacuum-hBN-graphene-hBN-superconductor heterostructure; as before, the coupling of GPs with the superconductor's Higgs mode manifests itself through the appearance of an avoided crossing in the vicinity of ω H , which occurs at successively larger wave vectors upon decreasing E gr F .Another source of tunability is the graphenesuperconductor distance, t, which corresponds to the thickness of the bottommost hBN slab.Strikingly, current experimental capabilities allow the latter to be controlled with atomic precision [19,20].We exploit this fact in Fig. 3b, where we have considered the same heterostructure, but now we have varied t instead, while keeping E gr F fixed.Naturally, the manifestation of the GP-Higgs mode interaction seems to be more pronounced for smaller t, reducing to faint feature at large t (cf. the result for t = 50 nm).Lastly, it should be noted that the net effect of decreasing the graphene-superconductor separation t is the outcome of two intertwined contributions: the graphene-superconductor interaction is evidently stronger when the materials lie close together, but equally important is the fact that the (group) velocity of plasmons in the graphene sheet gets continuously reduced as t diminishes due to the screening exercised by the nearby superconductor (and, consequently, the GP's dispersion shifts toward higher wave vectors, even reaching the nonlocal regime) [18,19,22].

III. HIGGS MODE VISIBILITY THROUGH THE PURCELL EFFECT
One way to overcome the momentum mismatch and investigate the presence of electromagnetic surface modes is to place a quantum emitter [17,[45][46][47] (herein modeled as a point-like electric dipole) in the proximity of an interface and study its decay rate as a function of the emitter-surface distance.With the advent of atomically-thin materials, and hBN in particular, all the relevant distances, i.e., emitter-superconductor, emittergraphene, and graphene-superconductor, can be tailored with sub-nanometer precision [e.g., by controlling the number of stacked hBN layers (each ∼ 0.7 nm-thick)].Although good emitters in the THz range are relatively rare, semiconductor quantum dots with intersublevel transitions in this range and long relaxation times have been demonstrated [48].The modification of the spontaneous decay rate of an emitter is a repercussion of a change in the electromagnetic LDOS, ρ(r), and it is known as the Purcell effect [18,27,28].Specifically, the Purcell factor-defined as the ratio ρ(r) ρ 0 (r) , where ρ 0 (r) is the LDOS experienced by an emitter in free-space-can be greatly enhanced by positioning the emitter near material interfaces supporting electromagnetic modes (which are responsible for augmenting the LDOS).
For p-polarized waves, the orientation-averaged Purcell factor-or, equivalently, the LDOS enhancement-can be determined via [17,28] where s z = √ 1 − s 2 and z is the vertical position of the emitter relative to the surface of the topmost hBN layer.
Figure 4 shows the LDOS enhancement experienced by an emitter in the proximity of a superconductor; Figs.4a-b,d-e refer to the case in the presence of graphene (located between the superconductor and the emitter), whereas Fig. 4c depicts a scenario where the graphene sheet is absent.The graphene sheet modifies the LDOS, affecting not only the absolute Purcell factor but also the peak/dip feature around the energy of the Higgs mode, ω h = 2∆.Such modification depends strongly on the emitter-graphene separation d (Figs.4a-b).Fig. 4d shows the LDOS enhancement for T > T c (i.e., above the superconductor's transition temperature) and thus the feature associated with the Higgs mode vanishes; all that remains is a relatively broad feature related to the excitation of graphene plasmons, which shifts to higher energies as the emitter-graphene distance is reduced.
Lastly, Fig. 5 depicts the LDOS enhancement for different values of graphene's Fermi energy (which can be tuned electrostatically), for two fixed emitter-graphene distances: d = 13 nm (top row of panels) and d = 2 nm (middle row of panels).For weakly doped graphene and the larger d the sharp feature associated with the hybrid GPs-Higgs mode dominates the Purcell factor, being eventually overtaken by the broader background with increasing E gr F .To unveil the nature of the peak seen on the Purcell factor, we plot in the bottom row of Fig. 5 the differential LDOS, which is the kernel of integration of Eq. ( 3).A first glance at that panel reveals that the line-shape of the (q -space) differential LDOS possesses a Fano-like line shape.This is most evident for the case where E gr F = 50 meV.To confirm this, we fit the numerical curve with an analytical expression using a Fano line shape resonance convoluted with the density of states in momentum space, and an exponential tail determining the behavior of the differential LDOS at large q .The simplest form of a Fano line shape reads [49] s Fano (q ) = s 0 q F + q gp −q ∆q where s 0 controls the amplitude of the resonance, q F is called the asymmetry parameter and determines the asymmetry of the line shape relative to the maximum of the resonance, whose position is given by q gp , the GP's wave vector, ∆q determines the width of the resonance, and Γ is a parameter characterizing the losses.For the q -space differential LDOS, we suggest the following convolution formula, dρ(q ) dq = s Fano (q )q 4 e −q /q 0 , ( where the term q 4 plays the role of a density of states in momentum space and e −q /q 0 gives the decay of the differential LDOS for large q .In the bottom panel of Fig. 5, we fit Eq. ( 5) to the numerical differential LDOS and find an excellent agreement (the fitting parameters are given in the Supplemental Material [36]).From the fit, we find that Γ and q 0 have the same value for both curves, each of them corresponding to two different Fermi energies.This is to be expected, as losses and q 0 are constant across the two devices which differ only on the Fermi energy.The values of the q gp between the two curves due to the two different Fermi energies are naturally different (smaller Fermi energies imply large wave vectors for the same frequency, which we have chosen to be the Higgs frequency).The Fano factor is q F = 1.0 for the curve corresponding to E gr F = 50 meV, while being 2.5 times larger for the curve representing the setup with E gr F = 250 meV (bottom panel of Fig. 5).It is clear from s Fano (q ) that the Fano line shape is symmetric when both q F = 0 and q F 1.This explains why the curve for E gr F = 50 meV is more asymmetric than that for E In a and b the graphene Fermi energy has been set at 0.25 eV and T = 1 K for the solid curves and T = 94 K (above T c ) for the dashed curves; the graphene sheet is placed 4 nm above the superconductor surface.We show results for two emitter-graphene distances: 13 nm (a) and 36 nm (b).In c the Fermi energy of graphene has been set to zero and T = 1 K.The red curve corresponds to an emitter-superconductor separation of 17 nm and the blue curve of 40 nm.In d we show results for the same distances as in a (red curve) and b (blue curve), but for T = 94 K.In e we show how the Purcell factor depends on the graphene-superconductor distance t at the frequency of the Higgs mode, ω h = 2∆ ≈ 28.32 meV.The other parameters are kept fixed: E gr F = 0.5 eV, T = 1 K, and emitter-graphene distance of 13 nm.(The nonlocal optical conductivity of graphene was used.)shape appears when a continuum is coupled to a single mode.
Here the single mode is the hybrid GP-Higgs mode with wave vector q gp , which, we stress, is indeed a hybridized mode between the graphene plasmon and the Higgs mode.For the case of E gr F = 50 meV the Fano line shape is well defined meaning that the coupling is strong thus leading to an enhanced LDOS, as seen in the top first panel of the topmost row in Fig. 5.That is, an additional decay channel is open.The coupling between the hybrid GP and the continuum of radiation is stronger for smaller distances.This is well perceivable in the middle row of Fig. 5, where a peak is always present in the LDOS.Let us conclude stressing that Eq. ( 4) holds only when the coupling between the continuum and the hybrid GP mode is sufficiently strong.For weak coupling, there will be deviations from the Fano line shape description.

IV. CONCLUSION
We have shown that signatures of a superconductor's Higgs mode can be detected by exploiting ultra-confined graphene plasmons supported by a graphene sheet placed in the superconductor's proximity.In particular, the presence of the Higgs mode for T < T c can be readily identified through an anticrossing feature that attests the coupling between graphene Here we show the effect of changing graphene's Fermi energy (indicated at the top of each column) while keeping all other parameters fixed: T = 1 K, emitter-graphene distance (d = 13 nm for the top line and d = 2 nm for the middle line), graphene-superconductor distance t = 4 nm.For d = 13 nm, the dependence of the decay rate on the emitter's frequency changes qualitatively from low (E gr F = 50 meV) to high (E gr F = 250 meV) graphene doping.In the bottom panel we depict the differential LDOS given by the integration kernel of Eq. ( 3), solid lines, and by the analytical expression given by Eq. ( 5), dashed lines.The energy is fixed at the value ω H . (The nonlocal optical conductivity of graphene was used.)plasmons and the superconductor's Higgs mode.Further, we suggest that the excitation of the Higgs mode of superconductors could also be detected through the emergence of a peak or dip in the near-field's Purcell factor, and whose shape (peak or dip) depends the coupling between the emitter continuum of radiation of the hybrid GP-Higgs mode.This coupling is most efficient for small Fermi energies and short distances between the superconductor and the emitter.
Experimentally, the interactions of GPs with the superconductor can be easily obtained patterning graphene into nanoribbons.In this geometry, the graphene plasmons can be excited by far-field methods, thereby generating the necessary optical near-field that strongly interacts with the superconductor's Higgs mode.In addition, the GP-Higgs interaction can also be experimentally investigated using state-of-the-art cryogenic scanning near-field optical microscopy (SNOM) [50].
In addition to quantum dots working as THz emitters, two other possibilities of emulating quantum emitters that are experimentally feasible are the following: (i) an antenna on the surface of a superconductor-graphene heterostructure illuminated with a THz beam; (ii) a variant of an ultra-fast Auston switch in a current carrying superconductor: The superconductor when illuminated with a fs laser pulse heats an it goes above the critical temperature.As superconductivity recovers, THz light is emitted.These are not quantum emitters per se, but these are emitters that one can relatively easily integrate with various types of superconductors.
Finally, there are a number of open questions that can spur from this work, e.g., if conductive thin films were added in direct electrical contact with the superconductor, then bound Andreev quasiparticle states inside the superconducting energy gap can form, being solutions to the Bogolubov-de Gennes equations [51].Another enticing outlook is the prospect of using highly-confined GPs for investigating Josephson plasma waves in layered high-T c superconductors [24,52,53].The present formalism could be extended to the coupling of the above-noted types of modes.
The work presented here sheds light on the fundamentals of collective excitations in novel architectures containing twodimensional materials and superconductors, paving the way for prospective experimental investigations on the electrodynamics of superconductors using ultra-confined graphene plasmons.

Figure 1 .
Figure 1.Schematic of the device studied in this paper.Schematic view of the heterostructure composed of a superconducting substrate, a few atomic layers of hexagonal boron nitride (hBN), a single sheet of graphene and a capping layer of hBN.The electric dipole is represented by the red-blue sphere.

Figure 2 .
Figure 2. Spectra of our devices as retrieved from the loss function.Spectra of surface electromagnetic waves in superconductors (a) and graphene-superconductor (b) structures, obtained from the calculation of the corresponding Im r p . a Dispersion diagram of SPPs supported by a vacuum-superconductor interface (the hatched area indicates the light-cone in vacuum).The inset shows a close-up of an extremely small region (notice the change of scale) where the SPP dispersion crosses the energy associated with the superconductor's Higgs mode; here, ∆E = E − ω H and ∆q = q − ω H /c. b Dispersion relation of GPs exhibiting an anticrossing feature that signals their interaction with the Higgs mode of the nearby superconductor; the graphene-superconductor separation is t = 5 nm.Setup parameters: We take T = 1 K; moreover, n = 6 × 10 21 cm −3 (so that E F ≈ 1.20 eV and ω p ≈ 2.88 eV), γ = 1 µeV, and T c = 93 K for the superconductor, and E gr F = 0.3 eV and γ gr = 1 meV, for graphene's Drude-like optical conductivity.

Figure 3 .
Figure 3. Tuning the hybridization of acoustic-like plasmons in graphene with the Higgs mode of a superconductor in air-hBNgraphene-hBN-superconductor heterostructures.The colormap indicates the loss function via Im r p .Spectral dependence upon varying the Fermi energy of graphene (a) and the graphene-superconductor distance (b).Setup parameters: the parameters of the superconductor are the same as in Fig.2, and the same goes for graphene's Drude damping.The thickness of the bottom hBN slab is given by t, whereas the thickness of the top hBN slab, t , has been kept constant (t = 10 nm).Here, we have modeled hBN's optical properties using a dielectric tensor of the form ↔ hBN = diag[ xx , yy , zz ] with xx = yy = 6.7 and zz = 3.6[19,43,44].

grFFigure 4 .
Figure 4. Purcell factor near a vacuum-hBN-graphene-hBNsuperconductor heterostructure.In a and b the graphene Fermi energy has been set at 0.25 eV and T = 1 K for the solid curves and T = 94 K (above T c ) for the dashed curves; the graphene sheet is placed 4 nm above the superconductor surface.We show results for two emitter-graphene distances: 13 nm (a) and 36 nm (b).In c the Fermi energy of graphene has been set to zero and T = 1 K.The red curve corresponds to an emitter-superconductor separation of 17 nm and the blue curve of 40 nm.In d we show results for the same distances as in a (red curve) and b (blue curve), but for T = 94 K.In e we show how the Purcell factor depends on the graphene-superconductor distance t at the frequency of the Higgs mode, ω h = 2∆ ≈ 28.32 meV.The other parameters are kept fixed: E gr F = 0.5 eV, T = 1 K, and emitter-graphene distance of 13 nm.(The nonlocal optical conductivity of graphene was used.)

Figure 5 .
Figure 5. Purcell factor as a function of graphene's Fermi energy.Here we show the effect of changing graphene's Fermi energy (indicated at the top of each column) while keeping all other parameters fixed: T = 1 K, emitter-graphene distance (d = 13 nm for the top line and d = 2 nm for the middle line), graphene-superconductor distance t = 4 nm.For d = 13 nm, the dependence of the decay rate on the emitter's frequency changes qualitatively from low (E gr F = 50 meV) to high (E gr F = 250 meV) graphene doping.In the bottom panel we depict the differential LDOS given by the integration kernel of Eq. (3), solid lines, and by the analytical expression given by Eq. (5), dashed lines.The energy is fixed at the value ω H . (The nonlocal optical conductivity of graphene was used.)