Second harmonic generation in graphene dressed by a strong terahertz field
Mikhail Tokman, Sergei B. Bodrov, Yuri A. Sergeev, Alexei I. Korytin,, Ivan Oladyshkin, Yongrui Wang, Alexey Belyanin, Andrei N. Stepanov
TL;DR
This paper reports the observation of enhanced second-harmonic generation in monolayer graphene under an ultra-strong terahertz pulse, revealing a nonperturbative light-matter interaction regime with significant nonlinear optical response.
Contribution
It introduces a novel experimental observation of enhanced second-harmonic generation in graphene driven by strong terahertz fields and develops an analytic theory to explain this effect.
Findings
Enhanced second-harmonic signal observed in graphene
Strong asymmetric carrier distribution causes nonlinear response
Analytic model accurately predicts intensity and polarization
Abstract
We observe enhanced second-harmonic generation in monolayer graphene in the presence of an ultra-strong terahertz field pulse with a peak amplitude of 250 kV/cm. This is a strongly nonperturbative regime of light-matter interaction in which particles get accelerated to energies exceeding the initial Fermi energy of 0.2 eV over a timescale of a few femtoseconds. The second-harmonic current is generated as electrons drift through the region of momenta corresponding to interband transition resonance at an optical frequency. The resulting strongly asymmetric distortion of carrier distribution in momentum space gives rise to an enhanced electric-dipole nonlinear response at the second harmonic. We develop an approximate analytic theory of this effect which accurately predicts observed intensity and polarization of the second-harmonic signal.
Peer Reviews
No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.
Videos
No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.
Second harmonic generation in graphene dressed by a strong terahertz field
Mikhail Tokman
Institute of Applied Physics, Russian Academy of Sciences, Nizhny Novgorod, 603950, Russia
Sergei B. Bodrov
Institute of Applied Physics, Russian Academy of Sciences, Nizhny Novgorod, 603950, Russia
Yuri A. Sergeev
Institute of Applied Physics, Russian Academy of Sciences, Nizhny Novgorod, 603950, Russia
Alexei I. Korytin
Institute of Applied Physics, Russian Academy of Sciences, Nizhny Novgorod, 603950, Russia
Ivan Oladyshkin
Institute of Applied Physics, Russian Academy of Sciences, Nizhny Novgorod, 603950, Russia
Yongrui Wang
Department of Physics and Astronomy, Texas A&M University, College Station, TX, 77843 USA
Alexey Belyanin
Department of Physics and Astronomy, Texas A&M University, College Station, TX, 77843 USA
Andrei N. Stepanov
Institute of Applied Physics, Russian Academy of Sciences, Nizhny Novgorod, 603950, Russia
Abstract
We observe enhanced second-harmonic generation in monolayer graphene in the presence of an ultra-strong terahertz field pulse with a peak amplitude of 250 kV/cm. This is a strongly nonperturbative regime of light-matter interaction in which particles get accelerated to energies exceeding the initial Fermi energy of 0.2 eV over a timescale of a few femtoseconds. The second-harmonic current is generated as electrons drift through the region of momenta corresponding to interband transition resonance at an optical frequency. The resulting strongly asymmetric distortion of carrier distribution in momentum space gives rise to an enhanced electric-dipole nonlinear response at the second harmonic. We develop an approximate analytic theory of this effect which accurately predicts observed intensity and polarization of the second-harmonic signal.
I Introduction
There has been a surge of interest in the nonlinear and quantum optics of graphene Bonaccorso et al. (2010); Kumar et al. (2013); Hong et al. (2013); Gu et al. (2012); Sun et al. (2010); Otsuji et al. (2012); Yao and Belyanin (2012); Tokman et al. (2013, 2016); Cheng et al. (2014a, b, 2017); An et al. (2013, 2014); Avetissian et al. (2012); Entin et al. (2010); Glazov and Ganichev (2014); Glazov (2011); Hendry et al. (2010); Mikhailov and Ziegler (2008); Mikhailov (2016, 2007, 2011, 2014); Yao et al. (2014); Constant et al. (2016); Wang et al. (2016, 2015); Smirnova et al. (2014); Winzer et al. (2013); Li et al. (2012); König-Otto et al. (2016); Dean and van Driel (2009, 2010); Lin et al. (2014); Brun and Pedersen (2015); Bykov et al. (2012); Wu et al. (2012). High speed of carriers in monolayer graphene, , and the resulting large dipole moment of the optical transition , give rise to high magnitudes of the nonlinear optical susceptibilities. A particularly intriguing question raised in recent studies is whether monolayer graphene can support second-order nonlinear processes such as second harmonic or difference frequency generation and other three-wave mixing processes. Since graphene has an in-plane inversion symmetry, the second-order nonlinear response is forbidden in the electric dipole approximation. Of course graphene, like any surface, has anisotropy between in-plane and out-of-plane electron motion, but the corresponding nonlinear response is very weak Dean and van Driel (2009, 2010). A much stronger second-order nonlinearity originates from nonlocal response beyond the electric dipole approximation, i.e. magnetic-dipole and electric-quadrupole response Tokman et al. (2016); Cheng et al. (2017); Glazov and Ganichev (2014); Glazov (2011); Mikhailov (2011); Yao et al. (2014); Constant et al. (2016); Wang et al. (2016); Smirnova et al. (2014). In Lin et al. (2014) second-harmonic generation (SHG) was observed in suspended graphene and attributed to the inversion symmetry breaking due to wrinkles and tears in a monolayer. In a bilayer graphene the inversion symmetry can be broken by applying a voltage bias in transverse direction Brun and Pedersen (2015).
The most obvious method of creating in-plane anisotropy is the anisotropic perturbation of carrier distribution in the k-space by a constant electric field Cheng et al. (2014b); An et al. (2013, 2014); Glazov and Ganichev (2014); Bykov et al. (2012); Wu et al. (2012). This can be also achieved by applying a low-frequency field, in particular at THz frequencies. This SHG mechanism corresponds effectively to a third-order nonlinearity. However, since the THz frequency is much lower than the optical frequencies, it is natural to introduce an effective second-order nonlinear response which depends on the THz field amplitude as a parameter Cheng et al. (2014b). For a weak low-frequency field the perturbation of an originally isotropic distribution is localized mainly near the Fermi surface. In this case the parameter of anisotropy is expressed through a low-frequency Drude current . The resulting nonlinear response obtained by perturbations turns out to be proportional to the magnitude of (see Cheng et al. (2014b)). The magnitude of the SH current scales as
[TABLE]
where , is the relaxation time of a DC current, is chemical potential, is the electric field amplitude of an optical pumping at frequency . Equation (1) is valid when the perturbation of electron distribution at Fermi surface is small.
The experimental results that we obtained below apparently agree with the scaling of Eq. (1). This is however very surprising for ultra-strong fields used in our experiment. Indeed, in recent experiments including our experiment the THz pulses of several hundreds fs duration and field amplitudes were used; see e.g. Tani et al. (2012); Oladyshkin et al. (2017). Such fields distort the original Fermi distribution far beyond small perturbations: an electron acceleration to energies of the order of 0.2 eV (which is a typical Fermi energy for CVD graphene on a glass substrate) happens over a timescale of only a few fs Oladyshkin et al. (2017). In such strong fields, a standard expression for the current using the Drude relaxation time of ps adopted in Cheng et al. (2014b) would yield a current amplitude much higher than the maximum possible value , where is the surface density of carriers! Furthermore, in ultra-strong fields the carrier density gets multiplied by a large factor during the THz pulse due to the electron-hole pair creation Tani et al. (2012); Oladyshkin et al. (2017).
Therefore, to interpret our experiments we need a theory which is not restricted to a standard perturbative model of a weak deformation of the Fermi surface. Within our analytic model the anisotropic deformation of the particle distribution is formed primarily in the vicinity of the interband transition resonance between the optical pump and particle states dressed by a low-frequency field. We consider the case when the energies of resonant particles are much higher than the Fermi energy. A dramatic change in the properties of a quantum system dressed by a strong field is a universal effect; see e.g. Kocharovskaya et al. (1999) and references therein. In our case this effect shows up as a broadening of the interband resonance due to particle acceleration by a strong low-frequency field in the process of an interband transition intiated by an optical field. The broadening of the resonant region in the k-space, , corresponds to the frequency bandwidth , which is equal to the inverse time of Schwinger pair creation in graphene by the field of magnitude Oladyshkin et al. (2017); Vajna et al. (2015); Lewkowicz et al. (2011). Under the action of an ultra-strong field the region of resonant perturbation of carriers in the k-space turns out to be asymmetric with respect to the resonant frequency given by . This asymmetry gives rise to an anisotropic electric-dipole nonlinear response. The case of an ultrafast electron scattering with characteristic time shorter than the Schwinger time is treated numerically in the Appendix.
Within our model the typical lifetime of the particle perturbation in the resonant region is determined by drift in a low frequency field: , i.e. it corresponds to the Schwinger time . Therefore, an increase in the magnitude of a THz field should broaden the applicability region of the perturbation theory with respect to the optical field. For the above parameters of the THz pulse the time is of the order of , which is smaller than, or of the same order as scattering times in graphene. In this case the saturation of the resonant absorption of the optical pump is weak up to the field amplitudes of order . This fact ensures the validity of the scaling despite the use of high-power femtosecond lasers.
Section 2 describes the experiment. Section 3 contains the basic set of equations describing SHG. Section 4 describes an approximate analytical model. In Section 5 the theory is compared with experiment. Derivation of certain formulas used in the analytic theory and numerical simulations for ultrafast scattering times can be found in the Appendix.
II Experiment
The schematic of the experimental setup is shown in Fig. 1. A Ti-Sapphire laser system (Spitfire, Spectra-Physics) generated pulses of energy 0.7 mJ, central wavelength of 0.795 nm, duration 70 fs and repetition rate 700 Hz. Optical radiation was injected through mirror PM3 parallel to the THz field. The diameter of the optical beam was m FWHM, the intensity of the optical beam on the sample was 3 .
THz pulses were generated in the LiNbO3 crystal as in Fülöp et al. (2010) and focused on the sample at an angle of . The diameter of the THz beam on the sample was m FWHM with respect to the field amplitude. The maximum THz electric field amplitude was 250 kV/cm. The maximum value of the P-polarized field was two times smaller than that of the S-polarized field.
The sample was a CVD graphene monolayer on borosilicate glass. Interaction with substrate led to p-doping to the level of Fermi energy eV Tani et al. (2012), which is much smaller than the particle energy of 0.75 eV corresponding to the interband resonance with pumping.
The SH signal for graphene on glass and for glass substrate only is shown in Fig. 2 together with theoretical results. Here p polarization corresponds to the field in the incidence plane, whereas s polarization is orthogonal to this plane. In the notations sss, pss etc. the first index is the polarization of the optical pump, the second index is the polarization of the THz field and the third index is the polarization of SH photons.
Clearly, significant SH signal existed only for SH photons with the same polarization as the THz field, in agreement with theory.
Fig. 3 shows the number of pss-polarized SH photons as a function of the delay time between the fs laser pulse and THz pulse, superimposed on the temporal profile of the THz pulse. Clearly, the SHG signal follows the THz field squared, . Since the number of SH photons is proportional to the current squared, this corresponds to the scaling Eq. (1). Inset to Fig. 3 shows the dependence of from the energy of a laser pulse, which also agrees with Eq. (1).
Note that the SH radiation propagated in the direction of mirror reflection of the incident optical pump, which proves the coherent nature of the SH signal.
III Theoretical model
III.1 Equations for the density matrix and current
Consider monolayer graphene located in the (x,y) plane. Not too far from the Dirac point, the effective Hamiltonian for carriers interacting with an electromagnetic field is Cheng et al. (2017); Glazov and Ganichev (2014); Wang et al. (2016):
[TABLE]
where , , are Pauli matrices, are unit vectors, and are vector and scalar potentials of an external EM field. The basis functions of the Hamiltonian Eq. (3), defined at , are
[TABLE]
Here eigenstates are determined on a unit area for periodic boundary conditions, indices numerate electron and hole states, is an angle between the quasimomentum and -axis. The eigenenergy corresponding to states in Eq. (3) is given by Katsnelson and Kat︠s︡nelʹson (2012).
In the electric dipole approximation an external electric field can be considered uniform, . It can be defined through either vector or scalar potential. The gauge invariance of observables calculated by solving the Schrödinger or master equation for the Hamiltonian (2) has been proven in Wang et al. (2016). Therefore, any EM gauge can be used. However, using a scalar potential leads to a slightly simpler derivation when solving a density matrix equation; see the comparison in Wang et al. (2016). Furthermore, according to Tokman (2009); Zhang et al. (2014), gauge invariance requires that the relaxation operator in the density matrix equation be made dependent on the vector potential, which would be an additional complication. Therefore, we define an external field through the scalar potential, assuming , in Eq. (2).
In the absence of relaxation (which we will take into account later) the Von Neumann equation in the basis of Eq. (3) is
[TABLE]
where
[TABLE]
For a uniform field the matrix element of the interaction Hamiltonian is diagonal with respect to , . Switching to a continuous k-spectrum, Eqs. (4) and (5) yield a closed set of equations for matrix elements . Denoting the quantum coherence as , population difference as , and the sum of populations as , we obtain
[TABLE]
[TABLE]
[TABLE]
where is the interband transition frequency, the Rabi frequency. Equation (8) is separated from the rest of the system and we won’t use it anymore.
Eqs. (6) and (7) make a closed set of equations, which is a version of semiconductor Bloch equations Cheng et al. (2014b); Al-Naib et al. (2015). Note that the “convective” terms in Eqs. (6)-(8) do not originate from some kind of phenomenological assumptions, e.g. an attempt to make it look like a Boltzmann-type equation or quasiclassical equations of motion with a Berry field. These terms rigorously follow from the matrix elements of the interaction Hamiltonian Eq. (5) in a given gauge. If we defined the same field through the vector potential, we would get the equations in a different form, see e.g. Malic et al. (2011). This would not change the observables of course.
The current operator is defined by . The observable current is , where are matrix elements of the current operator. Using and performing integration in k-space, we obtain
[TABLE]
[TABLE]
where g = 4 is spin and valley degeneracy. The first term in square brackets in Eqs. (9),(10) is the intraband current and the second term is an interband current.
III.2 Equations for slow variables
Now we make the following ansatz which separates the slow THz dynamics (subscript 0) from fast optical frequencies:
[TABLE]
[TABLE]
[TABLE]
[TABLE]
It also introduces envelopes of the optical fields at the fundamental and second harmonics.
In addition, we add the relaxation operator in its simplest form to the right-hand side of Eqs. (6), (7). This leads to the replacement in Eq. (6) and the population relaxation term in Eq. (7). Here are relaxation rates of the coherence and population difference and is the equilibrium Fermi distribution.
The resulting equations below contain only these slow-varying envelopes. However, we still keep counter-rotating terms, as shown below.
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
In Eqs. (11)-(15) the contribution of the terms oscillating at 2 is neglected since the SH field is very weak. The perturbation at the second harmonics is described by
[TABLE]
[TABLE]
[TABLE]
In numerical modeling the solution to Eqs. (11) - (15) was substituted into Eqs. (16) - (18), whereas the solution to Eqs. (16) - (18) was substituted into Eqs. (9) - (10) to find the SH current.
IV An approximate analytic solution
IV.1 Equations in rotating wave approximation
Within the rotating wave approximation (RWA) we can solve for the dynamics of carriers in the vicinity of an interband transition at the fundamental frequency of the optical field. The RWA corresponds to the following inequalities:
[TABLE]
In this case one can neglect counter-rotating terms and in terms of Eqs. 11, 12, 13, 14, 15, 16, 17 and 18. Furthermore, if the following conditions are satisfied,
[TABLE]
One can assume that , in Eqs. 11, 12, 13, 14, 15, 16, 17 and 18 and obtain approximate equations,
[TABLE]
[TABLE]
[TABLE]
where .
IV.2 The stationary phase solution of density matrix equations
For simplicity we consider the case when optical and THz fields are polarized along x; other polarizations are considered in Sec. IV D. Using and , the formal solution of Eq. (20) can be written as
[TABLE]
where
[TABLE]
The functions , and are obtained from functions , and as
[TABLE]
[TABLE]
[TABLE]
It is easy to check that Eqs. 23 and 24 is an exact solution to Eq. (20), satisfying the initial condition , where is the time when the optical field is turned on.
The integral in Eq. (23) can be evaluated by the method of stationary phase. Within this approach one assumes that the main contribution to the integral comes from the vicinity of a stationary point given by the condition
[TABLE]
Expanding near the stationary point , we obtain
[TABLE]
where ,
[TABLE]
The interval in Eq. (26) which makes the main contribution to the integral is given by
[TABLE]
The following hierarchy of timescales ensures the validity of the stationary phase method:
[TABLE]
Eq. (26) can be written in the form (see Appendix A)
[TABLE]
where
[TABLE]
, . The top and bottom signs in Eq. (30) correspond to and respectively.
The function is normalized as . In the limit (i.e. when ) Eqs. (29), (30) give
[TABLE]
where indicates a principal value of the integral, .
One can see from the expression in Eq. (30) that the scattering-induced broadening of the resonance is replaced by the nonlinear field-induced broadening: . The absorption line described by Eq. (30) is asymmetric, namely it is shifted towards for and towards for . An asymmetric shape is due to the drift of carriers in the k-space in the presence of a THz field.
The functions and in Eq. (29) are Lagrangian variables and , determined at the current moment of time t shifted with respect to a stationary point given by Eq. (25). The shift amount is , see Appendix. The possibility to use a “local” approximation and is discussed in Appendix A. It turns out that the first and the second inequality in Eqs. (28a) always ensure the validity the local approximation for the Rabi frequency. At the same time, the local approximation for the population difference requires in addition that the perturbation of the initial value in the resonant region be small enough. The perturbation of the population difference can be estimated using Eqs. (21), (29). The analysis in Appendix A yields
[TABLE]
Equation (31) looks like a standard result for a two-level system coupled to an optical field if we replace the population relaxation time with in the density matrix equations. This is to be expected, since the lifetime of the perturbation in the presence of a strong low-frequency field is determined by the drift of carriers out of the resonant region: . Using , and longer than 5 fs (i.e. for smaller than 300 kV/cm) one can obtain for . Here we also assumed , which is the case for high enough electron energies corresponding to the interband resonance .
IV.3 The nonlinear current
For the electric field polarized along x, Eqs. (9) and (22) give the following expression for the complex amplitude of the current at SH:
[TABLE]
Substituting Eq. (29) into Eq. (32) and assuming and leads to
[TABLE]
where is defined in Eq. (30). Limiting ourselves to the resonant region, we can expand when evaluating the inegral in Eq. (33) (here ) and change from to infinite integration limits when integrating over :
[TABLE]
The first two terms in the expansion in Eq. (34) give zero after angle integration . For the first term this is obvious due to the normalization condition ; the proof for the second term is in the Appendix. To calculate the integral for the third term in the expansion we take into account that the function is a linear combination of the even function and odd function , see Eq. (30). As a result, we obtain
[TABLE]
where . Recalling that the upper and lower signs correspond to and , we get
[TABLE]
IV.4 Polarization dependence
Eq. (35) has been obtained for collinear vectors . Consider now a different orientation when . Let and , but . The expression for the Rabi frequency takes the form . It is straightforward to obtain that the relationship between the SH current and has the same form as in Eq. (32). Therefore, instead of Eq. (35) we obtain
[TABLE]
Comparing Eqs. (35) with (36) one can see that the power of the SHG signal at is 9 times higher than at for the same magnitudes of the fields.
In the case we always get for any polarization of the optical field. This is expected, since for this geometry the low-frequency field cannot break the inversion symmetry along the direction of see the Appendix.
V Comparison with experiment
The field dependence and polarization dependence of the observed SH signal coincides with those predicted by the model. To compare absolute numbers of SH photons, one needs to take into account that
(i) Eqs. (35), (36) contain the components of the fields tangential to the monolayer. Using Fresnel formulas Landau and Lifshitz (1984), one can get
[TABLE]
Here is the amplitude of the transverse field of an incident S- or P-polarized wave, is its tangential component on the monolayer, the dielectric constant of the substrate, the incidence angle from air with respect to the normal to the surface of graphene.
(ii) The surface current at frequency and the amplitude of the S- or P-polarized field radiated by this current are related by
[TABLE]
where * *is the reflection angle.
(iii) The powers of the S- and P-polarized THz radiation in the experiment differ by a factor of 4.
For numerical estimates we assume , , .
Next we calculate the number of photons in the PSS polarization configuration. From Eqs. (36) (38) for the above parameters , , and we can estimate the field of the SH signal as
[TABLE]
which is equal to . For the pulse duration and beam cross-section used in the experiment, and for a 7% experimental efficiency of the detection system, the SH field in Eq. (39) corresponds to about 850 SH photons per series of 60,000 laser pulses. This agrees with an experimentally measured SH signal within the experimental accuracy.
Relative SHG efficiency for other polarizations can be obtained from Eqs. (35)(38):
[TABLE]
where
[TABLE]
For the same parameters we get , , which gives , , .
These theoretical values are compared with experiment in Fig. 2. As is clear from the figure, there is a good agreement between theory and experiment.
Acknowledgements.
This work has been supported by the RFBR grant No. 18-29-19091. Y.W. and A.B. acknowledge the support by Air Force Office for Scientific Research through Grants No. FA9550-17-1-0341 and FA9550-14-1-0376.
Appendix A The stationary phase solution for the density matrix equation
According to the first inequality in Eq. (28a) the optical pulse is longer than the integration interval in Eq. (26). Therefore the lower limit of integration can be taken as ; therefore, Eq. (26) yields
[TABLE]
where the upper and lower signs in correspond to and respectively. Furthermore, from Eqs. (25), (24) one can get
[TABLE]
Since the THz pulse is much longer than the optical pulse, one can write the integral as
[TABLE]
where can be treated as a constant during the optical pulse.
The calculations are greatly simplified if we assume
[TABLE]
we will check its validity later.
Substituting Eq. (43) into Eq. (42) and using Eq. (44), we obtain:
[TABLE]
Using Eq. (45) in Eqs. (24), (27), gives
[TABLE]
Substituting Eqs. (46), (47) into Eq. (41) results in
[TABLE]
where the upper and lower signs in are for and respectively. Taking into account that the integration interval which makes the main contribution near the resonance is and using Eq. 45, we obtain
[TABLE]
which means that the last inequality in Eq. 28b ensures the validity of the approximation Eq. 44. In the region the factor in Appendix A cannot be greater than , so if the last inequality in Eq. 28a is satisfied, it can be taken as 1. As a result, we obtain Eqs. (29), (30). We also provide here useful asymptotics of the function in Eq. (30) at :
[TABLE]
where the upper and lower signs are taken for and respectively.
Appendix B Approximate expressions for and
The functions can be written for the same conditions Eq. 44 and Eq. 45 as
[TABLE]
Taking into account that , and , we can drop the terms with time derivative in Eq. 52 and Eq. 53. The factor is small too. At the same time the term in Eq. 53 is small only if the perturbation of the population difference in the resonant region is small too. Therefore, we obtain
[TABLE]
Appendix C Quasistationary perturbation of populations
At high carrier energies one can take in Eq. (21). Consider a stationary solution of this equation:
[TABLE]
where the notation means . The boundary C of the integration limit in Eq. 56 is chosen where the effective source approaches zero. This choice depends also on the sign of , which determines the direction of particle drift in the k-space: for we get , whereas for we get . Next, we substitute expressions Eqs. (29), (30) into Eq. 56 and take into account that the characteristic “size” of the source in the k-space . As a result, under the condition we arrive at Eq.(31) for the perturbation of populations. Note that the inequality is equivalent to .
Appendix D The second term in the nonlinear current expansion
Here we prove that the second term in the expansion in Eq.(34) is equal to zero. Its expression is given by
[TABLE]
where and is a constant. We introduced the factor which makes the proof easier, and we will take the limit in the end. Taking into account that the function is a linear combination of the even function and odd function (see Eq. (30)), we obtain
[TABLE]
where
[TABLE]
In the last expression upper and lower signs are given by the signs of , as usual. Therefore we get
[TABLE]
Appendix E Polarization selection rules for THz field-induced SHG
Consider the orientation , and . In this case one should take , and in Eqs. (20)(22). The dependence of on has the same form as in Eq. (31), whereas in Eq. (29) one has to replace with . Upper and lower signs in all coefficients in Eq. (29) correspond to the signs of . As a result, instead of Eq.(35) we get
[TABLE]
If , and , similar considerations lead to
[TABLE]
Appendix F Numerical simulations
We used Eqs. 11, 12, 13, 14 and 15 to simulate the SHG in graphene illuminated with a strong THz pulse and an optical field beyond the stationary phase approximation. To derive these equations, the amplitudes of optical-frequency coherences and population differences and were assumed to be slow-varying. When the optical field is far off-resonance from an interband transition, Rabi oscillations can have a frequency comparable to the frequency detuning of the optical field. However, the slow-varying assumption is still valid if Rabi oscillations are strongly damped by ultrafast dephasing processes. So, we use ultrafast population relaxation and dephasing times for hot photoexcited electrons, fs and fs, which is consistent in order of magnitude with results from related studies Winzer et al. (2010); Xing et al. (2010); Tan et al. (2017) and can be attributed to strong Coulomb interaction between carriers in graphene which results in ultrafast carrier-carrier scattering through interband and intraband Auger recombination and impact ionization Winzer et al. (2010); Xing et al. (2010); Tan et al. (2017); Tani et al. (2012); Huang et al. (2018). (Auger recombination could be enhanced due to lattice imperfections in CVD graphene Tani et al. (2012).) Then we can still assume that the optical-frequency populations and coherences follow the source terms adiabatically, namely, we can put the to be zero in all equations except those for . Also, for reasons already discussed above we can assume the optical field to be weak enough to treat it in a perturbative way.
Another technical difficulty is that has a singularity at , which can lead to divergence in numerical simulations. To avoid this problem, we replaced in the denominator of by . We also assumed the chemical potential meV and electron temperature at equilibrium K. The THz field is chosen to be polarized in -direction, and the optical field is polarized in -direction.
The simulation shows that the SH current is generated predominantly in -direction, i.e. along the direction of the THz field. In Fig. 4 we plot the SH current calculated from the simulation, together with the profile of the THz pulse. The SH current generally follows the THz field, except for some small variations originated from the time evolution of the carrier distribution . These variations are likely beyond the detector resolution in the experiment.
Dressing of electron states by the THz field is expected to play an important role in the SHG process, even in the presence of ultrafast scattering. As an illustration, we assume that the carrier distribution is in thermal equilibrium, and calculate the SH current for two cases, with and without terms included. The latter would be similar to four-wave mixing in a two-level medium. Figure 5 shows the dependence of the SH current on the THz field amplituce for the two cases. It indicates that the dressed-state effect can enhance the SH current by one order of magnitude, and therefore the signal intensity by two orders. We can also see that the SH current at the highest THz field of 250 kV/cm is about three times higher in Fig. 5 as compared to Fig. 4 where the thermal equilibrium distribution for was assumed. One could say that roughly 1/3 of the SH current comes from direct parametric interaction between a THz field and an optical field, whereas the nonequilibrium distortion of carrier distribution contributes the remaining 2/3.
The reference list from the paper itself. Each links out to its DOI / PubMed record.
- 1Bonaccorso et al. (2010) F. Bonaccorso, Z. Sun, T. Hasan, and A. Ferrari, Nature photonics 4 , 611 (2010).
- 2Kumar et al. (2013) N. Kumar, J. Kumar, C. Gerstenkorn, R. Wang, H.-Y. Chiu, A. L. Smirl, and H. Zhao, Phys. Rev. B 87 , 121406 (2013) . · doi ↗
- 3Hong et al. (2013) S.-Y. Hong, J. I. Dadap, N. Petrone, P.-C. Yeh, J. Hone, and R. M. Osgood, Phys. Rev. X 3 , 021014 (2013) . · doi ↗
- 4Gu et al. (2012) T. Gu, N. Petrone, J. F. Mc Millan, A. van der Zande, M. Yu, G.-Q. Lo, D.-L. Kwong, J. Hone, and C. W. Wong, Nature Photonics 6 , 554 (2012).
- 5Sun et al. (2010) D. Sun, C. Divin, J. Rioux, J. E. Sipe, C. Berger, W. A. de Heer, P. N. First, and T. B. Norris, Nano Letters 10 , 1293 (2010) , p MID: 20210362, https://doi.org/10.1021/nl 9040737 . · doi ↗
- 6Otsuji et al. (2012) T. Otsuji, S. A. B. Tombet, A. Satou, H. Fukidome, M. Suemitsu, E. Sano, V. Popov, M. Ryzhii, and V. Ryzhii, Journal of Physics D: Applied Physics 45 , 303001 (2012) .
- 7Yao and Belyanin (2012) X. Yao and A. Belyanin, Phys. Rev. Lett. 108 , 255503 (2012) . · doi ↗
- 8Tokman et al. (2013) M. Tokman, X. Yao, and A. Belyanin, Phys. Rev. Lett. 110 , 077404 (2013) . · doi ↗