Transport through vertical graphene contacts under intense laser fields

We theoretically study the electronic and transport properties of two graphene layers vertically coupled by an insulating layer under the influence of a time-periodic external light field. The non-adiabatic driving induces excitations of electrons and a redistribution of the occupied states which is manifested in the opening of gaps in the quasienergy spectrum of graphene. When a voltage is applied between the top and bottom graphene layers, the photo-induced nonequilibrium occupation modifies the transport properties of the contact. We investigate the electronic and transport properties of the contact by using the nonequilibrium Green's function formalism. To illustrate the behavior of the differential conductance of the vertical contact under the light illumination, we consider two cases. First, we assume that both the bottom and top layers consist of graphene and second we consider a finite mass term in the bottom layer. We obtain that the differential conductance is strongly suppressed due to the opening of gaps in the quasienergy spectrum in graphene. Additionally, the conductance shows features corresponding to the tunneling of photoexcited electrons at energies of the van Hove singularity for both the top and bottom layers. In the case of a finite mass term in the bottom layer, the differential conductance can be directly related to the tunneling of photoexcited electrons.


I. INTRODUCTION
The outstanding mechanical, optical and electronic properties of graphene make it an attractive material for next-generation technology 1 . Prominent examples are applications of graphene in optoelectronic devices such as photodetectors 2-4 and sensors [5][6][7] . The core of the improved functionality in graphene-based optoelectronic devices lies in the interaction of light and matter in low dimensions. Although a remarkable absorption of 2.3% of incident light in monolayer graphene 8 is too weak to realize high-performance devices, several methods have been proposed to tune and enhance light-matter interaction in graphene. Among them, the absorption of graphene is enhanced by exciting surface plasmonic resonances [9][10][11][12] and by utilizing resonant structures to couple light with graphene [13][14][15] .
Apart from employing light-matter interaction in graphene to enhance the functionality of optoelectronic devices, graphene reveals several fundamental lightinduced phenomena 16 . A time-periodic external perturbation can turn a topologically trivial bandstructure into a topological one 17,18 . Depending on the polarization of the light field, gaps appear both at the Dirac points and other momenta in the Brillouin zone [19][20][21][22][23][24][25][26][27][28][29] . The most prominent gap opens at energies ε = ω/2 corresponding to a resonant absorption/emission of a photon with frequency ω between the valence and conduction band of graphene. Accompanied by the opening of a gap in bulk graphene, are chiral edge states. In particular, it has been discussed that nonlinear driving in graphene induces chiral edge photocurrents 30 and the light-induced Hall effect without a magnetic field 21,31 .
The periodic driving is in general studied within the Floquet formalism which is equivalent to Bloch theorem for periodic perturbation in time 32 . Besides the open-ing of a gap, the temporal periodicity induces copies of the original electronic bands spaced by integer multiples nω, the so-called Floquet sidebands. Floquet bands have been observed by time-and angle-resolved photoemission spectroscopy of surface Dirac fermions in a topological insulator under circular light irradiation 33,34 . An alternative way to probe Floquet bands and the related opening of a gap in the band structure, is by conductance measurement. Signatures of light-matter interaction in the conductance have been studied in planar two-terminal 19,21,23,[27][28][29] and multiterminal contacts [35][36][37] . The most dramatic effects of light-matter interation are manifested in a suppression of the conductance. In contrast to a time-periodic external light field, time-dependent modulation of gate or contact potentials and the accompanied photo-assisted tunneling of electrons lead to a great variety of interesting phenomena [38][39][40][41][42][43] .
Motivated by the growing interest in nonequilibrium driving of graphene devices, we study a contact consisting of two graphene layers vertically separated by an insulating layer. A voltage is applied independently between the top and bottom layer allowing for tunneling between the layers. A similar setup consisting of a graphene/boron nitride/graphene heterostructure was discussed in Ref. 44 and it was demonstrated that resonant tunneling with the conservation of energy and momentum can be achieved by aligning the crystallographic orientation of the top and bottom graphene layers. In addition to the applied voltage, we assume that the top layer is driven to a nonequilibrium state by a time-periodic monochromatic light. We discuss two different scenarios. First, we assume that the vertical contact consists of two graphene layers. Second, we assume that a finite gap opens in the bottom layer. In our approach, we focus on noninteracting electrons in which Coulomb interaction is neglected.
We study effects of both linear and circular polarization on the density of states, the occupation, and the differential conductance.
The main results are the following. The differential conductance is strongly suppressed in a broad range of voltages. In a qualitative picture, the suppression of the differential conductance is related to opening of gaps in the quasienergy spectrum due to time-dependent driving. Various aspects of the suppression of the conductance due to light-matter interaction have been discussed in Refs. [19,21,23,[27][28][29][35][36][37]. Here, we extend these previous results of the differential conductance for vertical graphene contacts and the regime beyond the linear Dirac approximation.
For strong driving and high frequencies, the momentum dependence of the optical matrix elements and the non-linearity of the graphene bandstructure becomes apparent in the differential conductance. Additional features occur at voltages that correspond to tunneling of photoexcited electrons to energies of the van Hove singularity in the top or bottom layer. For the case of a finite mass term in the bottom layer, the signatures in the differential conductance can be directly related to the tunneling of photo-excited electrons. Although we restrict the discussion to light absorption in graphene, the approach can be extended to other two-dimensional materials 45 .
The article is structured as follows. In Sec. II, we introduce the continuous Hamiltonian and discuss the individual contribution of the Hamiltonian. From the continuous Hamiltonian, we derive a tight-binding Hamiltonian describing the light-matter interaction and the tunneling between the contacts. In Sec. III we describe the nonequilibrium Green's function formalism and the relations of the Green's functions to the physical properties. Sections IV-VI contain the results of our article. In Sec. IV and V we study the density of states and the occupation of graphene under light illumination, respectively. Finally, in Sec. VI, we discuss the differential conductance and show that the light-matter interaction induces suppression and enhancement of the conductance. In Sec. VII, we summarize and conclude.

II. MODEL
We consider a vertical contact consisting of a top (t) and a bottom layer (b) which are coupled by an insulating layer described by a tunneling element t as shown in Fig. 1. A voltage V is applied between the top and the bottom layer and an external electric field perpendicular to the layers induces a nonequilibrium occupation of the electronic states. For simplicity, we assume that the electric field interacts only with the top layer and that the sample size is much smaller than the wavelength of the light such that we can regard the light field as being uniform in space.
The Hamiltonian of the vertical junction is given bŷ with the HamiltoniansĤ t (t) of the top layers,Ĥ b of the bottom layer and the tunneling HamiltonianĤ tun . In the following we discuss the contributions toĤ(t) separately. We introduce a continuous time-dependent Hamiltonian H t (t) of the top layer under light illumination and discretize it to obtain a tight-binding Hamiltonian. Similar to the top layer, we then introduce the Hamiltonian of the bottom layer and the tunneling Hamiltonian.

A. Hamiltonian of graphene under light illumination
The Hamiltonian describing the motion of electrons with mass m in the periodic lattice potential of graphene under light-radiation can be written aŝ with the field operators of the upper layerΨ † t (x) and Ψ t (x). The periodic lattice potential of graphene is composed of two triangular sublattices potentials A and B, respectively. The light-matter coupling is obtained from the minimal substitutionp →p − eA(t) and the vector potential A(t) is modeled classically in the dipole approximation. In the following studies we consider the case of linear polarized light A(t) = (A x cos(ωt), 0, 0) and circlar polarized light A(t) = (A x cos(ωt), A y sin(ωt), 0) with the frequency ω of the photons.
We then transform the continuous HamiltonianĤ t (t) into a tight-binding Hamiltonian and expand the field op- erators in terms of the carbon 2p z wave functions Φ t (x), with the annihilation operatorsâ kt andb kt of a particle with momentum k t on the sublattice A and B in the top layer, respectively. The tight-binding Hamiltonian is separated into a bare HamiltonianĤ 0 and a HamiltonianĤ I (t) descibing the interaction of graphene with the optical field,Ĥ t (t) = H 0 +Ĥ I (t) . It is convenient to transform the operators in sublatttice space to the space of the conduction and valence electrons. We write the creation operator of a quasiparticle in the conduction (c) and valence (v) band asψ † kt = (â c kt †âv kt † ). The bare Hamiltonian can then be written asĤ with the eigenenergieŝ and The Dirac point energy ε D sets the energy of the Dirac point relative to the chemical potential µ. The interaction energy between the next-nearest neighbor carbon atoms is given by γ 0 and the chemical potential is µ t . In the above equation with the vectors δ i connecting neighboring carbon atoms. These vectors are given by δ 1 = a 2 (1, √ 3), δ 2 = a 2 (1, − √ 3) and δ 3 = (−1, 0) with the nearest neighbor distance between the atoms a = 1.42Å. The vectors δ i together with the lattice of graphene are shown in Fig.2(a).
Following Refs. [46] and [47], we work in the Coulomb gauge ∇·A = 0 and apply the dipole approximation since the momentum of the photon is negligible compared to the momentum of the electrons. The light-matter interaction can then be written aŝ with the light-matter interaction matrix and the elements The time independent elements M k are with φ k = arg[f (k)] and a constant factor The light-matter interaction can be divided into intraband coupling, M cc k l (t) and M vv k l (t), and interband coupling, M cv k l (t) and M vc k l (t). The intraband coupling describes processes in which an electron in the valence (conduction)-band absorbs/emits a photon with frequency ω and stays in the valence (conduction)-band. In a process of interband coupling, an electron scatters from the valence (conduction)-band to the conduction (valence) band and exchanges energy of the photon.
In the following, we combine all the constant values of Eq. (11) in the amplitude of the vector potential A(t) and introduce the coupling strength in units of energy. To get an estimate of the magnitude of the vector potential, we consider the case of linear polarized light and convert the vector potential A(t) to an electric field by E(t) = ∂A(t)/∂t. In Ref. 31, a graphene sample was illuminated with an electric field amplitude E x = 4.0 × 10 7 V/m 31 and a photon energy ω ≈ 191 meV. Introducing the interaction energy of graphene γ 0 = 2.6 eV of the next-nearest neighbor carbon atoms, the photon energy corresponds to ω ≈ 0.07γ 0 . The coupling strenght in Eq. (12) then is of the order α x,y ≈ 0.02γ 0 with the factor M 3.0 nm −147 .
For the purpose of presentation, we show in following sections the results for slightly larger electric fields and photon energies. For example, in Secs. IV-VI, we assume an electric field amplitude of the order of E x = 0.1 V/Å and a photon energy ω ≈ 1.7 eV corresponding to α x,y ≈ 0.1γ 0 and ω ≈ 0.6γ 0 . We remark that high frequencies pulses with energies ω ≈ 1.7 eV have been for instance applied in dielectrics to control the electric current by the electric field 48 . Although in our approach, we assume a constant external field, high electric field amplitude can only be achieved by laser pulses. However, as pointed out in Ref. 31, a constant field amplitude is still a valid approach if the envelope of the Gaussian pulse varies slowly in time 31,49,50 .

B. Hamiltonian of the bottom lead and tunneling Hamiltonian
The electrons can tunnel to the bottom layer which we model as graphene with a finite mass term. The Hamil- tonian can be written aŝ with the excitation energies The tunneling Hamiltonian is given bŷ with the tunneling elements t ktk b . In the following we assume momentum conservation at tunneling and set t kk = δ kk t.

III. METHOD
Since we are interested in nonequilibrium transport properties, we apply the Keldysh Green's function approach to the HamiltonianĤ. In this method we express the physical properties in terms of Green's functions and calculate the corresponding Green's function. 51,52 To this end, we define the contour-ordered Green's function with the indices λ = (c, v) and λ = (c, v) indicating either the conduction band (c) or the valence band (v). The times τ and τ are located on the Keldysh contour and the contour-ordering operator is denoted as T .
To simplify the notation, we collect the Green's functions of the valence and conduction electrons in a matrix G k (τ, τ ) defined bŷ The general procedure then is to derive the Dyson equation for the Green's function in Eq. (16). This Green's function is given bŷ corresponding to the HamiltonianĤ 0 . In the next step, the Dyson equation has to be transformed from the contour time to the real-time. 51,52 After the transformation to the real time, every element in Eq. (16) becomes a matrix in Keldysh space defined by In the above expression, the upper indexes 1 or 2 refer to the position of the times t and t on the Keldysh contour. A definition of all the Green's functions in Eq. (19) and their relations is given in Appendix A. From the elements in Eq. (19) we can derive the retarded Green's

The advanced Green's function is related to the retarded Greens function by its hermitian conjugateĜ
The retarded (R), advanced (A) and the lesser (12) Greens functions constitute the building blocks by which we will express the charge current in the following section.

A. Green's functions without light irradiation
Following the procedure that we describe, we derive the Green's function for the bottom lead. Since we assume that the light only interacts with the upper layer, these Green's functions correspond to the Hamiltonian in Eq. (13) with momentum k b . In Fourier space, the retarded Green's function is given by , with a infinitesimal small real part η. The lesser Green's function iŝ with the chemical potential on the bottom layer µ b .

B. Green's functions under light irradiation
In this subsection, we derive the Dyson equation for the retarded and lesser Green's functions under lightirradiation. Because of the time-dependence in the optical field, we apply a discrete Fourier transformation to the Green's functions. To keep the notation simple, we replace the momentum k t on the top layer with k.
From Eq. (17) and the relations between the Green's functions in Keldsyh space we can derive the Dyson equation for the retarded Green's function The Green's functions in Eq. (20) depends only on momentum k due to the dipole approximation in the lightmatter interaction and the accompanied momentum conservation. Similarly, we can derive the lesser Green's function To compute the Green's functions in Eqs. (20) and (21) we have to apply a discrete Fourier transformation due to the time-dependence of the light-matter interaction. We can solve the Dyson equations by first transforming the time arguments to the Wigner representation with the relative time τ = t − t and the center-ofmass time T = (t + t )/2. 52 Second, we apply the Fourier transformation 51 (22) and then numerically solve the Dyson equation by matrix inversion. To simplify later equations, we introduce the notation for the Fourier transformed Green's function G R,12 k,nm (ε) =Ĝ R,12 k,n−m (ε + (n + m)ω/2). After the Fourier transformation of Eq. (20) we obtain In Eq. (24) A ± refers to the components of the vector potential A(t) = A + e iωt + A − e −iωt which are proportional to positive (+) and negative (−) frequencies, respectively. For the case of linear polarization the components are A ± = (A x /2, 0, 0) and for circular polarization A + = (A x /2, 0, 0) and A − = (0, −iA y /2, 0). The unperturbed retarded Green's functions in Eq. (23) of the conduction and valence electrons are given by Applying the Fourier transform of the lesser Green's function in Eq. (21) giveŝ The unperturbed lesser Green's function is given bŷ

C. Charge current
With the Green's function derived in the previous sections, we can calculate the charge current. The current on the bottom lead can be expressed as with the trace over the Green's functions in the space of the conduction and valence electrons. In Eq. (26) G 12 ktk b (t, t) is the Green's functions corresponding to the full HamiltonianĤ(t) taking the coupling between the leads into account to infinite order in the tunneling. However, to simplify the calculations, we restrict the discussion to the tunneling limit for which we can write the current as with with t kk = δ kk t and the momentum-independent tunneling element t. Since the summations over k t and k b separate in the tunneling limit, we can compute the summation over k b in the bottom layer analytically 53,54 . In Appendix B, we present a recursive method to compute the Green's functionsĜ R,A,12 k,nm .

IV. DENSITY OF STATES UNDER LIGHT IRRADIATION
In this section, we discuss the density of states and quasienergy spectrum of the graphene under light irradiation. We derive an effective Hamiltonian to study the splitting of quasi-eigenstates of graphene due to the lightmatter interaction. The interaction of graphene with the light induces sidebands of the conduction and valence band, which are shifted by ±nω from the eigenenergies ε k of the bare graphene Hamiltonian. In the following, we restrict the discussion to a case of n = 0, ±1, ±2 sidebands, since the effect of light on the quasienergy spectrum decreases with the number of sidebands.
The quasispectrum of graphene under light irradiation has been discussed in literature using Floquet theory 19,21,23,[27][28][29]32,[55][56][57][58] . For completness and discussions in later sections, we present here the results of the density of states of the band structure for momenta in the first Brillouin zone using the approach discussed in Sec. III.
The momentum-dependent density of states ρ k (ε) of graphene under light irradiation is given by with the retarded Green's function in Eq. (23). Figure 3 shows the density of states of graphene along the path (K − 0.7ŷ) → K → (K − 0.7K) in the first Brillouin zone shown in Fig. 2(b). The unit vectors in y and Kdirection are written asŷ andK, respectively. We discuss both the case of linear polarized light [ Fig. 3(a)] and circularly polarized light [ Fig. 3(b)]. Since a typical graphene sample is doped, we set the Dirac point energy to ε D = −0.2γ 0 such that states both above and below the Dirac point energy are occupied in an equilibrium state. The Dirac point energy in Fig. 3 is indicated by the dashed line. The interaction of light with graphene induces opening of gaps, the most prominent one is called the dynamical gap and occurs at energies ε k = ω/2 + ε D . The transition corresponds to a resonant absorption/emission of a photon with energy ω between the conduction and valence band. Importantly, the matrix elementsM ± k in Eq. (24) depend on momentum and on the polarization of the light. For certain momenta and polarization, the matrix element vanishes and a gap opening is suppressed. In the case of circularly polarized light, the quasienergies become fully gapped at the Dirac point energy. We remark that the bandgap opening in graphene due to circularly polarized light infers a topological phase transition from a trivial to a topological band structure 17,18,21,35 .
To further study the effect of the light-matter interaction on graphene, we restrict the discussion to a minimal model, in which we consider the emission or absorption of one photon such that n = 0, ±1. An effective Hamiltonian can be derived by comparing the retarded Green's function in Eq. (23) with the definition of a retarded Green's functionĜ R k,nm (ε) = (ε + iη −Ĥ eff k,nm ) −1 . We   Fig. 2(b). The unit vectors in y and K-direction are written asŷ andK, respectively. The dashed line shows that Dirac point energy which is set to εD = −0.2γ0. The frequency is ω = 0.6γ0, n = 0, ±1, ±2 and η = 10 −4 γ0. In (a) we consider linear polarized light with A(t) = (Axcos(ωt), 0, 0) and αx = 0.1γ0. In this case the element ∆ k vanishes along the path (K −0.7ŷ) → K and the spectrum remains ungapped. In (b), we consider circular polarized light with A(t) = (Axcos(ωt), Aysin(ωt), 0) and αx = αy = 0.1γ0. The density of state shows gaps at finite energies and becomes fully gapped at the Dirac point energy. The splitting close to the energies ε k = ω/2 + εD is given by ∆ k .
can then write the effective Hamiltonian aŝ Assuming thatM ± k = 0, the eigenvalues ofĤ eff k,nm are given by E c,v = ±ε k and the sidebands E c,v + = ±ε k + ω and E c,v − = ±ε k − ω which are shifted by ±ω from the original conduction and valence band, respectively.
To study the dynamical gap, we simplify the effective Hamiltonian and only consider the conduction band coupled with the first sideband of the valence band. Such an effective model corresponds to a 2x2 matrix inĤ eff k and describes the gap at finite energies. The eigenvalues for   Fig. 3. The red dashes line shows the boundary of the first Brillouin zone. In (a), the element vanishes along the path M → K and the density states remains ungapped as shown in Fig. 3(a). In (b) the the optical matrix element ∆ k is finite when approaching the K point and the density of states is fully gapped as shown in Fig. 3(b).
such a system are given by with the off-diagonal matrix element ∆ k = |Im[M k ]A ± | being proportional to the absolute value of the imaginary part of the light-matter matrix element connecting the valence and conduction band. Close to energies ε k ω/2 + ε D , the splitting between the eigenvalues is therefore 2∆ k and hence the gap is proportional to the coupling strength α x,y . The general behavior of gap openings under light illumination has been discussed for instance in Refs. 27 and 59.
Since the dynamical gap is proportional to ∆ k we can understand the gap opening at finite energies in Fig. 3 by consider the momentum dependence of ∆ k . Fig. 4 shows the momentum dependence ∆ k for (a) linear and (b) circular polarized light with the dashed line denoting the first Brillouin zone. When changing the momentum from (K − 0.7ŷ) → K [see Fig. 3(a)], the element ∆ k stays zero and states will be ungapped along the path (K − 0.7ŷ) → K. However, the element ∆ k is finite when changing the momentum along the path K → (K− 0.7K) and the spectrum is gapped. In Fig. 4(b), the element ∆ k is finite when approaching the K point from both (K − 0.7ŷ) and (K − 0.7K), and the spectrum is gapped at energies ε k ω/2 + ε D close to the K point.

V. ELECTRONIC OCCUPATION
In this section, we discuss the modification of the occupation on the upper lead due to light irradiation. Direct measurement of the nonequilibrium occupation has been performed on the surface states of a topological insulator by time and angle-resolved photoemission spectroscopy 33 . From the definition of the Green's functions, we can write the energy and k-dependent occupation of the top layer as with the lesser Green's function given by Eq. (25). In a simple picture, the occupation is given by the density of states multiplied by the Fermi functions at different energies f t (ε) and f t (ε ± ω). Figure 5 shows the occupation along the path (K − 0.7ŷ) → K → (K − 0.7K) in the first Brillouin zone and as a function of energy ε for the same parameters as in Fig. 4. We consider the case of linear [ Fig. 5(a)] and circular polarized light [ Fig. 5(b)] with a frequency of ω = 0.6γ 0 . We set the chemical potential to µ t = 0.15γ 0 , the Dirac point energy to ε D = −0.2γ 0 and we consider n = 0, ±1, ±2 sidebands. Without light-matter interaction and at zero temperature, states can be occupied to the chemical potential µ t . A finite coupling strength excites electrons and occupies states above the Fermi level µ t . The highest accessible states at zero temperature are at energies ε = µ t + 2ω corresponding to absorptions of two photons with energy ω by one electron.
It is interesting to note that the gaps in Fig. 5(b) close to energies ε ≈ 0.1γ 0 are different in size. These gaps correspond to the resonant absorption or emission of a single photon with energy ω. The different gap sizes occur due to the matrix elements shown in Fig. 4(b) whose magnitude dependent on the direction from which the Dirac point is approached. The different sizes of the gap at finite energy vanish within the Dirac approximation of the band structure and the optical matrix elements. The effect is a result of the consideration of the full momentum dependence of energies and optical matrix elements in the first Brillouin zone.

VI. DIFFERENTIAL CONDUCTANCE
With the current derived in Eq. (27), we can discuss different scenarios for transport through a graphene contact and identify signatures of the light-matter interaction in transport measurements. We focus on the differ-ential conductance which can be written as with the dc-current given by Eq. (27). Throughout the discussion, we assume that circularly polarized photons are absorbed in the top layer and the photo-excited electrons tunnel to the bottom layer. Motivated by the fact that an underlying substrate 60 or Coulomb interaction 46 opens a gap in the band structure of graphene, we discuss two scenarios: (i) we assume that both the top and the bottom layer consist of graphene. (ii) We discuss the case that the top layer consists of graphene and the bottom layer consists of graphene with a finite mass term opening a gap in the band structure.
A. Tunneling between graphene contacts Figure 6(a) shows the differential conductance for circularly polarized light and different values of the coupling strength α x,y . The frequency of the incoming light is set to ω = 0.6γ 0 and the temperature is T = 0. Since a typical graphene sample is doped we assume a finite Dirac point energy of ε D = −0.2γ 0 relative to the chemical potential µ t .
We first consider the case without light-matter interaction (α x = α y = 0). At small voltages, eV γ 0 , the differential conductance increases with the applied voltage and depends nonlinearly on the voltage for eV → γ 0 . Since we consider the tunnelling current in Eq. (27) this behavior is reflected by the overlap between the density of states of the top and bottom layer. The density of states of single-layer graphene is proportional to ε for ε γ 0 and approaches the van Hove singularity at ε = 0.8γ 0 assuming finite doping ε D = −0.2γ 0 61 .
As discussed in the previous section, a finite coupling strength opens a gap in the quasienergy spectrum of graphene. Even for moderate coupling strength α x = α y = 0.05γ 0 , the opening of a gap in the quasienergy spectrum is manifested in a strong suppression of the differential conductance at certain voltages.
The first suppression of the current close to eV = 0.1γ 0 corresponds to resonant absorption of a single photon between the conduction and valence band. Since the Dirac point energy is at ε = −0.2γ 0 and the frequency of the light is ω = 0.6γ 0 , the suppression of the conductance occurs at eV ≈ ω/2 + ε D = 0.1γ 0 . The second strong suppression of the differential conductance occurs close to eV = 0.4γ 0 and corresponds to resonant absorption of two photons between the conduction and valence band.
At strong coupling α x = α y = 0.1γ 0 , the differential conductance is suppressed over a larger range of voltages compared to the moderate coupling strength α x = α y = 0.05γ 0 since the gap in the quasienergy spectrum increases with increasing coupling strength. In-terestingly, the suppression and the enhancement of the conductance close to voltages eV ≈ 0.1γ 0 gradually decreases and increases at strong coupling α x = α y = 0.1γ 0 . The steps are a signature of the nonlinear band structure of graphene and the momentum dependence of the optical matrix elements for momenta k (K, K ). The gradual increase and decrease can be understood by considering for instance the conductance at voltage eV = 0.15γ 0 in Fig. 6(a). Figure 5(b) shows the occupation corresponding to the voltage eV = 0.15γ 0 and the same parameters as in Fig. 6(a). Approaching the K-point from (K − 0.7ŷ), states are occupied close to the chemical potential µ t . Starting from eV = 0.15γ 0 and increasing the voltage, first states close to k-values (K − 0.2ŷ) contribute to the current and give rise to the first step-like increase of the conductance. Second states close the (K −0.2K) are occupied by increasing the voltage and the conductance gradually increases a second time. A similar argument holds for the suppression of the conductance close to voltages eV = 0.4γ 0 which is related to the resonant absorption of two photons.
In addition to a suppression of the differential conductance, the light-matter interaction induces a peak at eV = 0.2γ 0 [black arrow in Fig. 6(a)]. The process corresponding to the peak is sketch in Fig. 6(b) and is related to absorption of a single photon at energy µ t = 0.2γ 0 . Figure 6(b) shows the density of states summed over all momenta in the first Brillouin zone of the top and bot-tom layer. In the top layer, the states are occupied to the chemical potential µ t . The light-matter interaction excites an electron to energies ε = µ t + ω matching the energy of the van Hove singularity of the top layer. The excited electron then tunnels to the van Hove singularity of the bottom layer providing a large density of unoccupied states. The peak at eV = 0.2γ 0 in the conductance is hence a result of the overlap of the van Hove singularities of both the top and the bottom layer.
B. Tunneling between graphene and graphene with finite mass term So far we discussed the differential conductance for the ideal case that light is absorbed in the top layer. Since the light is only partially absorbed in graphene 8 , also the lower layer will interact with the light. However, because of a substrate or interaction in graphene 46,60 , a band gap will open. In this section, we assume that the bandgap of the bottom layer is larger than the frequency of the light such that the lower layer becomes essentially transparent. In a more general setting, the bottom layer can be made of any two-dimensional material with a bandgap larger than the frequency of the light. In Fig. 7(a) and 8, we consider a finite mass m b in the bottom layer such that the energy of the bare graphene Hamiltonian is given by and the bandgap becomes 2m b . Figure 7(a) shows the differential conductance with the same parameters as in Fig. 6(a) except that now the bottom layer has a mass term m b = 0.4γ 0 and the bandgap is 0.8γ 0 . Without the light α x = α y = 0, the current is completely suppressed at voltages eV < m b because of the bandgap in the bottom layer. For eV > 0.2γ 0 , the voltage is larger than the bandgap on the bottom graphene and electrons can tunnel between the layers. We remark that the Dirac point energy is at ε D = −0.2γ 0 such that the gap edge on the bottom layer is reached for voltages eV = 0.2γ 0 . When the voltage is larger than the gap of the bottom layer, the behaviour is similar to the differential conductance in Fig. 6(a).
As we have seen in previous sections, a finite lightmatter interaction can occupy states above the chemical potential. We expect that the photoexcited electrons contribute to a current although the voltage is smaller than the gap on the bottom lead. This behavior is seen in Fig. 7(a) at finite light-matter coupling strength. The differential conductance for eV < 0.2γ 0 is solely due to the tunnelling of photoexcited electrons. However, since the energies of electrons that are excited by a single photon are smaller than the bandgap edge at eV = 0.2γ 0 of the bottom layer, the current is strongly suppressed for eV < 0.2γ 0 . In this case, the two-photon absorption is the dominant processes to the conductance.
Additionally to the finite conductance at voltages eV < 0.2γ 0 , the differential conductance shows a peak at eV ≈ 0.28γ 0 indicated by the black arrow in Fig. 7(a). ] α x = α y = 0 α x = α y = 0.05γ 0 α x = α y = 0.1γ 0 Figure 8. Differential conductance for circular polarized light and different coupling strength αx and αy. Compared to Fig. 7(a), the frequeny is ω = 0.7γ0. The other parameters are ω = 0.8γ0, n = 0, ±1, ±2, T = 0, η = 0.001γ0. The Dirac point energy is set to εD = −0.2γ0 such that the differential conductance is suppressed close to voltages eV ω/2 + εD = 0.2γ0. The finite conductance in the gap depends on the η = 0.001γ0. The steps in the case of αx = αy = 0.1γ0 are a signature of the non-linear band structure of graphene for k (K, K ). The two black arrows indicate the tunneling processes which are associated with the tunneling of photo-excited electrons to energies of the van Hove singularity of either the top or bottom layer.
The process corresponding to such a peak is sketch in Fig. 7(b) and its origin is similar to the peak at eV = 0.2γ 0 in Fig. 6(a), namely electrons at the chemical potential are photoexcited and can tunnel into a large number of empty states due to the van Hove singularity of the bottom layer. The finite mass term in the bottom layer shifts the van Hove singularity to higher energies such that the peak in Fig. 7(a) appears at eV ≈ 0.28γ 0 . In comparison to Fig. 6(b), the energies of the photoexcited electrons at energy ε ≈ 0.88γ 0 are larger than the van Hove singularity of the top layer but match the energy of the van Hove singularity of the bottom layer.
In Fig. 8 we discuss the differential conductance at higher frequency ω = 0.7γ 0 and the same parameters as in Fig. 7(a). In this case, the conductance shows complex features even below the bandgap eV < 0.2γ 0 . Because of the high frequency ω = 0.7γ 0 and the opening of the gap in the quasienergy spectrum of graphene, the absorption of a single photon is sufficient to excite an electron above the bandgap edge of the bottom layer. The resonant absorption of a single photon gives then rise to the increase of the conductance at voltages eV ≈ 0.08γ 0 .
Additionally, two smaller peaks occur at voltages eV = 0.1γ 0 and eV = 0.18γ 0 indicated by the black arrows in Fig. 8. The origin of both peaks can be associated with the van Hove singularities of both the top and the bottom layer. The van Hove singularity of the bottom layer occurs at energy ε = 0.8γ 0 due to the Dirac point energy ε D = −0.2γ 0 . However, similar to Fig. 7(b), the van Hove singularity of the bottom layer appears at energies ε ≈ 0.88γ 0 because of the finite mass term m b . When an electron at voltage eV = 0.1γ 0 absorbs a photon, it is excited to energies of the van Hove singularity of the top layer and then tunnels to the bottom layer. The opposite process is also possible corresponding to first the excitation of an electron at voltages eV = 0.18γ 0 and second, the tunneling to the energy of the van Hove singularity of the bottom layer.

VII. CONCLUSION
We studied the transport properties of two graphene layers vertically coupled with an insulating layer. A voltage is applied between the top and the bottom layer giving rise to a tunnelling current. In such a contact, we discussed the differential conductance when the top graphene layer is driven to a nonequilibrium state by an external light field.
In agreement with previous result, multiple gaps open in the quasienergy spectrum of graphene depending on the polarization. We obtained that the gaps corresponding to a resonant absorptions/emission of a single photon between the conduction and valence band have different sizes due to the momentum dependence of the optical matrix element. This effect becomes in particular important at frequencies ω 0.5γ 0 where the non-linearity of the bandstructure of graphene can not be neglected.
We studied the differential conductance in two kinds of contacts. First, we considered the differential conductance between two graphene layers. Second, we studied the modification of the conductance when a gap opens in the bottom layer. In both cases, the effects of the light-matter interaction are manifested in a strong suppression of the conductance at voltages eV = nω/2 + ε D due to the opening of gaps in the quasienergy spectrum of graphene. In the case of strong driving, the different gap sizes at energies ε = nω/2 + ε D become apparent in the differential conductance. Additional peaks in the differential conductance can be related to the tunneling of photo-excited electrons to the energies of the van Hove singularity of the top or bottom layer. For the case of a finite bandgap in the bottom layer, the differential conductance is solely due to the photo-excited electrons and is hence a direct signature of the light-matter interaction. Finally, we remark that although we studied the tunneling between two graphene layers, our approach can be extended to tunneling between arbitrary two-dimensional materials.
We conclude with neglected interaction and an outlook about extension to our approach. The light-matter interaction in graphene has been subject of much research and the nonequilibrium driving of graphene has been investigated in various approaches. Most of these approaches focus on the noninteracting picture in which the Coulomb interaction is neglected. However, it has been shown that the electron-electron interaction in two-dimensional ma-terials can be tuned by engineering the surrounding dielectric environment. The electron-electron interaction has a significant effect on the band structure. For instance, Coulomb interaction leads to an open of a gap in the bandstructure and the formation of excitons 46,62 .
The continuous applications of a laser pulse or light field to graphene results in heating and eventually damage of the sample. Under an intense laser pulse, the electron-phonon coupling in graphene will excite phonons and relax the electronic system. The phonon dissipation in the graphene under light irradiation and the effect on the transport properties can be studied by coupling the electronic system to a bath of phonons similar to the approaches in Refs. 63 and 64.
The Green's functions satisfy the following relations which are useful to derive the retarded and advanced Green's functions. The Green's functions are related by