Theory of resonantly enhanced photo-induced superconductivity

Optical driving of materials has emerged as a versatile tool to control their properties, with photo-induced superconductivity being among the most fascinating examples. In this work, we show that light or lattice vibrations coupled to an electronic interband transition naturally give rise to electron-electron attraction that may be enhanced when the underlying boson is driven into a non-thermal state. We find this phenomenon to be resonantly amplified when tuning the boson’s frequency close to the energy difference between the two electronic bands. This result offers a simple microscopic mechanism for photo-induced superconductivity and provides a recipe for designing new platforms in which light-induced superconductivity can be realized. We discuss two-dimensional heterostructures as a potential test ground for light-induced superconductivity concretely proposing a setup consisting of a graphene-hBN-SrTiO3 heterostructure, for which we estimate a superconducting Tc that may be achieved upon driving the system.


INTRODUCTION
2][3] Laser light has been shown to be a versatile tool in a variety of systems, capable of photoinducing ferroelectricity 4,5 , switching between charge density wave states [6][7][8] and even optically-stabilizing ferromagnetism. 9Arguably, one of the most exciting prospects is to use light to engineer high-temperature superconductivity.Experimentally, evidence for creating superconducting-like states in K 3 C 60 [10][11][12] and certain organic compounds 13 through laser driving in the THz range have made this prospect all the more tangible.
5][36][37] In this paper, however, we highlight the non-equilibrium properties of this model, which can potentially elucidate the microscopic mechanism behind photo-induced superconductivity.Furthermore, this model offers a simple, microscopic prescription for finding new nonthermal pathways to superconductivity.
In the following, we present a generic model of a boson coupled to the interband transition of two electronic bands and investigate its properties, focusing on driven superconductivity.We demonstrate that driving bosons enhances electron-electron attraction, particularly when the boson frequency matches the electronic band gap, resulting in a resonant enhancement.To support this claim we employ both perturbative methods and exact diagonalization for a two-site model, confirming the significant enhancement of electron-electron attraction on resonance.We then explore the consequences of this mechanism on superconductivity in extended systems by deriving the gap equation for our mechanism using manybody field theory.Consistent with the two site model results, we find that the resonance is absent in equilibrium, but the boson-mediated electron-electron attraction can be resonantly enhanced when the bosons are driven into a nonthermal state.
Finally we discuss a concrete system where our model can be tested experimentally with currently available techniques.We consider two-dimensional materials coupled to surface phonon-polariton modes of a substrate.As a practically realizable example, we explore a heterostructure consisting of graphene aligned with hexagonal boron nitride (hBN) [38][39][40][41] on top of the surface of SrTiO 3 .Graphene aligned with hBN hosts two electronic bands with a small gap, 42 while the SrTiO 3 surface provides surface phonon polaritons at a similar frequency as the gap. 43,44We estimate that photo-induced pairing in graphene could be found up to a critical temperature of around 15K.

Model for photo-induced superconductivity -
The generic model we propose for resonantly enhanced photoinduced superconductivity, has as basic ingredients two dispersive electronic bands: where c † k,σ creates an electron with momentum k and spin σ in the semi-filled conduction band with dispersion creates an electron with momentum k and spin σ in a higher lying conduction band with dispersion d (k) > c (k). Photo-induced superconductivity manifests when a boson is coupled to the interband transition via the Hamiltonian: where b q annihilates; b † q creates a boson with momentum q and frequency Ω(q) while g(q) parametrizes the cou-pling strength to the electrons, and X q = b † −q + b q is proportional to the coordinate of the oscillator.
Local approximation -To gain intuition for the microscopic mechanism that we propose, we first interrogate a local model consisting of lattice electrons where each site has two orbitals that are coupled locally via a boson.
Here c j,σ and c † j,σ are annihilation and creation operators of electrons in the lower level with spin σ at site j while d j,σ and d † j,σ are annihilation and creation operators of electrons in the upper level that is separated from the lower one by the energy gap ∆E.The local density operators for the two bands are defined as n c j,σ = c † j,σ c j,σ and n d j,σ = d † j,σ d j,σ .b † j and b j are bosonic creators and annihilators of the bosonic mode at site j that has eigenfrequency Ω and couples the two electronic levels with coupling strength g through the operator X j = b † j + b j .At this point we do not consider tunneling between the two sites.
Ground state attraction in the local model -We focus on the most simple two-site version, j ∈ {1, 2}, of the model in equation ( 3) and consider a half filled lower level with two electrons in the system overall.At vanishing coupling g = 0 the ground state of the system with no bosons is degenerate between two unpaired elec-trons on different sites and a singlet on the same site (E GS ↑,↓ = E GS ↑↓,− ).Upon turning on the coupling g > 0 this degeneracy is lifted: The state with a doubly occupied site obtains a slightly lower energy than two singly occupied sites which, to leading order in g reads We interpret this energy as an attractive interaction between the two electrons.Such an interaction has previously been noted in Ref. 34.
Attraction out of equilibrium -To account perturbatively for the out of equilibrium behaviour of the model, we perform a Schrieffer-Wolf transformation of the Hamiltonian in Eq. ( 3) for an arbitrary number of lattice sites to eliminate the coupling to leading order in g.As shown in the Methods section we find ( The occurring denominators (∆E 2 −Ω 2 ) −1 in the prefactors indicate resonant behaviour when ∆E ≈ Ω.We note that the resonance between the last two terms exactly cancels in the ground state, which can be seen by normal ordering the operators.We further analyze the second term Eq. ( 5) which constitutes an electronic density coupled to the squared boson displacement. 25In Ref. 14 such an interaction was considered on symmetry grounds and it was shown that it gives rise to a boson-number dependent attraction: This result indicates that driving the bosons out of equilibrium will enhance this attraction transiently and overcome the GS cancellation leading to a resonantly enhanced interaction.However, at the resonance, ∆E ≈ Ω, perturbation theory breaks down and we instead turn to exact numerical methods to investigate the resonance.
To explore the induced attraction away from the ground state, we perform exact diagonalisation calculations on the two-site model of Eq. (3) including a small hopping term parametrized by the hopping amplitude t First we compute the lower level doublon correlations in the ground state, as an indicator of an induced electron-electron interaction.The result is shown in Fig. 1(b).Consistent with our perturbative analysis in equation ( 4), in the ground state, we find positive correlations indicating an effective attraction that increases towards lower frequencies but show no specific features at ∆E ≈ Ω.
To analyse the out of equilibrium case, we coherently drive the bosons by adding a time-dependent term to the Hamiltonian according to Here F (t) is a Gaussian envelope and ω D the driving frequency that we always set resonant with the eigenfrequency of the bosons Ω = ω D .We compute the timeaveraged density correlations Eq. ( 8) as well as the timeaveraged boson number t1 t0 dt b † j b j (t) = N bos , where t 0 is a time shortly after the driving pulse and t 1 a later time after many driving periods.Details on how we perform the time evolution can be found in the methods section.We find an overall amplification of the electron-electron attraction compared to the equilibrium result.Moreover, while in equilibrium no enhancement of the interaction close to resonance Ω ≈ ∆E was found, such an enhancement can indeed be accessed out of equilibrium as can be seen in Fig. 1(b).This matches the intuition based on Eq. ( 6).The resonance is also evident in the timeaveraged number of bosons, which has a dip indicating that electrons are being excited into the upper level by absorbing bosons which would be detrimental to pairing.For completeness we note that for the two site model we observe a peak splitting at the resonance.The reason is that a finite hopping and coupling g leads to level splitting potentially enabling more than one resonance.Such features are expected to be smeared out in extended systems or due to dissipation.Dynamic gap equation.-Weshow that the resonantly enhanced photo-induced attraction persists in extended systems by using many-body field theory techniques.Starting from the Hamiltonian presented in Eq. ( 1) and (2), we integrate out the bosons exactly using the pathintegral approach yielding a dynamic electron-electron interaction.Then we perform a saddle point approximation to derive BCS type equations for the equilibrium phonon mediated superconductivity, and then discuss how they can be extended to the nonequilibrium case.
Within saddle point approximation, the intra-band pairing fields ∆ c = ψ k,↑,c ψ k,↓,c and ∆ d = ψ k,↑,d ψ k,↓,d are found self-consistently through two coupled gap equa-FIG.2. a.) Feynman diagrams for the coupled gap equations Eq. ( 10) and ( 11).b.) Critical temperature Tc as a function of frequency Ω for the sawtooth chain.The dispersion of the sawtooth chain is shown in the inset with ξ d = ε d − µ indicated.In thermal equilibrium (red line) no significant Tc is found while driving the bosons increasing their occupation by NDrive Eq. ( 16) yields a resonantly enhanced Tc for Ω ξ d where both the full expression Eq. ( 26) (yellow line) and the approximate expression Eq. ( 12) (blue dashed line) agree.c.) Nonequilibrium phase diagram as function of temperature T and driving strength nDrive (Eq.( 16)) in units of ξ d .At a critival driving strength the system changes from the normal state (NS) to a super conducting state (SC).The inset illustrates the non-thermal boson distribution we assume according to Eq. ( 16).The coupling g has been set to g = 0.2ξ d throughout while the chemical potential µ is set such that ξ d = 3t, where t is the nearest neighbour hopping.
tions shown diagramatically in Fig. 2(a): where we have introduced the Bogoliubov dispersion The dependence of the interaction on the state of the bosons is encoded in the frequency structure of the interaction.
In order to estimate T c we assume a local coupling g(q) = g and a flat bosonic dispersion Ω(q) = Ω.We next insert the equation for ∆ d into that of ∆ c arriving at a single self-consistent gap equation for ∆ c .We find that ∆ c (ω) exhibits only a weak frequency dependence for ω < Ω; hence, we disregard this dependence when estimating T c .This allows us to evaluate all frequency sums analytically and then continue the functions to real frequencies, introducing bosonic and fermionic occupation functions.The details of the calculation and the full gap equation are given in the methods section.
Similar to the standard BCS approach for supercon-ductivity in a single band, we only keep the most divergent terms contributing to the integrand of the gap equation.We identify two small parameters close to the Fermi surface when |ξ c | |ξ d | ≈ Ω, namely the detuning between the upper band and the phonon frequency,

and the dispersion relation of the lower band, |ξc|
Ω , and expand the gap equation to leading and next-to-leading order in these parameters.In this situation T c can be determined by a gap equation of the form: where the first term on the right hand side corresponds to a BCS type term, where ν(ξ) the density of states at energy ξ, with an effective electron-electron attraction given by U BCS = , which matches the attraction found from perturbation theory in the ground state (see equation ( 4)).Here n B is the Bose distribution function and n F is the Fermi distribution function and the integrals have a natural cut-off set by the bosonic frequency Ω, since the approximate expression for the gap equation is only valid within the energy window, |ξ c | < Ω and |ξ d − Ω| < Ω The second term, L Res , is given by: Intriguingly, L Res , has no BCS analogue and diverges even more strongly as ∆E(ξ d ) Ω , |ξc| Ω → 0. Similar to the BCS term the divergent integrand is cut-off by temperature.The appearance of a more strongly divergent term is a manifestation of the resonance found in our model when the frequency matches the energetic distance from the Fermi level to the upper band.However, at low temperatures T Ω, |ξ d |, pairing from L Res is exponentially suppressed because it is proportional to either the number of phonons or the number of electrons in the upper band through the factor n B (Ω) + n F (ξ d ).This leads to the absence of a resonance effect in equilibrium consistent with the two-site model.
We now conjecture that Eq. ( 12)-( 14) can be used to discuss non-equilibirum systems.This is accomplished by replacing the equilibrium, thermal distribution function n B (Ω) by a nonequilibrium distribution function.The main limitation of this procedure is that it does not include enhanced decoherence of electrons in the presence of photoexcited bosons. 30This decoherence provides a mechanism for pairbreaking and may result in a finite lifetime of photoinduced superconductivity. 15A nonequilibrium boson distribution activates L Res and reproduces a resonantly enhanced attraction consistent with our findings in the two-site model.
Model calculation in 1D -We use the sawtooth chain 45,46 to illustrate the interband phonon mechanism on a concrete example.This one-dimensional model consists of a dispersive lower band and a flat upper band: where a is the lattice constant, t the hopping between sites, and µ the chemical potential.The flat upper band in the sawtooth chain gives rise to a constant detuning, providing an idealised system to identify the resonance present in our model.We explore the possibility of using the nonequilibrium boson distribution as a resource by switching on the resonant contribution in the gap equation Eq. ( 12) through adjusting the equilibrium boson distribution to a driven one as illustrated in Fig. 2(c): We compute T c as a function of the boson frequency Ω and explore the resonance by tuning Ω across ξ d .The result is shown in Fig. 2(b).In equilibrium we do not find a significant T c for the coupling of g = 0.2ξ d and the selected frequencies.However, a significant T c appears close to the resonant region as soon as we assume a nonthermal boson distribution according to Eq. ( 16).Both the full expression of the gap equation shown in the methods section and the approximate expression in equation ( 12) match precisely close to the resonance.
Unlike in the two-site model, in the extended system the driven-boson-induced interaction turns into repulsion for blue-detuned frequencies Ω ≥ ξ d .8][49] We note that in the blue detuned region, real transitions from the lower to the upper band, that are not captured within our approach, become possible which are expected to be detrimental to pairing. 15The resonant enhancement is found only for a red detuned phonon where direct electron transitions are prohibited by energy conservation, hinting that a more involved nonequilibrium treatment of the problem within the Keldysh formalism, that takes the change of Fermi functions upon driving into account, is expected to exhibit the same qualitative behaviour.
In order to summarize our findings for the driven sawtooth chain, we compute an out-of-equilibrium phase diagram as a function of temperature and driving strength (Fig. 2c), fixing Ω slightly below resonance, Ω = 0.99ξ d .While the system is in the normal state down to very low temperatures in equilibrium (N Drive = 0) for g = 0.2ξ d , populating the bosonic mode gives rise to a finite T c .
Graphene-hBN-SrT iO 3 heterostructure.-Weconclude by proposing an experimentally realizable platform for resonantly induced superconductivity, namely a heterostructure of graphene aligned with hBN placed on top of bulk SrTiO 3 , as shown in Figure 1(c).Graphene aligned on hBN leads to a gap at the Dirac points of graphene, with size ∆E = 14meV. 42The bosonic mode within our pairing model is provided by a low energy surface-phonon-polariton 43,44 that exists in the quantum paraelectric SrTiO 3 and couples to the interband transition in gapped graphene through its dipole moment. 50,51e compute the relevant light-matter coupling using the mode functions obtained in Ref. 52 for the surface polaritons, and approximate the resulting coupling by a box function g(q) = g θ(|q| − 1 d ), where d is the distance between the graphene and the SrTiO 3 surface which provides a natural short-wavelength cut-off for surfacephonon-polariton induced interactions.We assume an hBN thickness of 5nm, corresponding to approximately 15 atomic layers. 53We include the local Coulomb repulsion that can be estimated via density functional theory 54 and take into account Morel-Anderson renormalisation processes that reduce its value to U * ≈ 1.1eV which is the value we use for our estimate.Using these inputs the T c estimate is obtained by solving the full gap equation (26).
Without any driving the Coulomb repulsion prevents pairing for reasonable temperatures.When driving the system into a state with a polaritonic occupation n Drive = 1, we obtain T c = 15.2K-a temperature that is readily accessible.Assuming that only polariton modes with momentum q < 1 d are excited, the required polaritonic occupation corresponds to exciting 3 × 10 −4 bosons per graphene unit cell.
In the Methods section, we estimate the change in temperature due to heating that the sample would undergo if all of the energy of the initial excitation was converted into heat close to the surface of the SrTiO 3 substrate.We find that at the critical temperature T c = 15.2K the expected temperature change would be less than 1K.This emphasizes that the driving we propose is of realistic strength.Surface modes cannot be directly excited by a laser due to their dispersion lying outside the lightcone.However, these modes can be driven by irradiating a tip, which generates excitations in the near-field that effectively drive the modes. 55 final comment is in order regarding the assumed pairing symmetry of resonantly induced superconductivity.In this work we considered only s-wave superconductivity.In reality, due to the electronic Coulomb repulsion, d-wave pairing might be favoured, and our mechanism can equally well lead to or contribute to d-wave superconductivity.In such a scenario, the negative effect of the local Coulomb repulsion would be further mitigated.
Discussion.-We introduced a simple microscopic mechanism based on a boson coupled to an electronic interband transition.We demonstrated both photoinduced enhancement of superconductivity, as well as resonant amplification of this enhancement when the phonon frequency is red-detuned but nearly aligned with the electronic excitation energy between the Fermi level in the partially filled lower band and the empty upper band.We demonstrated the existence of the resonantly amplified electron attraction by both exact diagonalization for a two-site model and by deriving a gap equation for extended systems.Tantalizingly, we showed that this interaction can be engineered in a two-dimensional material on a surface-phonon-polariton platform.For an example heterostructure of graphene-hBN-SrTiO 3 , we computed that superconductivity can arise at temperatures around T c ≈ 15K.
We believe that our simplified treatment of the nonequilibrium problem correctly captures the renormalization of the effective electron-electron interaction due to a non-thermal distribution of bosons.However, we do not exclude that heating and pair-breaking effects counteracting superconductivity 15,17 might be underestimated here, in particular close to the resonance, which might create a trade-off that needs to be taken into account when tuning the frequency.Moreover, dynamical generation of disorder that was theoretically shown to limit the possibility to realize light induced super-conducting like states in 1D systems, was not explored here either. 56t is worth noting, however, that in 2D systems, several studies have argued that disorder leads to multifractality that actually increases the transition temperature of superconductivity. 57The interplay between disorder and resonantly photo-induced superconductivity is thus left for future research.
Interestingly, even though the enhancement of the effective electron-electron attraction is likely a transient phenomenon, the superconducting state arising from it might be longer lived than the initial drive, as indicated by the accompanying study of the authors on the same model away from the resonance. 58n view of the resonant nature of our mechanism, an interesting future direction is to explore its applicability in the context of cavity materials engineering, [59][60][61][62][63] for which several proposals for cavity-induced superconductivity already exist. 64 using as transformation matrix 69 with α and β to be determined by the condition This condition is satisfied when The effective Schrieffer-Wolff Hamiltonian to order g 2 are calculated via 1 2 [S, V ] which gives rise to equation ( 5) in the results section.
Attractive correlations in 2-site model.-Inorder to explore the resonant regime Ω ≈ ∆ for the model in Eq. ( 7) we perform full diagonalization calculations.We limit the bosonic part of the Hilbert space to a maximumm of N max = 6 bosons for which all quantities converge.To explore the out of equilibrium properties, we add a driving term to the Hamiltonian according to Eq. ( 9).In Eq. ( 9), F (t) has a normalized Gaussian shape centered around t pump = 16π ωD and standard deviation s = 4π ωD multiplied by a bare driving strength F 0 that is set to F 0 = ∆E 2 .We always drive the bosons on resonance setting ω D = Ω.We evolve the system in time with the time evolution operator U (t) = T exp i t1 t0 H(t)dt , where T is the time-ordering operator, using finite time steps of width ∆t = π 10ωD , i.e. 20 time steps per driving period.We compute the time-averaged local doublon correlations where t 0 = 30π ωD and t 1 = 130π ωD , reported in Fig. 1b together with the time-averaged boson number.We show a concrete time-evolution for longer times in Fig. 3 for the resonant case Ω = ∆E and a slightly red-detuned case ω D = 0.9∆E.One can see that all quantities undergo heavy oscillations that are not decaying attributed to quantum revivals for this very small system that we don't expect to be present in extended systems.However, the time-averaged correlations quickly converge to a finite value.We additionally mark the time interval over which we computed the time average in the main part [t 0 , t 1 ].We average over 50 driving periods after the pump has decayed.The reason we limit the averaging time is that in a real system with dissipation the excitation would decay -presumably not reaching 50 oscillation cycles.Reducing the number of cycles for the averaging process does not qualitatively change our results.
Effective electron-electron interaction.-Weuse the path-integral formalism to derive a dynamic electronelectron interaction from the model Eq. ( 2).The bare action of the system, including bosonic modes, can be written as where k = ( k, iω) and q = ( q, iΩ) are composite indices of a lattice momentum k or q and a fermionic(bosonic) Matsubara frequency ω(Ω).The b q are complex field, bq their complex conjugate and D(q) −1 = iΩ − Ω( q) is the inverse bosonic propagator.ψ k,σ,n and ψk,σ,n are Grassmann valued fields with spin index σ ∈ {↑↓} and band-index n ∈ {c, d} and is the inverse bare electronic Green's function.We introduce an inter-band coupling via the bosonic fields as Assuming a 2-band model we have introduced n as the band complementing that labelled by the band index n.g( q) is the q dependent inter-band coupling.Now we can integrate out the bosons with a suitable substitution of the bosonic fields in the path integral to arrive at an effevtive interaction ψk −q,σ ,m ψ k ,σ , m  .
(24) Due to the sum over q and since all other terms are even in q one may perform the replacement taking only the part of the propagator even in q.Equations (10) and (11) are derived using the saddle point approximation in the Cooper channel, leading to a selfconsistent equation shown schematically in Fig. 2(a).A key aspect of our work is performing the Frequency sums analytically and expressing the result in terms of occupation distributions of the bosons and fermions.In the main text we show an approximate formula of the most divergent terms, for completeness we report here the full gap equation.We take the full ω and k dependence of ∆ d by replacing the expression for ∆ d in equation (11)  into equation (10), to arrive at a self-consistent equation for ∆ c only.We then find that the expression for ∆ c is weakly frequency dependent for ω Ω which allows to approximate ∆ c as indpendent of ω and k.The resulting gap equation that determines T c is found to be:  9).The times t0 and t1 used for the time average in Eq. ( 21) are marked in the bottom panel where the pump envelope according to Eq. ( 9) is shown.Additionally we show the time averaged correlations as a function of the upper time of the averaging.This quantity, which underlies Fig. 1 b.) in the main part, quickly converges.
Graphene coupled to surface phonon polaritons.-Inthis part we outline the estimate of the pairing temperature of our proposed heterostructure consisting of Graphene aligned with hBN on top of bulk SrTiO 3 .We use the q-dependent surface phonon-polariton modes given in Ref. 52 that read where the upper half-space z > 0 is assumed to be filled with hBN with dielectric constant ε hBN = 5 70 while the lower half-space is made up of the SrTiO 3 with ε TO −ω 2 .We use ε ∞ = 6.3 71 , ω TO = 2π×1.26THzand ω LO = 2π ×5.1THz. 72The normalization N of the modes can be determined from the normalization condition 52 with v(ω) = 2+ω∂ ω ε(ω) and ε 0 the vacuum permittivity.The limiting surface phonon frequency when c 2 q 2 ω 2

LO
is Ω s = 2π × 3.9THz which we will use throughout for our estimate.In order to estimate the light-matter coupling we focus on the Dirac cone at which the Hamiltonian, including the gap opening due to an effective onsite potential induced through the proximity of hBN to the graphene layer, reads where V 0 is the onsite potential, v F the Fermi velocity and A the vector potential due to the surface polaritons.At the K-point the inter-band coupling is equal to g inter = 1 while far away from the K-point it can be estimated from the Dirac Hamiltonian without including an onsite potential which yields g inter = i vF(k×A)

|k|
. Using the coupling away from the K-point as a conservative estimate, we get for the inter-band coupling where f (ω) = 4 ω 2 (ε∞+εHBN) 2 , S cell is the area of the unit cell in graphene, α the fine-structure constant and c the speed of light.In order to obtain a simple estimate we approximate the coupling by a box function g k,q → gθ |q| − 1 d .The cutoff 1 d is naturally provided by the mode-function Eq. ( 27) since ε(ω) → −ε hBN for q 2 c 2 ω 2 LO and therefore |θ(q)| 2 ∼ e −2qd .We determine g by fixing the integral g2 where we use g 2 since this appears in the interaction.Due to the box-function the coupling has now attained a q dependence in principle leading to a k-dependence of the gap itself.We compute the gap at k = K and assume it to be zero outside of the coupling radius θ(k).Additionally, the momentum transfer k − k is limited due to the box coupling by a factor θ(k − k ).We set the chemical potential µ = −15meV slightly into the lower band and otherwise the parameters from the main text.We include the local Coulomb repulsion which was estimated to be U ≈ 9.3eV based on density functional theory calculations. 54The bare local repulsion will be reduced by Morel-Anderson renormalisation according to 73 which we take as the value of the effective Hubbard U for our T c estimate.Importantly, the resulting T c does not depend very sensitively on the precise value of U * .Due to the exponential dependence of the polariton mode field strength on the distance d from the surface, T c does, however, vary significantly with d.In particular, short distances are very desirable to achieve a higher T c .We evaluate Eq. ( 26) on a k-grid of size 500 × 500.For the driven case we assume a non-thermal boson occupation according to Eq. ( 16), setting n Drive = 1 for all modes with q < 1 d .Estimated heating due to the drive In order to quantify the strength of the drive we propose to induce superconductivity in the Graphene-HBN-STO heterostructure, we provide an order of magnitude estimate of the heating of the sample that would result from the drive.The surface modes of STO couple to bulk modes through phonon non-linearities and disorder and therefore heat up the substrate.For a first estimate, we assume this to be the dominant heating effect and hence estimate the heating of the STO substrate assuming all energy of the initial excitation to be converted into heat.Shortly after the drive, the surface mode will only heat up a thin layer of the substrate close to the surface corresponding to the penetration depth which sets the length scale over which the surface modes couple to the bulk.At later times the heat will dissipate.We take ξ = 5nm as the penetration depth corresponding to the penetration depth of the mode with the largest q of all modes that we assume to have a non-zero coupling to graphene (see previous section).This constitutes a conservative estimate since modes with smaller q have a larger penetration depth.To get an estimate of the energy per volume, we compute expected number of excitations within the area of one unit-cell of graphene n ex = BZ b † q b q dq BZ 1 dq = BZ n Drive Θ 1 d − q dq BZ 1 dq ≈ 3 × 10 −4 (33)  where in the second step we have assumed that only modes with small wave-vectors are excited for which we chose the same cutoff as for the momentum coupling and integrals run over the Brillouin zone.The expected change in temperature is then where V STO is the unit-cell volume of STO, A C the area of the unit-cell of graphene, N av the Avogadro number and c P the specific heat of STO at constant pressure.We take the specific heat of STO at 15K from Ref. 74 as c p /T 3 ≈ 0.75×10 −4 J mol K 4 and use as the lattice constant a STO = 0.3905nm.Since STO has cubic lattice structure we have V STO = a 3 STO .Put together we obtain a change in temperature of ∆T = 0.7K.(35)  On the other hand, for an estimate of an upper bound on the induced heating, one may assume that all energy of the initial drive is converted into heat solely in the graphene flake.For this estimate we use the heat capacity of graphene c P = 0.24 J mol K at T = 10K taken from Ref. 75 and compute ∆T = Ω s n ex c P ≈ 1.9K.(36)   This is likely an overestimate, but given that we assume a temperature of the sample of 10K and computed T c ≈ 15K, the proposed drive would not destroy the superconductivity solely by heating the sample.

FIG. 3 .
FIG.3.Time evolution of N bos (blue line, right axis) and C (yellow line, left axis) red detuned from the resonance Ω = 0.9∆E and at resonance Ω = ∆E for long times of the system defined by Eq. (9).The times t0 and t1 used for the time average in Eq. (21) are marked in the bottom panel where the pump envelope according to Eq. (9) is shown.Additionally we show the time averaged correlations as a function of the upper time of the averaging.This quantity, which underlies Fig.1b.) in the main part, quickly converges.
Angel Rubio, John Sous, Andy Millis, Martin Eckstein and Hope Bretscher.M.M. is grateful for the financial support received from the Alex von Humboldt postdoctoral fellowship.S.C. is grateful for support from the NSF under Grant No. DGE-1845298 & for the hospitality of the Max Planck Institute for the Structure and Dynamics of Matter.ED acknowledges support from the ARO grant "Control of Many-Body States Using Strong Coherent Light-Matter Coupling in Terahertz Cavities" and the SNSF project 200021 212899.