# Ballistic transport through irradiated graphene

Jonathan Atteia, Jens H Bardarson, Jérôme Cayssol

► To cite this version:

Jonathan Atteia, Jens H Bardarson, Jérôme Cayssol. Ballistic transport through irradiated graphene. Physical Review B: Condensed Matter and Materials Physics (1998-2015), 2017, 96 (24), pp.245404. 10.1103/PhysRevB.96.245404. hal-01679895

HAL Id: hal-01679895

https://hal.science/hal-01679895

Submitted on 10 Jan 2018

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Ballistic transport through irradiated graphene

Jonathan Atteia,¹,* Jens H. Bardarson,²,³,† and Jérôme Cayssol¹,‡

¹LOMA (UMR-5798), CNRS and Bordeaux University, F-33045 Talence, France

²Department of Physics, KTH Royal Institute of Technology, Stockholm, SE-106 91 Sweden

³Max-Planck-Institut für Physik komplexer Systeme, Nöthnitzer Str. 38, 01187 Dresden, Germany

The coherent charge transport through an illuminated graphene ribbon is studied as a function of electronic doping and characteristics of the electromagnetic driving, for monochromatic circularly polarized light. We focus on the DC current carried by 2D bulk carriers which is dominant (over edge transport) for short and wide enough samples. We investigate how the ballistic conductance suppression, due to photon resonances between the valence and conduction bands, evolves when the experimentally tunable parameters are varied. The residual conductance can be associated with evanescent states and related to dynamical gaps in the Floquet quasienergy spectrum.

## I. INTRODUCTION

Topological insulators [1,2] are a novel class of materials which have attracted a lot of attention in the last few years. These materials possess a bulk band gap, and due to a specific topological order of the bands, robust metallic states appear at their boundary. The topologically protected edge or surface states arise from the particular crystalline structure of the material, which is unfortunately not easy to engineer or manipulate. Recently, it has been predicted that topologically trivial materials can be turned into topological insulators by driving them with a time-periodic external perturbation [3–7]. The formalism used to study those periodically driven systems is based on the Floquet theorem [8,9], which is the analog of the Bloch theorem for periodic potentials in time. These so-called Floquet topological insulators [5,10] (FTI’s) have been classified using new types of topological invariants [6,7,11].

An important candidate for the realization of a condensed matter FTI is irradiated graphene. Two distinctive situations may be encountered depending on whether the photon energy $\hbar\omega$ is larger or smaller than the bandwidth of the crystal. In both cases, a gap opens at quasienergy $\epsilon = 0$ realizing a Haldane topological insulator with photoinduced chiral edge states [3,5,6,12]. In the latter case, an additional dynamical gap at $\epsilon = \hbar\omega/2$ may also open due to one-photon resonance between the two bands of graphene [13–15]. Photoinduced chiral edge states appear in the gaps, thereby linking the distinct Floquet bands [16–19]. The topological nature of these edge states and their connection to transport properties (quantized Hall conductance) has been extensively studied [12,18,20–22] mostly in the absence of dissipation. Dissipation, occurring by coupling to phonons, leads to nontrivial occupation probabilities which usually gives rise to a nonquantized Hall conductance [23–27]. The Floquet-Bloch states have been observed by irradiating the surface states of a 3D topological insulator with circularly polarized midinfrared light and measuring the Floquet spectrum through

time resolved ARPES pump-probe experiments [28]. Note that FTIs have also been realized in cold atom systems [29] and in photonic crystals [30]. Besides, irradiated graphene has also been studied extensively in connection with pumping effects [31–33] and irradiated pn junctions for linear polarization [34–37]. More recently driven conductance of graphene and Weyl semimetal thin films under linear polarized THz radiation has been investigated in the regime of strong driving and low doping [38].

Our main motivation is to suggest the use of DC transport measurements to investigate the features of the Floquet quasienergy spectrum in coherent driven electronic systems. We therefore compute the conductance of graphene transistors driven by an electromagnetic wave, typically in the terahertz (THz) or infrared (IR) range. Indeed photon energies in the range 10–100 meV can be matched by chemical potential variations in typical graphene samples. Recently, photon-assisted shot noise has been studied experimentally in coherent diffusive graphene samples irradiated in the THz range [39].

In this paper, we consider a graphene based field effect transistor (gFET), whose conduction channel is irradiated by a circularly polarized electromagnetic wave at normal incidence (Fig. 1). The carrier density (in the absence of irradiation) can be tuned using a DC electrostatic backgate. The source and drain leads are heavily doped and are not irradiated. The two-terminal differential conductance is studied as a function of the irradiation strength, chemical potential $\mu$ in the central region, and photon energy $\hbar\omega$ for various lengths $L$ and widths $W$ of the graphene ribbon (Fig. 1). We consider the ballistic regime which has been reached experimentally in high-mobility suspended or encapsulated graphene samples [40–42]. The photon energies are in the range 10–100 meV, therefore much smaller than the electronic bandwidth of graphene. In this regime, during their time of flight across the transistor, the electrons are dressed coherently by the photons. The Floquet-Landauer-Büttiker scattering theory allows us to evaluate the source-drain conductance of the gFET. This approach implicitly assumes that the dissipation takes place only in the leads and not in the central irradiated region. This assumption puts a restriction in the length $L$ of our sample. Due to the low electron-phonon scattering [43] in graphene, especially below the optical phonon energies (200 meV), we expect that dissipation can be neglected in submicrons samples

<footnote>*jonathan.atteia@u-bordeaux.fr</footnote>

<footnote>†bardarson@kth.se</footnote>

<footnote>‡jerome.cayssol@u-bordeaux.fr</footnote>
<figure><img src="image_2.png" /><figcaption>FIG. 1. Field effect transistor consisting in a rectangular graphene ribbon contacted between source (S) and drain (D) leads. A remote electrostatic backgate (not represented) allows us to tune the chemical potential (and carrier density) in the central graphene region. A circularly polarized electromagnetic wave is shone at normal incidence, corresponding to an in-plane rotation vector potential $A(t) = A_0(\cos \omega t, \sin \omega t)$.</figcaption></figure>

at low enough temperatures. Experimentally, the conductance in a coherent ballistic periodically driven conductor has not been investigated yet for all $\mu$, and most theoretical studies focused so far on undoped graphene ($\mu = 0$) [3,5,12,20] and edge state transport. For doped graphene samples, the conduction by bulk states is always important and often exceeds the edge state contribution. For this reason, this paper mainly focuses on the bulk conduction through 2D evanescent and propagating electronic states. Finally the response of a graphene-based transistor to electromagnetic driving remains a crucial issue for optoelectronic applications [44].

The main results are the following. First, the conductance can be strongly suppressed in broad domains of chemical potential (or equivalently carrier density). The qualitative picture is that time-dependent driving opens gaps in the quasienergy spectrum, inducing evanescent states characterized by finite penetration lengths. Hence, the conductance-chemical potential curves depend drastically on the length of the central conducting channel. Such suppression of the graphene DC conductance by irradiation was first predicted by Calvo et al. [13–15]. In this paper, we study how these features evolve when experimentally tunable parameters are varied. Second, for undoped graphene, the conductance exhibits nonmonotonic variations as a function of the irradiation strength. At strong driving and high frequency, the conductance at the Dirac point has been interpreted as resulting from an infinite set of evanescent states [12]. Here we provide the full evolution of the conductance as a function of the irradiation strength and sample length, while showing that only a finite number of evanescent states contribute for moderate driving.

The paper is organized as follows. In Sec. II, the Floquet quasienergy spectrum of irradiated graphene is reviewed [3,12,16,17,45]. In Sec. III, the scattering problem is formulated and the expression for the differential conductance is derived within the Floquet-Landauer-Büttiker formalism. Section IV is devoted to the dependence of the conductance as a function of $\mu$ and irradiation, while we focus on the case of undoped graphene ($\mu = 0$) in Sec. V. Suitable parameters for experimental realizations are discussed in Sec. VI.

## II. FLOQUET SPECTRUM IN GRAPHENE

In this section, we review the Floquet quasienergy spectrum of an infinite graphene sheet irradiated by a circularly polarized electromagnetic wave [3,12,16,17,45]. The purpose of this section is to set the stage for the transport calculations (Secs. III, IV, and V) and make the paper self-contained.

### A. Geometry of the system and Hamiltonian

The graphene sheet is irradiated by a circularly polarized radiation of frequency $\omega = 2\pi/T$. Photon energies $\hbar\omega$ are much smaller than the bandwidth of graphene, which allows us to use the low-energy model for graphene. The massless Dirac fermions come in two flavors located around the two corners (valleys) of the Brillouin zone. Moreover, we consider ballistic graphene and only vertical (namely conserving momentum) electronic transitions (absorption or emission of photons) are allowed. Therefore intervalley scattering can be neglected, and the decoupled ($K$ and $K'$) valleys can be studied separately.

The time-dependent Hamiltonian describing the Dirac fermions in valley $\xi$ is written as:

$$H^{\xi}(t) = H_{0}^{\xi} + V^{\xi}(t), \qquad (1)$$

$\xi = \pm 1$ being the valley index. The time-independent Hamiltonian $H_0^\xi$, corresponding to nonirradiated graphene, reads [46,47]:

$$H_0^{\xi} = v(\xi\sigma_x p_x + \sigma_y p_y) - \mu, \qquad (2)$$

where $\sigma_x$ and $\sigma_y$ are the Pauli matrices associated with the sublattice isospin, $(p_x, p_y) = -i\hbar(\partial_x, \partial_y)$ are the components of the wave-vector operator $\mathbf{p}$, $v$ is the Fermi velocity, and $\mu$ is the chemical potential.

The vector potential $A(t) = A_0(\cos \omega t, \sin \omega t)$ couples to the electric charge via the Peierls substitution $\mathbf{p} \rightarrow \mathbf{p} - q\mathbf{A}(t)$, where $q = -e$ is the electron charge. The Hamiltonian for the coupling between electrons and the electromagnetic field can therefore be expressed as:

$$V^{\xi}(t) = (V_{1}^{\xi} e^{i\omega t} + V_{-1}^{\xi} e^{-i\omega t}), \qquad (3)$$

where the matrices $V_1^\xi$ and $V_{-1}^\xi$ read:

$$V_{1}^{\xi} = \frac{evA_{0}}{2}(\xi\sigma_{x} - i\sigma_{y}), \quad V_{-1}^{\xi} = \frac{evA_{0}}{2}(\xi\sigma_{x} + i\sigma_{y}). \qquad (4)$$

An important dimensionless parameter is

$$\beta = \frac{evA_0}{\hbar\omega}, \qquad (5)$$

which quantifies the electromagnetic driving strength. Following a simple quasiclassical argument, $evA_0$ is simply the energy gained by an electron traveling at speed $v$ in an electric field $E_0 = \omega A_0$ during a period of the electromagnetic wave $(1/\omega)$, while $\hbar\omega$ is the minimal energy quantum which can be absorbed by the electron. The typical electric field can be expressed as:

$$E_0 = \beta \frac{\hbar\omega}{el_{\omega}}, \qquad (6)$$

where $l_\omega = v/\omega$ is the distance traveled by the electron during a period of the electromagnetic field (divided by $2\pi$).
## B. Floquet spectrum

In this section, we assume translational invariance along the x and y axes and seek for solutions of the form $Ψ(t)e^{ik_x x+ik_y y}$. According to the Floquet theorem, the time-dependent part of the wave functions can be written as:

$$ \Psi(t) = \Phi_{\varepsilon}(t)e^{-i\varepsilon t/\hbar}, \qquad (7) $$

where $Φ_ε(t)$ is a T-periodic function and $ε$ is the corresponding quasienergy. The quasienergy $ε$ is defined modulo $\hbar\omega$, thereby making it possible to redefine $ε$ such that it always belongs to the interval $[-\hbar\omega/2, \hbar\omega/2]$.

The wave function obeys the Schrödinger equation

$$ i\hbar\partial_t\Psi(t) = H^\xi(t)\Psi(t). \qquad (8) $$

After injecting Eq. (7) in the Schrödinger equation, the T-periodic function $Φ_ε(t)$ appears to be an eigenstate of the Floquet Hamiltonian $H_F(t) = H(t) - i\partial_t$ with eigenvalue $ε$, namely:

$$ (H^\xi(t) - i\hbar\partial_t)\Phi_\varepsilon(t) = \varepsilon\Phi_\varepsilon(t). \qquad (9) $$

Using the following Fourier expansion for $Φ(t)$:

$$ \Phi_{\varepsilon}(t) = \sum_{m \in \mathbb{Z}} \Phi_{m} e^{im\omega t}, \qquad (10) $$

allows us to transform the time-dependent differential equation Eq. (9) into a time-independent eigenstate problem [8,9]:

$$ \sum_{n} H_{F,mn}^{\xi} \Phi_{n} = \sum_{n} (H_{0F,mn}^{\xi} + V_{F,mn}^{\xi}) \Phi_{n} = \varepsilon \Phi_{m}, \qquad (11) $$

where the quasienergy index $ε$ is omitted in all Fourier components $Φ_n$, and $(m,n) \in \mathbb{Z}^2$. The Floquet Hamiltonian $H_{0F,mn}^\xi$, defined as

$$ H_{0F,mn}^{\xi} = (\hbar v(\xi \sigma_x k_x + k_y \sigma_y) - \mu + m\hbar\omega)\delta_{mn}, \qquad (12) $$

is diagonal in the Floquet basis. In contrast, the driving Hamiltonian $V_{F,mn}^\xi$,

$$ V_{F,mn}^{\xi} = V_{1}^{\xi} \delta_{m,n-1} + V_{-1}^{\xi} \delta_{m,n+1}, \qquad (13) $$

couples distinct Fourier components $Φ_n$. Matrices $V_1^\xi$ and $V_{-1}^\xi$ are defined by Eq. (4).

In this Floquet-Fourier representation, the wave function is an infinite vector containing all the different Fourier harmonics of the wave function. Since the Floquet matrix is infinite, we need to set a cutoff $N$ in the number of Floquet replicas considered. The size of the matrix is thus $2(2N+1) \times 2(2N+1)$ (taking into account the sublattice isospin index). Afterwards, it will be necessary to check the stability of the results by increasing $N$.

Once the cutoff is made, the Floquet matrix $H_F^\xi$ can be easily diagonalized numerically to obtain the dispersion relation. Due to rotational invariance around $k=0$, the spectrum only depends on the norm $k=|k|$ of the wave vector. The spectrum consists of $2(2N+1)$ bands having the dispersion relation $ε_α(k)$ for $α \in [1,2(2N+1)]$. The Floquet bands centered around $ε=0$ are plotted as functions of $k$ in Fig. 2 for various driving strengths, using seven Fourier components ($N=3$).

<figure><img src="image_3.png" /></figure>

<figure><img src="image_4.png" /></figure>

<figure><img src="image_5.png" /><figcaption>FIG. 2. Dispersion relation $ε_α(k)$ ($k=|k|$) of an infinite graphene sheet irradiated by an electromagnetic wave at various driving strengths: (a) $β=0.1$, (b) $β=0.3$, and (c) $β=0.6$ using $N=3$ Floquet replicas. The dashed lines correspond to the nonirradiated case in the Floquet representation with dispersion $ε = \pm\hbar vk + n\hbar\omega$. The circles on (a) and (b) indicate the nonequivalent gaps $Δ_m$, located at momenta $k=m\omega/2v$, and at quasienergies $ε=0$ and $ε=\hbar\omega/2$, alternatively.</figcaption></figure>

For weak driving, $β=0.1$ [Fig. 2(a)], the bands corresponding to the irradiated case (solid curves) follow closely the nonirradiated ones (dashed curves) except around avoided crossings $Δ_1$ located at $k = \pm\hbar\omega/2v$ and centered around $ε = \pm\hbar\omega/2$. Upon increasing the driving, more gaps open at momenta $k_m = \pm m\hbar\omega/2v$, with $m$ integer [Fig. 2(b)]. These gaps correspond to the anticrossing of two Dirac bands dressed with $n$ and $n'$ photons such that $|n-n'| = m$. For example, at the gap $Δ_3$ in Fig. 2(b), the electron is in a coherent superposition of a valence band electron dressed with 2 photons (having absorbed two quanta of the electromagnetic field) and a conduction band electron dressed with $-1$ photon
(meaning having emitted one quantum of the field). For strong driving, multiphoton processes are more likely, leading to a strong modification of the quasienergy spectrum [Fig. 2(c)].

In the next sections, we study with more detail the gap structure as a function of $\beta$. The spectrum is periodic in energy (with period $\omega$), so we restrict ourselves to the study of the reduced bands where $\epsilon \in [-\hbar\omega/2, \hbar\omega/2]$. The gaps are located at different momenta because of the folding of the dispersion relation in the reduced Floquet zone $[-\hbar\omega/2, \hbar\omega/2]$. In this representation, there are two nonequivalent set of gaps centered at $\epsilon = 0$ and at $\epsilon = \hbar\omega/2$ [Fig. 2(b)]. The gaps of even order ($\Delta_0, \Delta_2, \Delta_4...$) are located in the center of the Floquet zone ($\epsilon = 0$, Sec. IIC) and the odd order ones ($\Delta_1, \Delta_3...$) are located at the edge ($\epsilon = \pm \hbar\omega/2$) [Fig. 2(b)].

### C. Photoinduced gaps vs $\beta$

The Floquet spectrum has a very rich gap structure [16] that depends strongly on the driving strength $\beta$. Understanding this structure is crucial in order to study the transport properties of irradiated graphene (see Secs. IV and V). The gap $\Delta_0$ and the gaps $\Delta_m$ with $m \ge 1$ have different origins. Multiphoton processes of order $m$ induce the gaps $\Delta_m$. For $m \ge 1$, the gap amplitudes have been derived within the generalized rotating wave approximation (RWA) in Ref. [45]:

$$ \frac{\Delta_m}{\hbar\omega} = \beta|J_{m+1}(2\beta) - J_{m-1}(2\beta)| \approx \frac{\beta^m}{(m-1)!}, \quad (14) $$

where $J_m$ is the $m$th order Bessel function, and the last approximation is only valid for $\beta \ll 1$. For weak driving, the size of the gap $\Delta_m$ ($m \ge 1$) is strongly reduced upon increasing $m$. Increasing the driving strength $\beta$, the gaps evolve in a nonmonotonic way and can also close. In this section, we will compare analytical expressions of the gap size obtained within the rotating wave approximation scheme and the results that we get from diagonalizing in Floquet space (Fig. 3). The RWA fits very well the numerics for low $\beta$, and describes qualitatively the size of the gaps for $\beta \lesssim 1$ (Fig. 3).

At energy $\epsilon = 0$, Oka and Aoki [3] first discussed the central gap $\Delta_0$ at $k = 0$ and found that:

$$ \frac{\Delta_0}{\hbar\omega} = \sqrt{1 + 4\beta^2} - 1 \approx 2\beta^2, \quad (15) $$

where the last approximation correspond to weak driving $\beta \ll 1$. The quadratic dependence in $\beta$ for weak driving points out that this gap originates from a second order process in the coupling: the emission and absorption of a virtual photon. The analytic expression Eq. (15) fits perfectly our numerics. For weak driving, the evolution is quadratic in $\beta$, and as the driving increases, it becomes linear until the gap size reaches the size of the Floquet zone $\epsilon = \hbar\omega$ at $\beta = \sqrt{3}/2$. At this particular value of $\beta$, the gap $\Delta_1$ vanishes as shown in Fig. 3(b). This is a case of a band touching induced by strong electromagnetic driving.

For weak driving, the gap $\Delta_1$ is located at wave vector $k = \pm \omega/2v$ and around energy $\epsilon = \hbar\omega/2$. This gap originates from the anticrossing of the bands $n = 1$ and $n' = 0$, and therefore it is a first order process in $\beta$ (exchange of one photon). This explains why $\Delta_1$ is the largest gap for weak driving [Figs. 2(a) and 2(b)]. The size of this gap first increases linearly with $\beta$,

<figure><img src="image_6.png" /></figure>

<figure><img src="image_7.png" /><figcaption>FIG. 3. Size of the gaps at (a) $\mu = 0$ and (b) $\mu = \hbar\omega/2$ as a function of the driving strength $\beta$ calculated using the numerical Floquet method (solid lines) and the rotating wave approximation (RWA), namely Eqs. (15) and (14) (dashed curves). The numerical calculations (solid lines) agrees with RWA estimations for weak driving. The gap $\Delta_1$ closes while the gap $\Delta_0$ tends to $\hbar\omega$ for $\beta = \sqrt{3}/2$.</figcaption></figure>

then it reaches a maximum value around $\beta \approx 0.5$, and starts decreasing until it reaches zero [Fig. 3(b)]. We see that at the same value of $\beta$, the gap $\Delta_0$ reaches $\hbar\omega$ (Fig. 3), and as we have seen earlier, this happens for $\beta = \sqrt{3}/2 \approx 0.87$. We conclude that this value corresponds to a band touching.

The gap $\Delta_2$ originates from the anticrossing of two Dirac bands dressed with $n = 1$ and $n' = -1$ photons. This is therefore a second order process in the coupling and is proportional to $\beta^2$ for weak driving. The idea is the same for the gap $\Delta_4$; it originates from the anticrossing of bands dressed with 2 and $-2$ photons, and is thus proportional to $\beta^4$ for weak driving.

### III. FORMALISM

In this section, we present the formalism used to compute the two-terminal conductance of a graphene-based transistor
whose central 2D conducting channel is irradiated by circularly polarized light at normal incidence. The transmission and reflection coefficients are evaluated within the Landauer-Büttiker formalism extended to driven Floquet systems [48]. Note that the electromagnetic field is uniform in the central gated graphene region, while the leads are not irradiated.

## A. Geometry and scattering problem

We consider an irradiated rectangular ribbon of width $W$ and length $L$ connected to two nonirradiated leads (Fig. 1). The Hamiltonian being time dependent, the single electron energy is not conserved. However, due to the invariance under translations of time $T$, the quasienergy $\epsilon$ is still a good quantum number. In the scattering problem, we consider electrons incoming at energy $\epsilon = 0$, which corresponds to the Fermi level of the leads. The electrons propagating inside the irradiated region can emit or absorb one or several photon(s) of energy $\hbar\omega$. These inelastic scattering events are described by a set of transmission (reflection) coefficients $t_{n0}$ ($r_{n0}$) corresponding to the amplitude of being transmitted (reflected) with a final energy $\epsilon + n\hbar\omega$.

In order to determine those Floquet scattering coefficients, periodic boundary conditions are used along the $y$ direction. Hence the transverse wave vector $k_y$ is quantized as:

$$k_y = \frac{2\pi n_y}{W}, \qquad (16)$$

with $n_y$ a relative integer. Note that with such periodic boundary conditions, only bulk 2D states are investigated. In contrast to $k_y$ and $\mu$, the longitudinal wave vector $k_x$ is not a good quantum number because the system is not translationally invariant along the $x$ axis.

In the following, all the two-component wave functions are written as $\Phi(x,y,t) = \Phi(x,t)e^{ik_y y}$, and we will work with the $\Phi(x,t)$ wave functions ($k_y$ being omitted) corresponding to the effective 1D transport problem at a given $k_y$. In order to calculate the two-terminal conductance, we need to match the wave functions at the interfaces $x = 0$ and $x = L$.

## B. Irradiated ribbon

In the irradiated region ($0 < x < L$), we need to find the Floquet eigenstates $\Phi_m$ and the longitudinal wave vector $k_x$ corresponding to a given set ($\mu, k_y$) (still having in mind the scattering problem for $\epsilon = 0$ corresponding to an electron incident from the lead Fermi level). Multiplying Eq. (11) by $\xi\sigma_x$ and rearranging it allows us to get both $k_x(\epsilon)$ and $\Phi_m$ by solving the following equations:

$$\sum_n (K_{0F,mn}^{\xi} + K_{VF,mn}^{\xi})\Phi_n = \hbar vk_x \Phi_m, \qquad (17)$$

where one has defined an infinite matrix $K_{0F}$ with matrix elements:

$$K_{0F,mn}^{\xi} = \xi((\mu - m\hbar\omega)\sigma_x - i\sigma_z\hbar vk_y)\delta_{mn}, \qquad (18)$$

and the infinite matrix $K_{VF}$ by the matrix elements:

$$K_{VF,mn}^{\xi} = -\xi\sigma_x V_1^{\xi} \delta_{m,n-1} - \xi\sigma_x V_{-1}^{\xi} \delta_{m,n+1} \qquad (19)$$

describing the effect of the electromagnetic field.

<figure><img src="image_8.png" /><figcaption>FIG. 4. Dispersion relation $k_x^\alpha(\epsilon)$ for $\beta = 0.3$ in the reduced Floquet zone with the real and imaginary parts of the wave vector $k_x$ for $k_y = 0$ and $N = 2$.</figcaption></figure>

Note that changing the valley index reverses the sign of $k_x$, and because the ribbon is space-inversion symmetric along the $x$ direction, the eigenvalues $k_x$ always come in pairs, therefore both valleys have the same contribution to the conductance. In practice, one can solve the scattering problem and compute the transport in a given valley, and multiply the single-valley result by 2 to get the total conductance.

Besides, $K_F = K_{0F} + K_{VF}$ is not Hermitian so $k_x$ can have an imaginary part, which corresponds to an evanescent state (Fig. 4). As we have seen before, we need to make a truncation in the number of Floquet side bands such that $m \in [-N, N]$. After numerical diagonalization, we obtain $2(2N + 1)$ eigenvalues $k_x^\alpha(\epsilon)$ and eigenvectors $\Phi_m^\alpha$ with $\alpha \in \{1, 2(2N + 1)\}$ labeling the different eigenmodes. The resulting two-component wave function is written as a superposition of the $2(2N + 1)$ eigenmodes labeled by index $\alpha$. The wave function of a given $\alpha$ eigenmode is expressed as:

$$\Phi^{\alpha}(x,t) = e^{ik_x^{\alpha}x} \sum_{m=-N}^{N} \Phi_m^{\alpha} e^{-im\omega t}, \qquad (20)$$

which is still a two-component spinor (attached to one valley, here $\xi = 1$).

When the Fermi level lies in the gap, $k_x^\alpha(\epsilon)$ has a nonzero imaginary part which signals an evanescent state. For a ribbon with a finite length, these states will allow the current to tunnel through the gapped central region. The penetration lengths $\xi$ are defined as the inverse of the imaginary part of the longitudinal wave vector $k_x(\epsilon)$ inside the gap (Fig. 4):

$$\xi_{\alpha} = \frac{1}{|\Im\{k_x^{\alpha}\}|}. \qquad (21)$$

For weak driving, $\beta \ll 1$, it is possible to associate each eigenmode ($\alpha$) with a given anticrossing between uncoupled Floquet replicas. Moreover, there is a direct relation between
the length $\xi$ and the size of the gaps:

$$ \xi_m = \frac{\hbar v}{\Delta_m}. \qquad (22) $$

As the size of the gap increases, the imaginary part of $k_x(\varepsilon)$ (where $\varepsilon$ is in the gap) increases, thereby the characteristic length decreases. According to Eqs. (15) and (14), and for $\beta \ll 1$, the decay lengths are expressed as:

$$ \xi_0 = \frac{l_\omega}{\sqrt{1 + 4\beta^2 - 1}} \approx \frac{l_\omega}{2\beta^2}, \qquad (23) $$

$$ \xi_m \approx \frac{(m-1)!}{\beta^m} l_\omega, \qquad (24) $$

in terms of the length:

$$ l_\omega = v/\omega, \qquad (25) $$

corresponding to the distance traveled by an electron during one cycle of driving (divided by $2\pi$). We compared the RWA method and the numerics and we see a good agreement between numerical diagonalization and RWA even for $\beta \approx 1$.

### C. Leads

To simplify the expression of the scattering states in the leads, we consider heavily n-doped leads, in the same manner as in Ref. [49] for the nonirradiated case. The chemical potential in the leads is $\mu_\infty > 0$ for $x < 0$ and $x > L$ such that $\mu_\infty \gg |\mu|$.

In order to account for electrons leaving the scattering region with energy $\varepsilon + n\hbar\omega$, we use the Floquet theorem in the leads. However, because there is no driving, the Floquet replicas are decoupled. In the leads, for an infinite system, Eq. (9) with $\varepsilon = 0$ becomes:

$$ \hbar v(\sigma_x k_x + \sigma_y k_y)\Phi_m(x) = (\mu_\infty - m\hbar\omega)\Phi_m(x). \qquad (26) $$

For a given 2D wave vector $(k_x, k_y)$, there is an infinite amount of plane wave solutions labeled by their Floquet index $m$ and their band index $s = \pm 1$:

$$ \Phi_{m,s}(x) = \frac{1}{\sqrt{2}} \begin{pmatrix} s \\ e^{i\varphi} \end{pmatrix} e^{ik_x x} \qquad (27) $$

with $\mu_\infty = m\hbar\omega + s\hbar v\sqrt{k_x^2 + k_y^2}$ and $\cos\varphi = k_x/\mu_\infty$, $\varphi$ being the angle of incidence of the electron. Since $\mu_\infty \to \infty$, we have $\varphi \to 0, \pi$ and $s = +1$ in both source and drain leads. Finally, the spinor in the leads is independent of $m$ and we obtain one left (+) and one right (-) going solution:

$$ \Phi_m^{\pm}(x) = \Phi^{\pm} e^{ik_x^{\pm}x}, \qquad (28) $$

with $k_x^{\pm} = \pm\mu_\infty$, and

$$ \Phi^{\pm} = \frac{1}{\sqrt{2}} \begin{pmatrix} 1 \\ \pm 1 \end{pmatrix}. \qquad (29) $$

### D. Matching at the interfaces

Using the expression of the eigenstates in the different regions, it is possible to construct the scattering states. We consider an incoming wave from the left lead at the quasienergy $\varepsilon = 0$ and Floquet index $n = 0$. In the left lead, there will be

$2N + 1$ waves reflected at energies $\varepsilon = n\hbar\omega$ with amplitude $r_{n0}$. The wave function in the left lead is then:

$$ \Psi^L(x,t) = \left( \Phi^+ e^{ik_x^+ x} + \sum_{n=-N}^{N} r_{n0} \Phi^- e^{ik_x^- x} e^{in\omega t} \right) e^{-i\varepsilon t}. \qquad (30) $$

In the irradiated region, the wave function is a superposition of the eigenstates $\Phi^\alpha$ with amplitudes $a_\alpha$:

$$ \begin{aligned} \Psi^L(x,t) &= \sum_{\alpha} a_{\alpha} \Phi^{\alpha}(x,t) e^{-i\varepsilon t} \\ &= \left( \sum_{\alpha} a_{\alpha} e^{ik_x^{\alpha} x} \sum_{n=-N}^{N} \Phi_n^{\alpha} e^{-in\omega t} \right) e^{-i\varepsilon t}. \end{aligned} \qquad (31) $$

In the right lead, only the right going states with amplitude $t_{n0}$ on the n's replica are chosen:

$$ \Psi^R(x,t) = \left( \sum_{n=-N}^{N} t_{n0} \Phi_n^+ e^{ik_x^+(x-L)} e^{in\omega t} \right) e^{-i\varepsilon t}. \qquad (32) $$

To calculate the reflected and transmitted amplitudes, we need to match the wave functions at the interfaces. The boundary condition at $x = 0$ is $\Psi^L(x = 0^-, t) = \Psi^I(x = 0^+, t)$:

$$ \Phi^+ + \sum_{n=-N}^{N} r_{n0} \Phi_n^- e^{-in\omega t} = \sum_{\alpha} \sum_{n=-N}^{N} a_{\alpha} \Phi_n^{\alpha} e^{-in\omega t}, \qquad (33) $$

and at $x = L$, the condition is $\Psi^I(x = L^-, y, t) = \Psi^R(x = L^+, t)$ so:

$$ \sum_{\alpha} \sum_{n=-N}^{N} a_{\alpha} \Phi_n^{\alpha} e^{ik_x^{\alpha} L} e^{-in\omega t} = \sum_{n=-N}^{N} t_{n0} \Phi_n^{+} e^{-in\omega t}. \qquad (34) $$

These boundary conditions, at $x = 0$ and $x = L$, must be valid at any time $t$, so we can project them on the different Fourier harmonics:

$$ \Phi^{\pm}\delta_{n0} + r_{n0}\Phi^{-} = \sum_{\alpha} a_{\alpha}\Phi_{n}^{\alpha}, \qquad (35) $$

$$ \sum_{\alpha} a_{\alpha} \Phi_{n}^{\alpha} e^{ik_{x}^{\alpha} L} = t_{n0} \Phi^{+}. \qquad (36) $$

Finally, on one hand, the number of unknown scattering parameters is $4(2N + 1)$ since we have $(2N + 1)$ reflection coefficients $r_{n0}$, $(2N + 1)$ transmission coefficients $t_{n0}$, and $2(2N + 1)$ coefficients $a_\alpha$ to determine. On the other hand, each matching condition represents $2N + 1$ spinor relations, hence there is a total of $4(2N + 1)$ linear relations between those coefficients.

The number of Floquet replicas $N$ defines the cutoff used to calculate the eigenstates and eigenvectors of the system and therefore to obtain the transmission coefficients. We have systematically checked the convergence of our results with respect to $N$.

### E. Conductance formula

To calculate the conductance of the sample, we use the scattering theory extended to Floquet systems [48]. The scattering matrix element $t_{n0,k_y,\mu}(\varepsilon)$ is the probability amplitude for an electron entering from the left lead at energy $\varepsilon$ and wave
vector $k_y$ to exit in the right lead with energy $\epsilon + n\hbar\omega$. The expression of the current through the sample is [48,50]:

$$I(\mu) = \frac{e}{h} \int_{-\infty}^{\infty} d\epsilon \sum_{n=-N}^{N} (T_n(\epsilon)f_L(\epsilon) - T'_n(\epsilon)f_R(\epsilon)), \quad (37)$$

where $f_L(\epsilon)$ and $f_R(\epsilon)$ are the Fermi-Dirac distributions in the left and right leads, respectively. The transmission $T_n(\epsilon)$ is the sum over all $k_y$ channels of the transmission probabilities of electrons from energy $\epsilon$ in the left lead towards energy $\epsilon + n\hbar\omega$ in the right lead:

$$T_n = \sum_{k_y} |t_{n0,k_y}|^2. \quad (38)$$

The current flowing in the scatterer is conserved so the scattering matrix is unitary. This implies the relation:

$$\sum_{n=-N}^{n=N} (T_n + R_n) = 1, \quad (39)$$

where $R_n$ is defined as the sum of the reflection probabilities $|r_{n0,k_y}|^2$ of the transverse modes $k_y$:

$$R_n = \sum_{k_y} |r_{n0,k_y}|^2. \quad (40)$$

We checked that the unitary relation Eq. (39) is fulfilled in the numerical implementation. We can also notice that $T_n(\epsilon) = T'_n(\epsilon)$, $T'_n(\epsilon)$ being the transmission probability defined similarly as $T_n(\epsilon)$, but for the reversed scattering process (incident electron going from the right lead and transmitted to the left lead). The equality $T_n(\epsilon) = T'_n(\epsilon)$ ensures that there is no pumped current for zero bias. Therefore, at zero temperature, for small bias so that the transmission coefficient varies weakly, we obtain the conductance formula:

$$G(\mu) = \frac{\partial I}{\partial V} = G_0 \sum_{n=-N}^{N} T_n(\mu), \quad (41)$$

where $G_0 = \frac{4e^2}{h}$, where the factor 4 accounts for valley and spin degeneracy.

## IV. CONDUCTANCE-CHEMICAL POTENTIAL CURVES

The main goal of this paper is to evaluate the two-terminal conductance $G$ of graphene-based transistors as a function of the chemical potential $\mu$, radiation strength $\beta$, frequency $\omega$, and geometrical parameters $L$ and $W$ (Fig. 1). We consider the ballistic regime relevant for currently achievable high-mobility samples [40–42]. In this section, the dependence of the conductance upon the chemical potential is discussed for various irradiation strengths $\beta$. After a short review of the nonirradiated case ($\beta = 0$) [49], the results are presented by increasing the driving parameter $\beta$ from low to strong driving. The main features consist of strong suppressions of the conductance in wide ranges of chemical potential, especially around $\mu = \pm \hbar\omega/2$. The suppression of the conductance originates from photoinduced gaps in the quasienergy spectrum which lead to evanescent states. A simple phenomenological model is introduced which accounts for the residual conductance around $\mu = \pm \hbar\omega/2$. A crossover between 2D transport through

bulk evanescent states and 1D edge transport is predicted and shown to depend on the geometrical parameters $L$ and $W$ of the graphene ribbon. Specific features at $\mu = 0$ (undoped graphene) are discussed more thoroughly in the next section V.

### A. Nonirradiated ribbon $\beta = 0$

The conductance of a graphene ribbon contacted by two heavily doped leads as a function of the chemical potential has been studied by Tworzydło et al. [49]. The system acts as an electronic Fabry-Perot interferometer whose transmission coefficient reads:

$$T_{k_y}^{(p)} = \frac{k_x^2}{k_x^2 \cos^2(k_x L) + \left(\frac{\mu}{\hbar v}\right)^2 \sin^2(k_x L)}, \quad (42)$$

for propagating modes with transverse wave vector $k_y = 2\pi n_y / W$. The longitudinal wave vector, which is given by

$$k_x = \pm \sqrt{\left(\frac{\mu}{\hbar v}\right)^2 - k_y^2}, \quad (43)$$

is real for low enough incident angles: $\hbar vk_y < \mu$.

For larger transverse wave vector ($\hbar vk_y > \mu$), the longitudinal wave vector $k_x$ becomes imaginary:

$$k_x = \pm i \kappa_x = \pm i \sqrt{k_y^2 - \left(\frac{\mu}{\hbar v}\right)^2}, \quad (44)$$

and the transmission is given by:

$$T_{k_y}^{(e)} = \frac{\kappa_x^2}{\kappa_x^2 \cosh^2(\kappa_x L) + \left(\frac{\mu}{\hbar v}\right)^2 \sinh^2(\kappa_x L)}. \quad (45)$$

Since the conductance of the ribbon is obtained by summing over all $k_y$ modes, the transport originates from two contributions: the evanescent modes and the propagating ones, such that the total conductance reads:

$$G(\mu) = G_0 \sum_{n_y} T_{k_y} = \frac{G_0 W}{2\pi} \int_{-k_\infty}^{k_\infty} dk_y T_{k_y} \quad (46)$$

$$= \frac{G_0 W}{\pi} \left[ \int_0^{\frac{\mu}{\hbar v}} dk_y T_{k_y}^{(p)} + \int_{\frac{\mu}{\hbar v}}^{k_\infty} dk_y T_{k_y}^{(e)} \right], \quad (47)$$

where $k_\infty$ is the (very large) Fermi wave vector in the leads. Note that the discrete sum has been replaced by an integral over $k_y$ since large width $W$ is assumed. There is a competition between two contributions: for small $\mu$, the transport will be dominated by evanescent modes, whereas the propagating ones dominate for large $\mu$.

At the Dirac point ($\mu = 0$), the conductance originates only from evanescent modes with transmission probability:

$$T_{k_y}^{(e)} = \frac{1}{\cosh^2(k_y L)}, \quad (48)$$

leading to the minimal conductance:

$$G(\mu = 0) = \frac{G_0 W}{\pi L} \int_0^{k_\infty L} \frac{dk_y L}{\cosh^2(k_y L)} = \frac{G_0 W}{\pi L}. \quad (49)$$
<table><caption>TABLE I. Table of the gap sizes and the characteristic length of the corresponding evanescent states for a driving strength of $\beta = 0.1$.</caption><thead><tr><th></th><th colspan="2">Gaps at ε = 0</th><th colspan="2">Gaps at ε = ħω/2</th></tr><tr><th></th><th>m = 0</th><th>m = 2</th><th>m = 1</th><th>m = 3</th></tr></thead><tbody><tr><th>Δ<sub>m</sub>/ħω</th><td>0.02</td><td>0.01</td><td>0.1</td><td>0.0005</td></tr><tr><th>ξ<sub>m</sub>/l<sub>ω</sub></th><td>50</td><td>100</td><td>10</td><td>2000</td></tr></tbody></table>

In this large W regime, the conductance at the Dirac point is proportional to W and $1/L$. It is therefore convenient to normalize the conductance by the ratio $W/L$, as in Fig. 5, green curves, which corresponds to the conductivity $\sigma = GL/W$ of the rectangular ribbon.

## B. Weak driving $\beta = 0.1$

In the presence of electromagnetic radiation, the conductance (blue curves in Fig. 5) exhibits broad dips around chemical potentials $\mu = \pm \hbar\omega/2$. These dips mainly correspond to gaps centered at quasienergy $\epsilon = \pm \hbar\omega/2$, and located near $k = \pm \omega/2v$, in the quasienergy dispersion relation, see Fig. 2(a). These gaps originate from the one-photon resonance between the valence and the conduction band. The electromagnetic coupling leads to an avoided crossing and the opening of a gap $\Delta_1$, associated to a typical decay length $\xi_1 = \hbar v/\Delta_1$. Besides this main gap, the quasienergy spectrum Fig. 2(a) also contains a set of very tiny gaps located at higher wave vectors, around $k = \pm 3\omega/2v$, $\pm 5\omega/2v$, ..., and also all nested around $\pm \hbar\omega/2$. However, for weak driving, the weight of the wave function on these states is negligible, and therefore their contribution to the conductance turns out to be far smaller (than the $\Delta_1$ contribution) although their decay lengths are much larger (nearly propagating states).

The conductance at $\mu = \pm \hbar\omega/2$ is therefore controlled by the ratio between the sample length $L$ and the decay length $\xi_1$ of the evanescent state associated to the gap $\Delta_1$, typically $\Delta_1 \simeq 0.1 \hbar\omega$ and $\xi_1 = 10 l_\omega$ for $\beta = 0.1$ according to Eq. (14) (see also Table I). Indeed the conductance dips become more pronounced and deeper as the length $L$ increases. For $L = 10 l_\omega = \xi_1$, the conductance is significantly suppressed (about 25 percent less than the nonirradiated sample conductance) around $\mu = \pm \hbar\omega/2$ [Fig. 5(a)]. For $L = 25 l_\omega = 2.5 \xi_1$, the conductance is strongly reduced, roughly by a factor 4, with respect to the nonirradiated value [Fig. 5(b)]. Finally when the length exceeds the penetration length by an order of magnitude, namely for $L = 100 l_\omega = 10 \xi_1$, the conductance dips are sharp and well defined [Fig. 5(c)]. Besides, the width of the dips [which is well defined only in Fig. 5(c)] corresponds to the value of the gap $\Delta_1$: $\Delta_1 = 0.1 \hbar\omega$ (see Table I). The remaining conductance originates from the states in the gap $\Delta_3$ that have a very long characteristic length.

Away from these dips, the conductance of the irradiated ribbon follows approximately the nonirradiated one. Here, for $\beta = 0.1$, one can notice that the gaps around $\epsilon = 0$ are very small so their characteristic length are much larger than at $\pm \hbar\omega/2$ (Table I), which means that even for $L = 100 l_\omega$, the deviation from the nonirradiated curve keeps rather small. The conductance near the Dirac point, $\mu = 0$, is discussed in more

<figure><img src="image_9.png" /></figure>

<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /><figcaption>FIG. 5. Conductance of a rectangular graphene ribbon as a function of chemical potential for the nonirradiated case ($\beta = 0$ in green) and for $\beta = 0.1$ in blue ($N = 2$) for a width $W = 250 l_\omega$ and various lengths: (a) $L = 10 l_\omega = \xi_1$, (b) $L = 25 l_\omega = 2.5 \xi_1$, and (c) $L = 100 l_\omega = 10 \xi_1$. For the irradiated case, some dips develop around filling $\mu = \pm \hbar\omega/2$ which corresponds to the gaps of order $m = 1$ at $k = \pm \omega/2v$ in Fig. 2.</figcaption></figure>
detail in Sec. V. We conclude that for weak driving, the effect of the driving on the conductance is mainly observable around $\mu = \pm \hbar\omega/2$.

### C. Bulk and edge conductance versus length $L$ for $\beta = 0.1$ near $\mu = \pm \hbar\omega/2$

#### 1. Bulk states

Inside the dips at $\mu = \pm \hbar\omega/2$, the current is mainly carried by evanescent modes, thus, the conductance of the sample is expected to decrease exponentially with the length. At short length, the conductance in the gap $\epsilon = \hbar\omega/2$ is dominated by the evanescent modes coming from the gap $\Delta_1$, but for larger length, the modes originating from the gap $\Delta_3$ will be dominant. Therefore, it is possible to express the conductivity in the gap $\epsilon = \hbar\omega/2$ with a simple model taking into account this interplay between states at gaps $\Delta_1$ and $\Delta_3$. First, we have checked that the conductance is proportional to $W$, and thus our ansatz reads:

$$G_{\text{bulk}} = G_0 W (ae^{-L/\xi_1} + be^{-L/\xi_3}), \quad (50)$$

where $\xi_1$ and $\xi_3$ are given by Table I, a and b being the only fitting parameters of our model that depend only on $\beta$. The latter parameters (expressed in units of $l_\omega^{-1}$) quantify the relative importance of the evanescent states in transport. Finally, using the rescaling by $G_0 W/L$, we have plotted in Fig. 6(b) the following ratio:

$$\frac{G}{G_0} \frac{L}{W} = L(ae^{-L/\xi_1} + be^{-L/\xi_3}), \quad (51)$$

and extracted the value of a and b. The peak with maximum around $L = 10l_\omega$ corresponds to the current carried by the evanescent modes carried mainly by the evanescent states originating from the first order side band ($\xi_1 = 10l_\omega$). The curve reaches a minimum around $L = 70l_\omega$ and then increases again as the states coming from the third order side band become the dominant source of current. The parameter b is two orders of magnitude smaller than a, which indicates that for $\beta = 0.1$, the weight of the wave function on the evanescent states corresponding to the gap $\Delta_3$ is very weak. This fact corroborates the approximation we made by considering only the evanescent states originating from the gaps $\Delta_1$ and $\Delta_3$.

#### 2. Edge states

So far only transport through 2D bulk states (propagating or evanescent) has been evaluated and discussed. One can evaluate the edge state contribution and compare it with the bulk contribution. The dotted lines correspond to the edge state contribution. There is one edge state per valley that contributes to the dc conductance [16,17,21], taking into account for valley and edge degeneracy, the maximal conductance is equal to:

$$G_{\text{edge}}(\mu = \hbar\omega/2) = \frac{4e^2}{h}. \quad (52)$$

We said “maximal” conductance because the irradiation may reduce the conductance of the edge states. This effect has been studied by Aaron et al. [51] using the Bernevig-Hugues-Zhang model [52]. They found that the conductance is reduced by a Bessel factor with argument $\beta$, in the same fashion as in the photon-assisted transport problem. In this paper, as

<figure><img src="image_12.png" /></figure>

<figure><img src="image_13.png" /><figcaption>FIG. 6. Conductance of a ribbon as a function of the length $L$ of the ribbon for driving strength $\beta = 0.1$ and chemical potential $\mu = \hbar\omega/2$ ($N = 3$) using two different representations: (a) logarithm of the conductance for $W = 100l_\omega$, and (b) conductance normalized by the ratio $W/L$ (independent of the ribbon width $W$). The green smooth curve corresponds to the fit in Eq. (51), and the resulting fitting parameters are a and b. For length $L \lesssim 70l_\omega$, the conductance is carried mainly by the modes with characteristic length $\xi_1$ and for larger length, the modes with length $\xi_3$ dominate.</figcaption></figure>

we consider weak driving, we will not take into account this reduced edge state conductance. Interedge scattering is neglected as only large sample widths $W$ are considered. We have plotted the ratio $G_{\text{edge}}L/G_0W$ for the edge states for different widths of the sample (Fig. 6). This demonstrates that the conductance is dominated by edge states only for long enough samples, namely when $L$ exceeds the typical decay length of the bulk evanescent states. This crossover between 2D transport by evanescent states and 1D edge transport occurs at a typical length $L_c$ which increases with the ribbon width $W$. Using the fit function (50) with the parameters a and b, it is possible to plot $L_c$ as a function of $W$ using the relation:

$$G_{\text{bulk}}(W, L_c) = G_{\text{edge}}. \quad (53)$$
<figure><img src="image_14.png" /><figcaption>FIG. 7. Critical length $L_c$ where bulk and edge have the same contribution to the current at $\mu = \hbar\omega/2$ for a driving strength $\beta = 0.1$. For lengths longer than $L_c$, the current is carried mainly by edge states whereas for shorter lengths, the current is dominated by the evanescent bulk states.</figcaption></figure>

The resulting $L_c(W)$ curve is plotted in Fig. 7. This “phase diagram” allows us to see the competition between bulk and edge current which depends on the shape of the sample and the frequency of irradiation (through $l_\omega$). This curve shows that bulk transport dominates for short and wide graphene ribbons which can be understood qualitatively. Quantitatively, it gives a criterium needed to separate bulk and edge conductance in a particular experiment. There is a clear break in the curve around $W = 500l_\omega$, where the edge conductance becomes smaller than the conductance carried by bulk evanescent states originating from the $\Delta_3$ gap.

### D. Electromagnetic driving $\beta = 0.3$

For $\beta = 0.3$, the scenario is similar to the case with $\beta = 0.1$, except that some dips develop at $\mu = 0$ and $\pm \hbar\omega$ in addition to the dips at $\pm \hbar\omega/2$. Those dips correspond to the second order processes in $\beta$ which are not negligible anymore. Table II shows the gap sizes and the length of the corresponding evanescent states. As the length increases, the dips develop and are better defined.

Each gap is larger at $\beta = 0.3$ than its value at $\beta = 0.1$ case (Fig. 3). Hence the corresponding evanescent state decay lengths have all decreased, thereby making the residual conductance larger, provided $L$ is kept constant. However, the weight of the wave function on the higher order sidebands has

<table><caption>TABLE II. Table of the gaps size and the characteristic length of the corresponding evanescent states for a driving strength of $\beta = 0.3$.</caption><thead><tr><th rowspan="2"></th><th colspan="3">Gaps at ε = 0</th><th colspan="2">Gaps at ε = ħω/2</th></tr><tr><th>m = 0</th><th>m = 2</th><th>m = 4</th><th>m = 1</th><th>m = 3</th></tr></thead><tbody><tr><td>Δ<sub>m</sub>/ħω</td><td>0.16</td><td>0.08</td><td>0.0013</td><td>0.3</td><td>0.0135</td></tr><tr><td>ξ<sub>m</sub>/l<sub>ω</sub></td><td>6.25</td><td>12.5</td><td>761</td><td>3.33</td><td>74.1</td></tr></tbody></table>

<figure><img src="image_15.png" /></figure>

<figure><img src="image_16.png" /></figure>

<figure><img src="image_17.png" /><figcaption>FIG. 8. Conductance of an irradiated graphene ribbon as a function of chemical potential for the nonirradiated case ($\beta = 0$, green curve) and for $\beta = 0.3$ ($N = 2$, blue curve) for a width $W = 250l_\omega$ and a length: (a) $L = 10l_\omega$, (b) $L = 25l_\omega$, and (c) $L = 100l_\omega$. In the irradiated case, some dips develop around filling $\mu = 0$, $\pm \hbar\omega/2$ and $\pm \hbar\omega$ which corresponds to the gaps of order $m = 1$ at $k = \pm \omega/2v$ and the gaps in Fig. 2.</figcaption></figure>
increased. In the gaps at $\mu = \pm \hbar\omega/2$, we can see in Fig. 8(c) that the conductance is not zero anymore because the third order sidebands now carry some current. We can also see the presence of the gap $\Delta_3$ at $\mu = \pm \hbar\omega/2$ inside the gap $\Delta_1$ because there is a small dip in the middle inside the larger dip. Nevertheless, the gap at $\mu = 0$ is perfectly defined for $L = 100l_\omega$ because the current carried by the fourth order sidebands is zero. In the gaps at $\mu = \pm \hbar\omega$, the current is not zero because the density of states is much higher than in the gap at $\mu = 0$.

## V. UNDOPED GRAPHENE

In this section, we analyze more thoroughly the contribution of bulk states to the conductance when the chemical potential is equal to zero. We discuss the evolution of the conductance at the Dirac point as a function of the irradiation strength $\beta$, and as a function of length $L$, comparing it with the edge contribution. We compare our results with those of Gu et al. [12].

### A. Minimal conductivity vs driving strength

At the Dirac point, for $\beta = 0$, the current is carried by evanescent states that belong to the elastic channel. Due to the semimetallic nature of graphene, their decay length is infinite for $k_y = 0$ and decreases as $|k_y|$ increases [Eq. (48)]. For $\beta \neq 0$, the Floquet bands are coupled, so the inelastic channels with $m \neq 0$ are opened and the whole set of evanescent states located at wave vector $k_m^x = \pm m\omega/2v$, with even $m$, will contribute to the current. The evanescent states in the gap at $k = 0$ that originate from the elastic channel $m = 0$ are now gapped, and their decay length is no longer infinite.

The characteristic length of the evanescent states in the gap of order $m$ is given by Eqs. (23) and (24). When $\beta$ increases, their decay length decreases and therefore the current carried by the evanescent state with number of photon $m$ should decrease. However, with increasing $\beta$, the amplitude of the wave function over states with higher number of photons $m$ inside the irradiated region increases. There is a competition between these two effects that creates oscillations of the conductance at the Dirac point when varying $\beta$ (Fig. 9).

To understand this competition, it is useful to plot the conductivity $G_m$ carried by each channel $m$, where $m$ corresponds to the even number $m$ of photons absorbed or emitted, defined as:

$$G_m = \frac{4e^2}{h}(T_m + T_{-m}), \qquad (54)$$

evaluated at $\mu = 0$. We normalize this quantity by $W/L$ in order to compare our results to the nonirradiated case. We realize that the shape of this curve is independent of $W$ which means that the conductance $G$ is a linear function of $W$. This is consistent with the fact that the density of states is proportional to $1/W$. However this curve depends on $L$ because $L$ has to be compared to the characteristic length $\xi$ of the evanescent states (which are functions of $\beta$). The effect of varying $L$ on this curve is to compress or expand the curve along the $\beta$ axis. In the next section, we study the evolution of the conductance as a function of $L$ for a fixed $\beta$.

<figure><img src="image_18.png" /><figcaption>FIG. 9. Conductance $G_m$ at the Dirac point carried by electrons having emitted or absorbed $m$ photons in a function of the driving strength for a length $L = 50l_\omega$ and $N = 5$ Floquet replicas. $G_m$ is defined by Eq. (54) in the text. The competition between the decreasing characteristic length of the evanescent states and the increasing weight of the wave function over replicas with a higher number of photons creates these oscillations.</figcaption></figure>

For $\beta = 0$, we recover the nonirradiated result minimal conductance:

$$G(\mu = 0, \beta = 0) = \frac{G_0 W}{\pi L}. \qquad (55)$$

For $\beta < 0.15$, the two main contributions come from the channels $m = 0$ and $m = 2$ (Fig. 9). As expected, the current carried by the channel $m = 0$ decreases with $\beta$ because: (i) the gap $\Delta_0$ increases and (ii) the weight of the wave function on this channel decreases. We now turn to the current carried by the channel $m = 2$: For $\beta = 0$ it is zero because there is no irradiation, thus it necessarily increases with $\beta$ as we open this channel. We can see in Fig. 9 that the conductance $G_2$ increases faster than $G_0$ decreases, therefore the total conductance increases.

At $\beta = 0.15$, the channel $m = 4$ starts to contribute although the main part of the current is being carried by channel $m = 2$. Then for $\beta > 0.2$, the conductance $G_2$ of this channel $m = 2$ starts to decrease while the conductance $G_4$ increases with $\beta$. Around $\beta = 0.3$, the total conductance $G$ reaches a minimum, while contributions from the $m = 2$ and $m = 4$ channels are almost equal. The same scenario repeats with increasing $\beta$.

### B. Minimal conductivity vs length

For $\beta = 0.3$, we obtain that the conductance can be fitted by a simple model including only three evanescent states [associated with three nested gaps around $\epsilon = 0$ in Fig. 2(a)]:

$$\frac{G}{G_0} \frac{L}{W} = L(ae^{-L/\xi_0} + be^{-L/\xi_2} + ce^{-L/\xi_4}), \qquad (56)$$

where the fitting parameters $a$, $b$, and $c$ are shown on the plot. For $\beta = 0.3$, the conductance originating from evanescent states that belong to the $\Delta_4$ gap is very weak. We checked
<figure><img src="image_19.png" /></figure>

<figure><img src="image_20.png" /><figcaption>FIG. 10. Conductance of a ribbon as a function of the length L of the ribbon for driving strength β = 0.3 and chemical potential µ = 0 (N = 3) using two different representation: (a) logarithm of the conductance for W = 100l_ω and (b) conductance normalized by the ratio W/L (this curve is independent of the ribbon width W). The smooth green curve corresponds to the fit in Eq. (56), and the resulting fitting parameters are a, b, and c. The dashed lines correspond to the edge state conductance given by Eq. (57), which are linear with respect to L because of their ballistic nature.</figcaption></figure>

that adding the conductance states originating from the gap Δ₆ in the fitting function (56) doesn't change the results (for lengths smaller than 800l_ω, according to Table II). Therefore the scaling of the bulk conductance is quite different than the one obtained by Gu et al. [12]. Our model suggests that the bulk conductance is described by Eq. (56) where only the evanescent states with small m (m is the order of the anticrossing) are required, whereas in Gu et al., the whole set of evanescent states originating from the gaps Δₘ including those with high m has to be taken into account, resulting in an approximate power law behavior [12]. The discrepancy is related to the fact that we have considered the regime of small β < 1 which corresponds to moderate laser powers for photon energy below 100 meV to avoid heating. In this regime we

<figure><img src="image_21.png" /><figcaption>FIG. 11. Critical length L_c where bulk and edge have the same contribution to the current at µ = 0 for a driving strength β = 0.3. For lengths longer than L_c, the current is carried mainly by edge states whereas for shorter lengths, the current is dominated by the evanescent bulk states.</figcaption></figure>

find that only a finite number of Floquet replica is needed to get the DC conductance.

## 1. Edge state conductance

So far, only the bulk conductance has been evaluated and it is necessary to explicit the conditions upon which edge state contributions dominates. The dotted lines in Fig. 10 correspond to the edge state conductance. There is one edge state linking the valleys that contributes to dc conductance [12,16,20,21], therefore the maximal edge conductance is equal to:

$$G_{\text{edge}}(\mu = 0) = \frac{2e^2}{h} = G_0/2. \qquad (57)$$

This situation is completely identical to the case studied in Sec. IVC except for the factor 1/2 in the edge conductance and the different expression for the bulk conductance (56). We don't take into account the nonquantization of the edge states due to the driving [51] but we know that the edge states are chiral and robust to disorder, so that their conductance is independent of the size of the sample. We apply the same procedure as in Sec. IVC to see the competition between edge and bulk states and obtain the curve L_c(W) for which the edge and the bulk states conductance is equal (Fig. 11).

# VI. EXPERIMENTAL PARAMETERS

Here we discuss the assumptions made in the previous sections and the conditions necessary to observe the various transport phenomena investigated in this paper.

## A. Dissipation

The Landauer-Büttiker formalism implicitly assumes that dissipation occurs only in the leads. In the irradiated region, the electrons are coherently dressed by one or several quanta of the electromagnetic radiation, and the corresponding excess energy is dissipated in the drain electrode of the transistor.
Clearly, this is only valid when the irradiated region length L is short enough. If one considers coupling to a phonon bath, the relevant time scales are the characteristic time for excited electrons to emit phonons, and the dwell time of the electron through the central irradiated graphene region. In the ballistic limit, this condition gives an upper bound on the sample length

$$L < v\tau_{ph}, \qquad (58)$$

where $\tau_{ph}$ is the average time between electron-phonon inelastic scattering processes. The dominant process for relaxation of excited electrons is the emission of optical phonons of energy 194 or 330 meV [53]. The characteristic time for this process is $\tau_{opt} \approx 1$ ps, which limits the sample length to $L < 1$ µm. For electrons with smaller energy, the only dissipating channel is the acoustic phonon [54,55], which is a much slower process that takes place on a timescale of the order of a nanosecond [55], which corresponds to length $L < 1$ mm. Finally, the assumption of a ballistic conductor requires the length of the sample to be smaller than the mean free path $l_e$ of scattering of electrons by impurities, which in high-quality graphene is $l_e \approx 1-10$ µm.

In conclusion, dissipation by (optical) phonons can be safely neglected for photon energies $\hbar\omega$ smaller than 200 meV. Besides, we note that such photon energies are much smaller than the bandwidth of graphene, so the low energy Dirac equation approximation is valid.

## B. Observability

According to the previous section, we consider photon energies ranging from 10 to 100 meV, namely terahertz frequencies from 2.5 to 25 THz. This frequency window corresponds to characteristic lengths $l_\omega = 6$ nm to 0.6 nm.

Experimentally, the typical electron densities that can be reached in graphene are around $5 \times 10^{12}$ cm$^{-2}$ which corresponds to Fermi energies up to 250 meV. With the electron energies considered, it is experimentally feasible to reach the energy dips at $\pm \hbar\omega/2$.

Finally, the laser intensity required to reach the driving strength $\beta = 1$ (considered as an upper limit for the scope of this paper) is 1 mW µm$^{-2}$ for photons with energy 10 meV and 10$^4$ mW µm$^{-2}$ for photons of 100 meV, which are very

high. This is why we have restricted our plots to values of $\beta$ much smaller than one (typically $\beta = 0.1$ or $\beta = 0.3$).

## VII. CONCLUSIONS

The two-terminal conductance of a rectangular graphene ribbon irradiated by an electromagnetic wave (frequency $\omega$) has been studied using the Floquet theory for experimentally realizable setups. In the ballistic regime our numerical calculations confirm that the coherent dressing of the original Dirac cone leads to the opening of a set of two nonequivalent photoinduced gaps in the Floquet zone, around quasienergies $\epsilon = 0$ and $\epsilon = \hbar\omega/2$. The size of the gaps are found to be in good agreement with RWA estimations [45] at low irradiation $\beta$. When the chemical potential is tuned inside one of these photoinduced gaps, the conductance of the sample decreases drastically because the current is now carried by evanescent states (rather than propagating states when $\beta = 0$). The conductance curve as a function of the chemical follows closely the nonirradiated one except at integer multiples of $\hbar\omega/2$ where broad dips appear. For weak driving, the main effect on the conductance is seen at doping $\mu = \hbar\omega/2$, while for stronger driving additional dips also develop around $\mu = \hbar\omega$. The effect of the irradiation is also effective at the Dirac point where it modulates the value of bulk minimal conductance. The widths of these dips (on the axis of chemical potential) correspond to the sizes of the gaps in the Floquet spectrum, while the depth of these dips (on the conductance axis) depends on the length of the sample. The residual conductance in dips has been fitted by assuming that the current is carried by a few evanescent states characterized by distinct decay lengths.

Besides the transport by 2D bulk states, electrons may also propagate through photoinduced 1D edge states. We have shown that the edge state contribution to the total conductance can be neglected for short enough samples, and a quantitative criterium, for being in the bulk transport regime, has been given.

## ACKNOWLEDGMENTS

We are grateful to Dario Bercioux and Bogusz Bujnowski for useful comments. J.H.B. was supported by the ERC Starting Grant No. 679722.

[1] M. Z. Hasan and C. L. Kane, Rev. Mod. Phys. 82, 3045 (2010).

[2] X.-L. Qi and S.-C. Zhang, Rev. Mod. Phys. 83, 1057 (2011).

[3] T. Oka and H. Aoki, Phys. Rev. B 79, 081406 (2009).

[4] J.-i. Inoue and A. Tanaka, Phys. Rev. Lett. 105, 017401 (2010).

[5] N. H. Lindner, G. Refael, and V. Galitski, Nat. Phys. 7, 490 (2011).

[6] T. Kitagawa, E. Berg, M. Rudner, and E. Demler, Phys. Rev. B 82, 235114 (2010).

[7] M. S. Rudner, N. H. Lindner, E. Berg, and M. Levin, Phys. Rev. X 3, 031005 (2013).

[8] J. H. Shirley, Phys. Rev. 138, B979 (1965).

[9] H. Sambe, Phys. Rev. A 7, 2203 (1973).

[10] J. Cayssol, B. Dóra, F. Simon, and R. Moessner, Phys. Status Solidi RRL 7, 101 (2013).

[11] A. Gómez-León and G. Platero, Phys. Rev. Lett. 110, 200403 (2013).

[12] Z. Gu, H. A. Fertig, D. P. Arovas, and A. Auerbach, Phys. Rev. Lett. 107, 216601 (2011).

[13] H. L. Calvo, H. M. Pastawski, S. Roche, and L. E. F. Foa Torres, Appl. Phys. Lett. 98, 232103 (2011).

[14] H. L. Calvo, P. M. Perez-Piskunow, S. Roche, and L. E. F. Foa Torres, Appl. Phys. Lett. 101, 253506 (2012).

[15] H. L. Calvo, P. M. Perez-Piskunow, H. M. Pastawski, S. Roche, and L. E. F. F. Torres, J. Phys.: Condens. Matter 25, 144202 (2013).
[16] P. M. Perez-Piskunow, G. Usaj, C. A. Balseiro, and L. E. F. Foa Torres, Phys. Rev. B 89, 121401 (2014).

[17] G. Usaj, P. M. Perez-Piskunow, L. E. F. Foa Torres, and C. A. Balseiro, Phys. Rev. B 90, 115423 (2014).

[18] A. Kundu, H. A. Fertig, and B. Seradjeh, Phys. Rev. Lett. 113, 236803 (2014).

[19] A. Quelle, M. Goerbig, and C. Morais Smith, New J. Phys. 18, 015006 (2016).

[20] T. Kitagawa, T. Oka, A. Brataas, L. Fu, and E. Demler, Phys. Rev. B 84, 235108 (2011).

[21] L. E. F. Foa Torres, P. M. Perez-Piskunow, C. A. Balseiro, and G. Usaj, Phys. Rev. Lett. 113, 266801 (2014).

[22] M. Fruchart, P. Delplace, J. Weston, X. Waintal, and D. Carpentier, Phys. E (Amsterdam, Neth.) 75, 287 (2016).

[23] H. Dehghani, T. Oka, and A. Mitra, Phys. Rev. B 90, 195429 (2014).

[24] H. Dehghani, T. Oka, and A. Mitra, Phys. Rev. B 91, 155422 (2015).

[25] H. Dehghani and A. Mitra, Phys. Rev. B 92, 165111 (2015).

[26] H. Dehghani and A. Mitra, Phys. Rev. B 93, 205437 (2016).

[27] K. I. Seetharam, C.-E. Bardyn, N. H. Lindner, M. S. Rudner, and G. Refael, Phys. Rev. X 5, 041050 (2015).

[28] Y. H. Wang, H. Steinberg, P. Jarillo-Herrero, and N. Gedik, Science 342, 453 (2013).

[29] G. Jotzu, M. Messer, R. Desbuquois, M. Lebrat, T. Uehlinger, D. Greif, and T. Esslinger, Nature (London) 515, 237 (2014).

[30] M. C. Rechtsman, J. M. Zeuner, Y. Plotnik, Y. Lumer, D. Podolsky, F. Dreisow, S. Nolte, M. Segev, and A. Szameit, Nature (London) 496, 196 (2013).

[31] E. Prada, P. San-Jose, and H. Schomerus, Phys. Rev. B 80, 245414 (2009).

[32] P. San-Jose, E. Prada, S. Kohler, and H. Schomerus, Phys. Rev. B 84, 155408 (2011).

[33] P. San-Jose, E. Prada, H. Schomerus, and S. Kohler, Appl. Phys. Lett. 101, 153506 (2012).

[34] S. V. Syzranov, M. V. Fistul, and K. B. Efetov, Phys. Rev. B 78, 045407 (2008).

[35] M. V. Fistul, S. V. Syzranov, A. M. Kadigrobov, and K. B. Efetov, Phys. Rev. B 82, 121409 (2010).

[36] S. V. Syzranov, Y. I. Rodionov, K. I. Kugel, and F. Nori, Phys. Rev. B 88, 241112 (2013).

[37] M. V. Fistul and K. B. Efetov, Phys. Rev. B 90, 125416 (2014).

[38] Y. I. Rodionov, K. I. Kugel, and F. Nori, Phys. Rev. B 94, 195108 (2016).

[39] F. D. Parmentier, L. N. Serkovic-Loli, P. Roulleau, and D. C. Glattli, Phys. Rev. Lett. 116, 227401 (2016).

[40] V. E. Calado, S.-E. Zhu, S. Goswami, Q. Xu, K. Watanabe, T. Taniguchi, G. C. A. M. Janssen, and L. M. K. Vandersypen, Appl. Phys. Lett. 104, 023103 (2014).

[41] L. Banszerus, M. Schmitz, S. Engels, M. Goldsche, K. Watanabe, T. Taniguchi, B. Beschoten, and C. Stampfer, Nano Lett. 16, 1387 (2016).

[42] L. Wang, I. Meric, P. Y. Huang, Q. Gao, Y. Gao, H. Tran, T. Taniguchi, K. Watanabe, L. M. Campos, D. A. Muller, J. Guo, P. Kim, J. Hone, K. L. Shepard, and C. R. Dean, Science 342, 614 (2013).

[43] A. A. Balandin, Nat. Mater. 10, 569 (2011).

[44] K. F. Mak, L. Ju, F. Wang, and T. F. Heinz, Solid State Commun. 152, 1341 (2012).

[45] Y. Zhou and M. W. Wu, Phys. Rev. B 83, 245436 (2011).

[46] A. H. Castro Neto, F. Guinea, N. M. R. Peres, K. S. Novoselov, and A. K. Geim, Rev. Mod. Phys. 81, 109 (2009).

[47] M. O. Goerbig, Rev. Mod. Phys. 83, 1193 (2011).

[48] M. Moskalets and M. Büttiker, Phys. Rev. B 66, 205320 (2002).

[49] J. Tworzydło, B. Trauzettel, M. Titov, A. Rycerz, and C. W. J. Beenakker, Phys. Rev. Lett. 96, 246802 (2006).

[50] S. Kohler, J. Lehmann, and P. Hänggi, Phys. Rep. 406, 379 (2005).

[51] A. Farrell and T. Pereg-Barnea, Phys. Rev. Lett. 115, 106403 (2015).

[52] B. A. Bernevig, T. L. Hughes, and S.-C. Zhang, Science 314, 1757 (2006).

[53] D. Sun, Z.-K. Wu, C. Divin, X. Li, C. Berger, W. A. de Heer, P. N. First, and T. B. Norris, Phys. Rev. Lett. 101, 157402 (2008).

[54] W.-K. Tse and S. Das Sarma, Phys. Rev. B 79, 235406 (2009).

[55] R. Bistritzer and A. H. MacDonald, Phys. Rev. Lett. 102, 206410 (2009).