Temporal scattering of a graphene plasmon by a rapid carrier density decrease

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 $ϵ$; see Fig. 1(a). In all numerical results we set $ϵ=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

 

Fig. 1. (a) Illustration of a surface plasmon propagating along a graphene sheet at $t<0$. (b) Time dependence of the graphene carrier density. (c) Dispersion diagram showing the frequency transformation of the initial plasmon when the carrier density decreases from $N_1$ to $N_2$. The lines labeled by 1 and 2 are the dispersion curves for the plasmons at densities $N_1$ and $N_2$, respectively. The shaded region shows the continuous spectrum for bulk waves; the hatched region defines the spectrum of the generated transient radiation. The dashed lines $h=\sqrt{ϵ}|\omega |/c$ are the light lines.

Download Full Size | PPT Slide | PDF

 

Fig. 2. (a) Phase $n_{\mathrm{ph}}$ and group $n_{\mathrm{gr}}$ refractive indices for graphene plasmons. The dashed line shows the asymptotics $n_{\mathrm{ph}}\approx \omega /(2\pi \mathrm{\Omega })$. (b) Decomposition of the plasmon energy $W_0$ into the electric $W_E$, magnetic $W_H$, and kinetic $W_k$ parts.

Download Full Size | PPT Slide | PDF

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

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).

 

Fig. 3. Frequency ${\omega }_2$ of the excited plasmon relative (a) to the initial frequency ${\omega }_1$ and (b) to the frequency parameter ${\mathrm{\Omega }}_2$ as a function of ${\mathrm{\Omega }}_2/{\mathrm{\Omega }}_1$.

Download Full Size | PPT Slide | PDF

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 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$.

 

Fig. 4. Energy distribution after temporal scattering as a function of ${\mathrm{\Omega }}_2/{\mathrm{\Omega }}_1$ for (a) ${\omega }_1/{\mathrm{\Omega }}_1=15$ and (b) ${\omega }_1/{\mathrm{\Omega }}_1=1$.

Download Full Size | PPT Slide | PDF

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.

 

Fig. 5. Far-field angular distribution $w_b(\theta )/W_0$ for ${\mathrm{\Omega }}_2/{\mathrm{\Omega }}_1=0.5$ and several values of ${\omega }_1/{\mathrm{\Omega }}_1$.

Download Full Size | PPT Slide | PDF

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 }/(ϵ\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

 

Fig. 6. Comparison of the plasmon transmission $t_H$ and reflection $r_H$ coefficients obtained here (general and quasi-static results) and available from the literature for ${\omega }_1/{\mathrm{\Omega }}_1=15$. Equation (11) is the general formula and Eq. (18) is the quasi-static formula derived here. Equation (19) is from Ref. [35]. Equation (20) is from Ref. [36].

Download Full Size | PPT Slide | PDF

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$.

 

Fig. 7. Snapshots of the electric field distribution $E_x(x,y,t)$ produced by a plasmon wavepacket with central frequency ${\omega }_1/{\mathrm{\Omega }}_1=10$ at three time moments $t=0$, $0.15/{\mathrm{\Omega }}_1$, $7/{\mathrm{\Omega }}_1$. Initially the wavepacket propagates along the graphene sheet at $y=0$ in the $+x$ direction (frame $t=0$). Soon after the carrier density jump ${\mathrm{\Omega }}_2/{\mathrm{\Omega }}_1=0.2$ at $t=0$, transient processes take place (frame $t=0.15/{\mathrm{\Omega }}_1$). Finally, two propagating plasmons are formed (frame $t=7/{\mathrm{\Omega }}_1$).

Download Full Size | PPT Slide | PDF

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.

 

Fig. 8. Obtaining the transmission coefficients $t_{E,H}$ from the FDTD simulations. (a) Absolute value of the incident $E_x^i(\omega )$, $H_z^i(\omega )$ and transmitted $E_x^t(\omega )$, $H_z^t(\omega )$ field spectra for ${\mathrm{\Omega }}_2/{\mathrm{\Omega }}_1=0.6$. (b) Analytically calculated frequency ${\omega }_2$ of the excited plasmon as a function of the frequency ${\omega }_1$ of the incident plasmon. (c) Differential spectral transormation $f=d{\omega }_2/d{\omega }_1$, which is the derivative of ${\omega }_2({\omega }_1)$ shown in frame (b). (d) Comparison of the transmission coefficients obtained numerically from the FDTD simulations and analytically using Eqs. (11)–(13).

Download Full Size | PPT Slide | PDF

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.