## Abstract

Recently, the temporal control of graphene carrier density has emerged as a viable means to create various frequency shifting, modulation, and sensing photonic devices. Here we describe a general theoretical approach to calculate the graphene plasmon transformation after rapid changes of the Fermi level and carrier density. The approach is based on solving the Maxwell equations supplemented by the microscopic current equation. We derive formulas for the amplitudes of the transmitted and reflected plasmons after a rapid carrier density drop. The relation of these amplitudes and the Fourier transformed finite-difference time-domain (FDTD) fields is also established by introducing the concept of differential spectral transformation of wavepackets. The results of the analytical and FDTD approaches refute the claims of plasmon amplification under rapid carrier changes that appeared in recent theoretical studies. The presented theoretical and computational approaches form a basis of time-varying electromagnetics of graphene plasmonics.

© 2018 Optical Society of America under the terms of the OSA Open Access Publishing Agreement

## 1. INTRODUCTION

The ability to control guided-mode propagation through the refractive index allows one to engineer various optoelectronic devices, for example, ultrafast modulators in telecommunications [1]. The refractive index can be changed, for example, by using an external electric field (electro-refraction) or by manipulating the carrier density (plasma-dispersion effect). The smallness of these effects in conventional materials stimulates researchers to look for novel materials with enhanced nonlinear properties. The appearance of graphene—a promising candidate for photonic applications in infrared and terahertz ranges—has attracted significant interest in both experimental studies and their theoretical understanding as reviewed recently [2,3].

Graphene sheets with finite carrier density are capable of supporting propagating modes—graphene plasmons [4–6]. Their propagation can be directed by creating various sheet corrugations [7], curved landscapes [8], or topological effects [9]. Graphene-based nanowaveguides hold promise as printed interconnects [10]. Strong nonlinearities in graphene can enable all-optical modulation [11,12]. Applying gate voltage one can also tune graphene conductivity and realize two-dimensional (2D) transformation optics using graphene [13], dynamically controlled devices [14], dynamically tunable and active hyperbolic metamaterials [15], and reconfigurable terahertz nanodevices [16]. The conductivity control is an apparent advantage of graphene over metal layers and metal interfaces [13] in addition to lower losses and stronger confinement [17], which can further be improved by making heterostructures [18–20]. Depending on the characteristic time scale, the conductivity changes can be slow or fast compared to the plasmon period. In particular, the variation of graphene properties on a scale of a few 10s of femtoseconds was demonstrated [21–23]. The instantaneous response of graphene to ultrafast optical fields, elucidating the role of hot carriers on sub-100 fs time scales, has also been recently studied [24]. In general, the control by external voltage in graphene extends other directions of nonlinear control of optical properties of 2D materials, for example, using terahertz fields for ultrafast optical modulation in semiconductor quantum wells [25,26].

The feasibility of controlling graphene properties in time necessitates the development of a theoretical understanding of plasmon transformation under time variations of graphene parameters, in particular, for rapid density changes. In general, the transformation of waves in time-varying media has a rather long history. Most of the studies considered the transformation of waves in bulk materials [27–29]. The scattering of surface waves under rapid carrier injection was also investigated for a plasma half-space [30,31] and a plasma layer [32]. The consideration of wave scattering, either bulk or surface, requires a need for careful consideration of material equations for the medium, for example, plasma [33] or polar molecules [34].

The problem of temporal graphene plasmon transformation was addressed very recently in Refs. [35,36]. Among other regimes, these studies considered step-like variations of the Fermi level and calculated the amplitudes of the surface plasmons excited at the temporal discontinuity. A particularly noticeable prediction is the appearance of plasmon amplification. In Ref. [35], the plasmon amplification occurs if the Fermi level rapidly decreases. In contrast, in Ref. [36], the plasmon amplification was predicted for rapid density increases. In our opinion, there are at least two issues with these predictions. First, solving the same problem, these studies provide different transmission and reflection coefficients for the plasmons for step-like variations of carrier density; compare Eq. (10) in Ref. [35] and Eq. (7) in Ref. [36]. Second, and more fundamental, the appearance or disappearance of carriers under a rapid Fermi level variation cannot inject electromagnetic energy into the plasmon. Energy injection is possible under several scenarios: carriers recombine producing photons (or in this case, plasmons) or external forces perform work on the system (parametric processes). To achieve plasmon emission one essentially needs negative conductivity, or population inversion. Although it is possible to obtain in graphene [37,38], this was not present in the problem considered in Refs. [35,36]. Parametric processes are claimed to be responsible for the amplification in Ref. [36]. The energy of a graphene plasmon consists of the energy of the electromagnetic field and the kinetic energy of the moving carriers (current). Since the change of density neither increases the energy of the electromagnetic field nor the kinetic energy of the carriers, there is no parametric mechanism that can inject energy into the plasmons. A resolution of these paradoxes is important both for theoretical advancement of ultrafast graphene plasmonics as well as for creation of novel optoelectronic devices.

In this paper we solve the problem of graphene plasmon scattering under a rapid time variation of carrier density using the Maxwell equations and the microscopic current equation. Our model is not limited by the quasi-static plasmon approximation adopted in Refs. [35,36]. While our general approach is applicable to either decrease or increase, we consider explicitly only the case of rapid decrease. The small time scale of the changes also justifies the neglect of collisions, which we adopted here for the sake of shedding some light on energy redistribution after the temporal discontinuity. We investigate the frequencies of the transformed plasmons and derive formulas for their amplitudes. A comparison of the energies of the excited modes shows that the plasmon amplification does not exist at temporal discontinuities, which is in contrast to the predictions of Refs. [35,36]. We also study the plasmon scattering using finite-difference time-domain (FDTD) simulations and present a correct procedure to relate the FDTD Fourier transforms to the transmission and reflection coefficients for monochromatic waves. Finally, we discuss the culprits behind the predictions of plasmon amplification in Refs. [35,36].

## 2. ANALYTICAL SOLUTION OF THE SCATTERING ON TIME DISCONTINUITY

#### A. Dispersion and Energy of Graphene Plasmons for Constant Parameters

We take a surface plasmon at frequency $\omega $ propagating along a graphene sheet surrounded by a dielectric with permittivity $\u03f5$; see Fig. 1(a). In all numerical results we set $\u03f5=1$. The electromagnetic field of a transverse magnetic (TM) plasmon has three components $\{{H}_{z},{E}_{x},{E}_{y}\}$. The tangential ${E}_{x}$ component can be written as

where $h$ is the propagation wavevector, $\varkappa $ is the decay constant, and ${E}_{0}$ is the amplitude, which can be assumed to be real. The magnetic field at $y={0}^{+}$ is ${H}_{0}=i\omega {E}_{0}/(c\varkappa )$. The equation for the graphene current in the Drude model is where $\mathrm{\Omega}$ is the frequency parameter characterizing the graphene response. It conveniently includes the carrier concentration, temperature, band structure, etc. Specific values for $\mathrm{\Omega}$ are not required since the solution will depend only on the ratio $\omega /\mathrm{\Omega}$. For graphene, this parameter is $\mathrm{\Omega}\sim {E}_{f}\sim \sqrt{N}$, where ${E}_{f}$ is the Fermi level and $N$ is the electron density. For a regular electron gas, $\mathrm{\Omega}\sim N$. In Eq. (2) we neglected collisions for the sake of simplicity. Using the standard boundary conditions across the current sheet, we arrive at the following dispersion for the TM plasmon:The time-averaged plasmon energy ${W}_{0}={W}_{f}+{W}_{k}$ consists of the energies of its electric ${W}_{E}$ and magnetic ${W}_{H}$ fields (${W}_{f}={W}_{E}+{W}_{H}$) outside of the sheet and the kinetic energy ${W}_{k}$ of the oscillating graphene carriers:

The dependence of the relative energy components on frequency is shown in Fig. 2(b). With increasing $\omega /\mathrm{\Omega}$, which corresponds to increasing plasmon confinement, the magnetic component decreases while the kinetic one increases.

#### B. Formulation of the Temporal Scattering Problem

We assume that initially, at $t<0$, there is a surface plasmon with properties described in Section 2.A propagating along a graphene sheet. We denote the parameters of the initial plasmon by subscript 1: ${\omega}_{1}$, ${\mathrm{\Omega}}_{1}$, ${h}_{1}$, and ${\varkappa}_{1}$. The amplitudes of ${E}_{x}$ and ${H}_{z}$ fields are ${E}_{0}$ and ${H}_{0}$, respectively. At $t=0$, the Fermi level rapidly decreases and, therefore, the carrier density drops from ${N}_{1}$ to ${N}_{2}$; see Fig. 1(b). This translates into a change from ${\mathrm{\Omega}}_{1}$ to ${\mathrm{\Omega}}_{2}$, which can be described by the jump parameter

as in Refs. [35,36]. We are interested in finding the fields created by the temporal discontinuity defined by Eq. (5).#### C. Frequencies of Scattered Waves

The scattering results in the creation of two plasmons and transient bulk radiation going to $y\to \pm \infty $. One plasmon propagates in the $+x$ direction (transmitted) and the other propagates in the $-x$ direction (reflected). Their frequencies can be found using the invariance of the spatial dependence ${e}^{i{h}_{1}x}$ set by the initial plasmon and their dispersion relation at ${\mathrm{\Omega}}_{2}$, as illustrated in Fig. 1(c). The initial dispersion at ${\mathrm{\Omega}}_{1}$ defines the relation between ${h}_{1}$ and ${\omega}_{1}$. After the density changes, the dispersion for ${\mathrm{\Omega}}_{2}$ gives two frequencies $\pm {\omega}_{2}$ that correspond to the fixed ${h}_{1}$. One frequency describes the transmitted plasmon, the other—reflected. For $\gamma <1$, the plasmons will be frequency downshifted while the bulk radiation will always be frequency upshifted. We will analyze the scattering at several values of ${\omega}_{1}/{\mathrm{\Omega}}_{1}$ corresponding to various degrees of plasmon localization; see the dispersion curves in Fig. 2(a).

Figure 3(a) shows the change in the relative plasmon frequency ${\omega}_{2}/{\omega}_{1}$ as a function of ${\mathrm{\Omega}}_{2}/{\mathrm{\Omega}}_{1}$. For all values of ${\omega}_{1}/{\mathrm{\Omega}}_{1}$, the frequency decreases with decreasing ${\mathrm{\Omega}}_{2}/{\mathrm{\Omega}}_{1}$; in the quasi-static regime ${\omega}_{1}/{\mathrm{\Omega}}_{1}\gg 1$, ${\omega}_{2}/{\omega}_{1}\approx \gamma $, as in Ref. [35]. It is quite interesting that as the carrier density ${N}_{2}$ decreases, the ratio ${\omega}_{2}/{\mathrm{\Omega}}_{2}$ increases, making the excited plasmon more and more confined; see Fig. 3(b).

#### D. Amplitudes of Scattered Waves

To find the scattered fields $\propto {e}^{i{h}_{1}x}$, we use the Maxwell equations (in Gaussian units)

To solve Eqs. (6) and (8) with the specified initial conditions, we use the Laplace transform technique [30,32]. From the system of linear equations, we express the transform ${\tilde{E}}_{x}(x,y,p)$ of ${E}_{x}(x,y,t)$, where $p$ is the Laplace variable, as

Evaluating the residues at $p=\mp i{\omega}_{2}$, which is the solution of $D(p)=0$, we obtain the amplitudes of the transmitted and reflected plasmons:

#### E. Energy Balance

The energy balance not only gives us physical insight into the temporal scattering but also provides justification for the solution. In the absence of collisions, the energy of the initial plasmon ${W}_{0}$ should transform into the energy of the excited (transmitted ${W}_{t}$ and reflected ${W}_{r}$) plasmons and transient bulk radiation ${W}_{b}$. Also, the removed carriers had some kinetic energy, which can be considered as loss ${W}_{l}$. The energy balance becomes

The plasmon energies can be directly calculated using Eqs. (4) and (11). The loss is, apparently, ${W}_{l}=(1-{\mathrm{\Omega}}_{2}/{\mathrm{\Omega}}_{1}){W}_{k}$. The bulk energy ${W}_{b}$ can be calculated from the branch cut integrals and represented in terms of its angular density ${w}_{b}(\theta )$ usingThe results for several values of ${\omega}_{1}/{\mathrm{\Omega}}_{1}$ are shown in Fig. 4. In general, the plasmon transformation depends strongly not only on the jump ${\mathrm{\Omega}}_{2}/{\mathrm{\Omega}}_{1}$ but also on the initial condition ${\omega}_{1}/{\mathrm{\Omega}}_{1}$. The lost energy ${W}_{l}$ increases linearly with decreasing ${\mathrm{\Omega}}_{2}/{\mathrm{\Omega}}_{1}$. Furthermore, the loss becomes higher with increasing confinement of the initial plasmon or ${\omega}_{1}/{\mathrm{\Omega}}_{1}$ as expected from the energy distribution in the initial plasmon; see Fig. 2(b). In contrast, the transient radiation ${W}_{b}$ decreases with increasing confinement. The transmitted plasmon energy ${W}_{t}$ always decreases with decreasing ${\mathrm{\Omega}}_{2}/{\mathrm{\Omega}}_{1}$. The reflected plasmon energy ${W}_{r}$ increases with decreasing ${\mathrm{\Omega}}_{2}/{\mathrm{\Omega}}_{1}$. Interestingly, the excitation of the plasmons takes place even when the carrier density becomes arbitrarily small but finite. This process is especially pronounced when the initial plasmon is strongly confined ${\omega}_{1}/{\mathrm{\Omega}}_{1}=15$. In this case, when most of the existing carriers disappear, the electromagnetic energy, which is ${W}_{f}/{W}_{0}=0.5$, seems to drive the remaining carriers so that both transmitted and reflected plasmons are excited almost equally, ${W}_{t}/{W}_{0}\approx {W}_{r}/{W}_{0}\approx 0.25$.

The transient radiation propagates from the graphene sheet, and its angular distribution is described by Eq. (15). Figure 5 shows typical angular distributions. The weakly localized plasmon, ${\omega}_{1}/{\mathrm{\Omega}}_{1}=1$, produces a rather narrow radiation peak at small angle. The strongly localized plasmon, ${\omega}_{1}/{\mathrm{\Omega}}_{1}=15$, emits almost vertically in a very broad angular range.

#### F. Temporal Scattering in Quasi-Static Approximation

Since the quasi-static approximation is often used, it is instructive to derive the corresponding transmission and reflection coefficients without using the general formulas in Eqs. (11)–(13). It was shown in Section 2.A that in this limit the plasmon energy consists of the electric and kinetic parts. The magnetic field is rather small $|{H}_{z}/{E}_{x}|=2\pi \mathrm{\Omega}/(\u03f5\omega )\ll 1$. According to Fig. 4(a), the temporal discontinuity produces scattered (transmitted and reflected) plasmons and very weak bulk radiation. Also, a significant part of the initial plasmon energy is taken away by the removed carriers. Thus, we can write the electric field of the scattered plasmons at $t>0$ as

#### G. Analysis of Transmission and Reflection Coefficients

Here we compare our results with those of Refs. [35,36] where the scattering of a graphene plasmon on a temporal discontinuity was considered in the quasi-static approximation, i.e., for a strongly confined plasmon. In Ref. [35], matching the field distributions of the incident plasmon and excited (reflected and transmitted) plasmons using the continuity of ${H}_{z}$ and its time derivative resulted in, see Eq. (10) in Ref. [35],

Having discussed the issues with Eqs. (19) and (20), let us study how our results, given by Eqs. (11)–(13), compare with them; see Fig. 6. One can clearly see that our results agree neither with Eq. (19) nor with Eq. (20). In contrast to increasing ${t}_{H}$ with decreasing $\gamma $ obtained in Ref. [35], our ${t}_{H}$ decreases discarding the claim of gain appearance; see also Fig. 4(a). Moreover, in contrast to ${r}_{H}\approx {t}_{H}\approx 1/2$ at ${\mathrm{\Omega}}_{2}/{\mathrm{\Omega}}_{1}\to 0$, our ${r}_{H}$ and ${t}_{H}$ vanish as expected. In general, it is quite possible to obtain $|{t}_{E,H}|>1$ since the coefficients define only the field values at $y=0$. Indeed, our simulations show that, for example, for ${\omega}_{1}/{\mathrm{\Omega}}_{1}=1$ we obtain $|{t}_{E,H}|>1$ at $\gamma =0.5$. However, due to the change of the mode profile, the plasmon energy always decreases.

Where is the origin of mistakes in Refs. [35,36]? Both of them considered the problem explicitly in the quasi-static approximation. In this approximation, the modes are defined by their electric fields and surface currents. Thus, one should use the boundary conditions for these components. However, the behavior of the microscopic current at the temporal discontinuity, see Eq. (7), was not employed in Refs. [35,36] as a necessary condition to assess the plasmon transformation.

## 3. EXTRACTING TRANSMISSION FROM FDTD SIMULATIONS WITH TEMPORAL DISCONTINUITIES

While analytical results can be derived in some limited geometries, numerical simulations are required in most practical situations related to device modeling. From this point of view, it is instructive to compare the analytical results with FDTD simulations. We applied an in-house developed 2D FDTD code [40] to simulate the temporal scattering of a graphene plasmon. Initially, a plasmon wavepacket propagating in the $+x$ direction along the graphene sheet is created; see Fig. 7, frame $t=0$. The current discontinuity at $t=0$, see Eq. (7), gives rise to its scattering. Figure 7, frame $t=0.15/{\mathrm{\Omega}}_{1}$shows a typical distribution of ${E}_{x}(x,y)$ after the jump when some transient radiation escapes from the sheet region. Subsequently, one observes the separation of the field into the transmitted and reflected plasmon wavepackets; see Fig. 7, frame $t=7/{\mathrm{\Omega}}_{1}$.

The spectra of the electric and magnetic fields of the transmitted plasmon were computed at some specific points near the sheet using the running Fourier transform and plotted together with the incident spectra in Fig. 8(a). According to the analytical dependence of the frequency conversion ${\omega}_{2}({\omega}_{1})$ shown in Fig. 8(b), the components at ${\omega}_{1}/{\mathrm{\Omega}}_{1}=15$ should be transformed to ${\omega}_{2}/{\mathrm{\Omega}}_{1}=11.8$; this agrees well with the simulation result in Fig. 8(a). Focusing on amplitudes, according to Eqs. (11)–(13), see also Fig. 6, at ${\omega}_{1}/{\mathrm{\Omega}}_{1}=15$ and ${\mathrm{\Omega}}_{2}/{\mathrm{\Omega}}_{1}=0.6$ we should have ${t}_{E}={E}_{x}^{t}/{E}_{0}=0.924$ and ${t}_{H}={H}_{z}^{t}/{H}_{0}=0.704$. None of these is seen in Fig. 8(a). Instead, Fig. 8(a) shows that the spectrum of the electric field of the transmitted plasmon is higher than of the incident; the magnetic spectrum is only slightly lower. Apparently, the textbook procedure of simply dividing the spectra does not work. Furthermore, the transmitted spectra are visibly narrower than the incident spectra. In fact, to find the transmission coefficient from the FDTD results one needs not only to apply the spectral shift but also the spectral compression or expansion.

Let us define the magnetic field spectrum of the incident plasmon wavepacket by ${H}_{z}^{i}(\omega )$ and the transmitted spectrum by ${H}_{z}^{t}(\omega )$. How can we find the transmission coefficient for a monochromatic plasmon at some fixed ${\omega}_{1}$ within the spectral range of ${H}_{z}^{i}(\omega )$? After a temporal jump, a monochromatic plasmon at ${\omega}_{1}$ will create a plasmon at ${\omega}_{2}$. Let us consider a small spectral interval $\mathrm{\Delta}{\omega}_{1}$ of the wavepacket spectrum. The corresponding spectral components will be transformed into $\mathrm{\Delta}{\omega}_{2}$ after the temporal discontinuity according to

Thus, according to derived Eq. (22) the transmission coefficients ${t}_{H}$ or ${t}_{E}$ for monocromatic waves can be found from the wavepacket spectra by taking the ratio of the spectral components at corresponding frequencies and multiplying by the spectral transformation factor $f$. Figure 8(c) shows the spectral transformation factor calculated analytically based on the frequency dependence ${\omega}_{2}({\omega}_{1})$ shown in Fig. 8(b). In this specific case, this factor differs significantly from 1 and therefore plays a crucial role in extracting the transmission coefficients from the FDTD simulations. Figure 8(d) compares the analytical results obtained from Eq. (11) and the numerical results calculated using Eq. (22). Accounting for the spectral transformation allowed us to match perfectly the analytical results in the whole spectral range of the incident wavepacket. We also note that Ref. [35] presented the spectra obtained from time-domain finite element method (FEM) modeling. However, the transmission was evaluated as the ratio of the peak heights. This may have led to the erroneous estimation of the transmission coefficient since the value of the spectral transformation was not taken into account.

The described procedure can also be applied to quasi-monochromatic wavepackets if the transmission or reflection coefficients at the central frequency, rather than a finite range, are needed. This allows one to integrate Eq. (21) over the frequencies and pull ${t}_{H}$ out of the integral

## 4. CONCLUSION

To conclude, the temporal scattering of a graphene plasmon by a rapid carrier density change was considered. By solving the Maxwell equations supplemented by the material equation for the sheet current using the Laplace transform technique we obtained general formulas for the temporal transmission and reflection coefficients as well as for the transient bulk radiation in the case of carrier density decrease. A procedure to derive the coefficients in the quasi-static limit using the boundary conditions for the electric field and current was proposed and its results were verified by a comparison with the general formulas. It was shown that the energy of the excited plasmons and bulk radiation is smaller than the initial plasmon energy by the amount of the kinetic energy of the removed carriers. This refutes the claim of plasmon amplification in Ref. [35]. The presence of carrier collisions, which are neglected in this study, will lead to the damping of the excited surface plasmons but should not significantly change the energy balance immediately after the rapid density change.

When new carriers are created, one has to solve separately for the dynamics of the existing and newly created carriers since they have different initial velocities. Among other things, this leads to the appearance of the so-called free-streaming mode, which consists of a self-consistent distribution of direct current and magnetic field [41]. This mode takes some energy of the initial plasmon. Thus, a rapid increase of carrier density always reduces the energy, in contrast to the predictions of Ref. [36].

The temporal scattering was also modeled using the FDTD approach. The procedure to find the temporal transmission and reflection coefficients was developed and its results agreed with the analytical ones. It was shown that finding the temporal transmission and reflection coefficients from the FDTD simulations of wavepacket propagation requires not only finding the Fourier transforms of the wavepackets but also accounting for their differential spectral compression or expansion. Besides scattering, rapid changes in graphene, as in plasma layers [42], can create a viable mechanism of launching graphene plasmons using optical beams without using various corrugations or coupling prisms.

## Funding

Ministry of Education and Science of the Russian Federation (Minobrnauka) (3.3854.2017/4.6).

## REFERENCES

**1. **G. T. Reed, G. Mashanovich, F. Y. Gardes, and D. J. Thomson, “Silicon optical modulators,” Nat. Photonics **4**, 518–526 (2010). [CrossRef]

**2. **F. J. García de Abajo, “Special issue “2D materials for nanophotonics”,” ACS Photon. **4**, 2959–2961 (2017). [CrossRef]

**3. **Q. Guo, C. Li, B. Deng, S. Yuan, F. Guinea, and F. Xia, “Infrared nanophotonics based on graphene plasmonics,” ACS Photon. **4**, 2989–2999 (2017). [CrossRef]

**4. **A. N. Grigorenko, M. Polini, and K. S. Novoselov, “Graphene plasmonics,” Nat. Photonics **6**, 749–758 (2012). [CrossRef]

**5. **J. Chen, M. Badioli, P. Alonso-González, S. Thongrattanasiri, F. Huth, J. Osmond, M. Spasenović, A. Centeno, A. Pesquera, N. Camara, F. J. García de Abajo, R. Hillenbrand, and F. H. L. Koppens, “Optical nano-imaging of gate-tunable graphene plasmons,” Nature **487**, 77–81 (2012). [CrossRef]

**6. **Z. Fei, A. S. Rodin, G. O. Andreev, W. Bao, A. S. McLeod, M. Wagner, L. M. Zhang, Z. Zhao, M. Thiemens, G. Dominguez, M. M. Fogler, A. H. C. Neto, C. N. Lau, F. Keilmann, and D. N. Basov, “Gate-tuning of graphene plasmons revealed by infrared nano-imaging,” Nature **487**, 82–85 (2012). [CrossRef]

**7. **T. M. Slipchenko, M. L. Nesterov, R. Hillenbrand, A. Y. Nikitin, and L. Martín-Moreno, “Graphene plasmon reflection by corrugations,” ACS Photon. **4**, 3081–3088 (2017). [CrossRef]

**8. **D. Smirnova, S. H. Mousavi, Z. Wang, Y. S. Kivshar, and A. B. Khanikaev, “Trapping and guiding surface plasmons in curved graphene landscapes,” ACS Photon. **3**, 875–880 (2016). [CrossRef]

**9. **D. A. Kuzmin, I. V. Bychkov, V. G. Shavrov, and V. V. Temnov, “Plasmonics of magnetic and topological graphene-based nanostructures,” Nanophotonics **7**, 597–611 (2018). [CrossRef]

**10. **G. Lovat, P. Burghignoli, and R. Araneo, “Low-frequency dominant-mode propagation in spatially dispersive graphene nanowaveguides,” IEEE Trans. Electromagn. Compat. **55**, 328–333 (2013). [CrossRef]

**11. **W. Li, B. Chen, C. Meng, W. Fang, Y. Xiao, X. Li, Z. Hu, Y. Xu, L. Tong, H. Wang, W. Liu, J. Bao, and Y. R. Shen, “Ultrafast all-optical graphene modulator,” Nano Lett. **14**, 955–959 (2014). [CrossRef]

**12. **J. D. Cox and F. J. García de Abajo, “Transient nonlinear plasmonics in nanostructured graphene,” Optica **5**, 429–433 (2018). [CrossRef]

**13. **A. Vakil and N. Engheta, “Transformation optics using graphene,” Science **332**, 1291–1294 (2011). [CrossRef]

**14. **T.-T. Kim, H.-D. Kim, R. Zhao, S. S. Oh, T. Ha, D. S. Chung, Y. H. Lee, B. Min, and S. Zhang, “Electrically tunable slow light using graphene metamaterials,” ACS Photon. **5**, 1800–1807 (2018). [CrossRef]

**15. **J. S. T. Smalley, F. Vallini, X. Zhang, and Y. Fainman, “Dynamically tunable and active hyperbolic metamaterials,” Adv. Opt. Photon. **10**, 354–408 (2018). [CrossRef]

**16. **P.-Y. Chen, H. Huang, D. Akinwande, and A. Alù, “Graphene-based plasmonic platform for reconfigurable terahertz nanodevices,” ACS Photon. **1**, 647–654 (2014). [CrossRef]

**17. **M. Jablan, M. Soljačić, and H. Buljan, “Plasmons in graphene: fundamental properties and potential applications,” Proc. IEEE **101**, 1689–1704 (2013). [CrossRef]

**18. **A. V. Kretinin, Y. Cao, J. S. Tu, G. L. Yu, R. Jalil, K. S. Novoselov, S. J. Haigh, A. Gholinia, A. Mishchenko, M. Lozada, T. Georgiou, C. R. Woods, F. Withers, P. Blake, G. Eda, A. Wirsig, C. Hucho, K. Watanabe, T. Taniguchi, A. K. Geim, and R. V. Gorbachev, “Electronic properties of graphene encapsulated with different two-dimensional atomic crystals,” Nano Lett. **14**, 3270–3276 (2014). [CrossRef]

**19. **A. Woessner, M. B. Lundeberg, Y. Gao, P. A.-G. A. Principi, M. Carrega, K. Watanabe, M. P. T. Taniguchi, G. Vignale, J. Hone, R. Hillenbrand, and F. H. L. Koppens, “Highly confined low-loss plasmons in graphene-boron nitride heterostructures,” Nat. Mater. **14**, 421–425 (2015). [CrossRef]

**20. **J. D. Caldwell and K. S. Novoselov, “Mid-infrared nanophotonics,” Nat. Mater. **14**, 364–366 (2015). [CrossRef]

**21. **I. Gierz, J. C. Petersen, M. Mitrano, C. Cacho, I. C. E. Turcu, E. Springate, A. Stöhr, A. Köhler, U. Starke, and A. Cavalleri, “Snapshots of non-equilibrium Dirac carrier distributions in graphene,” Nat. Mater. **12**, 1119–1124 (2013). [CrossRef]

**22. **S. Ulstrup, J. C. Johannsen, F. Cilento, J. A. Miwa, A. Crepaldi, M. Zacchigna, C. Cacho, R. Chapman, E. Springate, S. Mammadov, F. Fromm, C. Raidel, T. Seyller, F. Parmigiani, M. Grioni, P. D. C. King, and P. Hofmann, “Ultrafast dynamics of massive Dirac fermions in bilayer graphene,” Phys. Rev. Lett. **112**, 257401 (2014). [CrossRef]

**23. **M. Trushin, A. Grupp, G. Soavi, A. Budweg, D. De Fazio, U. Sassi, A. Lombardo, A. C. Ferrari, W. Belzig, A. Leitenstorfer, and D. Brida, “Ultrafast pseudospin dynamics in graphene,” Phys. Rev. B **92**, 165429 (2015). [CrossRef]

**24. **M. Baudisch, A. Marini, J. D. Cox, T. Zhu, F. Silva, S. Teichmann, M. Massicotte, F. Koppens, L. S. Levitov, F. J. García de Abajo, and J. Biegert, “Ultrafast nonlinear optical response of Dirac fermions in graphene,” Nat. Commun. **9**, 1018 (2018). [CrossRef]

**25. **A. V. Maslov and D. S. Citrin, “Numerical calculation of the terahertz field-induced changes in the optical absorption in quantum wells,” IEEE J. Sel. Top. Quantum Electron. **8**, 457–463 (2002). [CrossRef]

**26. **S. G. Carter, V. Birkedal, C. S. Wang, L. A. Coldren, A. V. Maslov, D. S. Citrin, and M. S. Sherwin, “Quantum coherence in an optical modulator,” Science **310**, 651–653 (2005). [CrossRef]

**27. **F. D. Morgenthaler, “Velocity modulation of electromagnetic waves,” IRE Trans. Microwave Theory Tech. **6**, 167–172 (1958). [CrossRef]

**28. **D. K. Kalluri, *Electromagnetics of Time Varying Complex Media: Frequency and Polarization Transformer*, 2nd ed. (CRC Press, 2010).

**29. **K. Qu, Q. Jia, M. R. Edwards, and N. J. Fisch, “Theory of electromagnetic wave frequency upconversion in dynamic media,” Phys. Rev. E **98**, 023202 (2018). [CrossRef]

**30. **M. I. Bakunov and S. N. Zhukov, “Conversion of a surface electromagnetic wave at the boundary of a time-varying plasma,” Plasma Phys. Rep. **22**, 649–658 (1996). [CrossRef]

**31. **M. I. Bakunov, A. V. Maslov, and S. N. Zhukov, “Time-dependent scattering of a standing surface plasmon by rapid ionization in a semiconductor,” Opt. Lett. **25**, 926–928 (2000). [CrossRef]

**32. **M. I. Bakunov, A. V. Maslov, and S. N. Zhukov, “Scattering of a surface plasmon polariton by rapid plasma creation in a semiconductor slab,” J. Opt. Soc. Am. B **16**, 1942–1950 (1999). [CrossRef]

**33. **M. I. Bakunov and A. V. Maslov, “Reflection and transmission of electromagnetic waves at a temporal boundary: comment,” Opt. Lett. **39**, 6029 (2014). [CrossRef]

**34. **K. Huang and T. Hong, “Dielectric polarization and electric displacement in polar-molecule reactions,” J. Phys. Chem. A **119**, 8898–8902 (2015). [CrossRef]

**35. **G. A. Menendez and B. Maes, “Time reflection and refraction of graphene plasmons at a temporal discontinuity,” Opt. Lett. **42**, 5006–5009 (2017). [CrossRef]

**36. **J. Wilson, F. Santosa, M. Min, and T. Low, “Temporal control of graphene plasmons,” Phys. Rev. B **98**, 081411 (2018). [CrossRef]

**37. **S. Boubanga-Tombet, S. Chan, T. Watanabe, A. Satou, V. Ryzhii, and T. Otsuji, “Ultrafast carrier dynamics and terahertz emission in optically pumped graphene at room temperature,” Phys. Rev. B **85**, 035443 (2012). [CrossRef]

**38. **T. Li, L. Luo, M. Hupalo, J. Zhang, M. C. Tringides, J. Schmalian, and J. Wang, “Femtosecond population inversion and stimulated emission of dense Dirac fermions in graphene,” Phys. Rev. Lett. **108**, 167401 (2012). [CrossRef]

**39. **Y. V. Bludov, A. Ferreira, N. M. R. Peres, and M. I. Vasilevskiy, “A primer on surface plasmon-polaritons in graphene,” Int. J. Mod. Phys. B **27**, 1341001 (2013). [CrossRef]

**40. **A. V. Maslov, “Levitation and propulsion of a Mie-resonance particle by a surface plasmon,” Opt. Lett. **42**, 3327–3330 (2017). [CrossRef]

**41. **M. I. Bakunov and A. V. Maslov, “Trapping of an electromagnetic wave by the boundary of a time-varying plasma,” Phys. Rev. E **57**, 5978–5987 (1998). [CrossRef]

**42. **M. I. Bakunov and A. V. Maslov, “Trapping of electromagnetic wave by nonstationary plasma layer,” Phys. Rev. Lett. **79**, 4585–4588 (1997). [CrossRef]