Nonlinear optical response of graphene in terahertz and near-infrared frequency regime

Abstract

In this review, we discuss our recent theoretical work on the nonlinear optical response of graphene and its sister structure in terahertz (THz) and near-infrared frequency regime. Due to Dirac-like linear energy-momentum dispersion, the third-order nonlinear current in graphene is much stronger than that in conventional semiconductors. The nonlinear current grows rapidly with increasing temperature and decreasing frequency. The third-order nonlinear current can be as strong as the linear current under moderate electric field strength of 10⁴ V/cm. In bilayer graphene (BLG) with low energy trigonal warping effect, not only the optical response is strongly nonlinear, the optical nonlinearity is well-preserved at elevated temperature. In the presence of a bandgap (such as semihydrogenated graphene (SHG)), there exists two well separated linear response and nonlinear response peaks. This suggests that SHG can have a unique potential as a two-color nonlinear material in the THz frequency regime where the relative intensity of the two colors can be tuned with the electric field. In a graphene superlattice structure of Kronig-Penney type periodic potential, the Dirac cone is elliptically deformed. We found that not only the optical nonlinearity is preserved in such a system, the total optical response is further enhanced by a factor proportional to the band anisotropy. This suggests that graphene superlattice is another potential candidate in THz device application.

Key words： graphene; terahertz (THz) response; nonlinear effect; photomixing

Introduction

Physical properties of graphene

Graphene is a one-atom thick, 2-dimensional honeycomb structure made up of carbon atoms (Fig. 1). The single layer was first isolated and systematically studied by Novoselov et al. in 2005 [ 1, 2]. The first theoretical study of graphene however dated back to 1947. In his pioneering work [ 3], Wallace presented that graphene is a gapless semiconductor, whose valence band touches the conduction band at K and K′ points of its Brillouin zone (Fig. 2). The most fascinating aspect of this touching point, or the ‘Dirac point’, is that the energy band is in a linear form of E_k = ±ħv_Fk and hence the electrons around the Dirac points behave like a massless ultra-relativistic fermions, but moving with a much reduced ‘speed of light’v_F≈c/300 (c = vacuum speed of light) [ 1]. This aspect is fundamentally different from the Schrodinger fermions E_k = ħ²k²/2m* in conventional semiconductors. Many unusual physical phenomena arise due to the relativistic quasiparticle dynamics. For example, the K point electrons exhibit anomalous perfect tunneling effect despite the potential barrier height and width [ 4]. This anomalous tunneling behavior is related to the Klein tunneling of massless spin-1/2 Dirac particles in quantum electrodynamics [ 5]. The Klein tunneling was thought to be a textbook example to illustrate the bizarre consequence of Dirac equation. The scale-down condensed matter version can now be realized in graphene [ 6, 7]. Moreover, geometrical asymmetry and applied magnetic field can result in strong enhancement of optical and magnetic properties [ 8- 12].

Because of the relativistic dynamics, electron scattering in graphene is strongly suppressed, and this results in unusually high electron mobility [ 13, 14]. It is believed that electron mobility of 100000 cm²/Vs can be achieved in high quality sample [ 15], suggesting a promising transistor application [ 16, 17]. The massless Dirac fermion exhibits another unusual behavior in the presence of a magnetic field. The quantum Hall conductivity follows half-integer steps {\sigma }_{xy}=\pm \mathrm{4}{e}^{\mathrm{2}}/h(N+\mathrm{1}/\mathrm{2}) , where h is the Planck constant, [ 1, 18, 19] instead of the conventional integer-multiple quantized conductivity in conventional semiconductor. This is again due to the relativistic spectrum of the K electrons where the E = 0 Landau level is shared by both electrons and holes. Furthermore, the large energy separation between the two lowest Landau levels allows the quantum Hall effect to survive even at room temperature [ 20]. The existence of a minimal direct current conductivity in the absence of charge carrier is another surprising result. The minimal conductivity has a well-established experimental value of {\sigma }_{\min }=\mathrm{4}{e}^{\mathrm{2}}/h [ 1] yet its physical origin is not well-understood since various theoretical models [ 4, 21- 24] yield very different values of σ_min. It has been suggested that the many-body interactions, wrinkling and ripples of the graphene sheet and the formation of electron-hole puddles [ 25] could be the possible underlying mechanisms. When interact with photon, the massless Dirac fermion manifests itself as a universal interband optical conductivity of e²/4ħ [ 26- 31].

Fig.1 Graphene, an atomically thin layer of carbon atoms arranged in honeycomb structure, a₁ and a₂ are the lattice unit vectors. Red and green dots are atoms from two sublattices

Full size|PPT slide

Apart from the unusual electronic and optical properties, electrons in graphene can be manipulated in completely different ways. In graphene, the electron not only has spin, but also possesses two additional degrees of freedom: valley and pseudospin. Such additional degrees of freedom arise from the fact that the low energy electron resides in two K and K′ valleys and their relativistic nature is described by a two-component pseudospinor wave function. Although still in the early conceptual stages, the valley and pseusospin degree of freedoms in graphene open up the possibilities of ‘valleytroics’ and ‘pseudospintronics’ devices. The concept of ‘valleytronics’ was first proposed in a nano-constricted device in which two graphene sheets are connected by a narrow zigzag-edge nanoribbon [ 32]. The electron transport across the junction becomes valley-dependent and the degree of valley polarization is tunable by a gate voltage. Many other strategies have since been proposed. For example, valley-dependent scattering by a line defect [ 33], spatial splitting of valley current by the trigonally warped band structure at high energy regime [ 34, 35], tunneling barriers based on gapped graphene and strain-engineered grapheme [ 36- 38], and the valley-dependent focusing and de-focusing effect in bilayer graphene (BLG) n-p junction [ 39] can all be utilized to produce valley polarization. In addition to the valley degree of freedom, pseudospin magnetization can be generated in graphene with a bandgap [ 40, 41] or spontaneously generated via electron-electron interaction [ 42, 43]. More importantly, the pseudospin magnetization can be optically probed [ 42]. Pseudospin valves in a graphene/superconductor/graphene heterostructure and in a BLG electrostatic tunneling barrier [ 41, 44] offer further possibilities to manipulate the transport of the pseudospins.

Fig.2 (a) Reciprocal lattice of graphene, b₁ = a₂ × e_z/A, b₂ = a₁ × e_z/A, A = (a₁ × a₂) ∙ e_z; (b) band structure near the Dirac point

Full size|PPT slide

In terms of device application, graphene is a ‘designer’ structure whose electronic properties can be tailor-made to meet any device requirement. Graphene can be cut into ribbons or be transformed into superlattices via electrostatic gating. The electronic properties of graphene nanoribbon can be tuned by varying the width and the type of its edges into armchair or zigzag configurations [ 45- 47]. For instance, a bandgap can be opened in armchair graphene nanoribbon and the size of the bandgap is tunable via the nanoribbon width [ 45, 46]. In the case of graphene superlattices [ 48- 50], elliptical deformation of the Dirac cone can be engineered without breaking the k-linearity of the band structure.

Although it has only been 10 years since the first isolation of graphene [ 2], myriads of unusual properties, such as the strong suppression of weak localization [ 51- 53], thermoelectric transport [ 54- 56], quantum spin Hall effect [ 57], chiral superconductivity [ 58], just to name a few, have been discovered and many more are still continually emerging. It is therefore impossible to fully cover all aspects of graphene in this brief overview. Broader discussions of graphene can be found in several classic review articles [ 15, 59- 62]. Finally, we remark that in a new class of materials, i.e. the topological insulator (TI), the surface states can also described by a Dirac cone. Although it is not the scope of this paper to discuss the physics of TI, it is worth-noting that many of the unusual properties of graphene can be directly translated into TI [ 63]. Together with the emerging single layer honeycomb structures of group IV atoms such as silicene [ 64, 65], germanene [ 66], and stanene [ 67] (single layer silicon, germanium and tin, respectively), it is not unreasonable to speculate that the physics of Dirac fermions shall play an important role in the upcoming developments of condensed matter physics.

Nonlinear optical properties of graphene

We theoretically study the nonlinear optical properties of graphene and its sister-structures in terahertz (THz) and far-infrared (FIR) frequencies [ 68- 71]. The motivation behind these studies arises from two factors. THz wave plays important role in the study of condensed matter since many dynamical processes occur in the THz frequency regime (approximately a few meV). THz is also an invaluable tool in the field of astrophysics, telecommunication, non-destructive imaging and chemical/bio-molecules identifications [ 72]. Unfortunately, THz frequency situated right in between the optics and electronic regimes. Efficient generation and detection of THz waves are problematic because it is too high of a frequency via electronic approach and too low of a frequency via photonics approach. The hunt for an efficient mean of THz generation and detection is therefore one of the ongoing primary objectives. Second, exceptionally strong optical response has been reported in graphene both theoretically and experimentally [ 73- 77]. The third-order nonlinear susceptibility χ⁽³⁾ in graphene is 10⁸ stronger than that in a bulk insulator [ 73- 75]. Furthermore, Wright et al. [ 77] has found that the THz/FIR interband optical conductivity can be significantly enhanced by 3-photon nonlinear interband optical processes under electric field strength in the order of 10³ V/cm. The rather weak 2.3% absorption (corresponding to the universal conductivity e²/4ħ) can hence be overcome by the nonlinear optical absorption. Although not directly observed in free standing single layer graphene in THz range, giant nonlinear transmittance has experimentally been observed in graphene dispersions [ 78, 79] and, recently, the third-harmonic generation in graphene on a substrate in near-infrared frequency has been experimentally demonstrated [ 80]. In BLG, second harmonic can be generated by breaking the symmetry using an in-plane electric field [ 81]. Although the 0.2 eV photon energy is well-beyond the THz regime, the unusually large χ⁽²⁾≈ 10⁵ pm/V highlights the potential of BLG in nonlinear photonics application.

The optical nonlinearity in graphene is directly related to its linear energy spectrum. The energy dispersion of the massless Dirac fermion around the K point is written in the form of ϵ_s = sħv_F|k| and the group velocity is v_s = sv_F \widehat{\mathit{k}} where \widehat{\mathit{k}} is the unit vector of the wavevector k. The group velocity is completely independent of wavevector k. From a pedagogical point of view, the Dirac fermions are expected to oscillate abruptly between the two values of +v_F and -v_F when driven by an external oscillating electric field. This gives rise to a series of square-wave optical response. Since a square function is rich in higher-order harmonics, the massless Dirac fermion is expected to exhibit strong nonlinear optical response. This is in contrast to the Schrodinger electrons of ϵ_k = ħ²k²/2m and group velocity v = ħk/m. The k-dependent group velocity allows the optical current response to oscillate continuously with the external electric field and hence the anharmonicity is absent. Although Ishikawa has shown that the highly anharmonic intraband current response (i.e., the square current response as discussed above) is reduced by a interband component [ 82], nonlinear optical response such as frequency up-conversion is still expected to be a significant optical process in graphene.

Nonlinear optical response of graphene and its related structures

Terahertz photon-mixing effect of gapless and gapped single layer graphene

The nonlinear intraband optical response of gapless graphene has been previously studied by Mikhailov et al. using the semiclassical electron transport equation for two limiting cases: (i) zero doping at finite temperature; and (ii) finite doping at zero temperature [ 74, 75]. The intermediate regime between (i) and (ii), i.e., doped graphene at finite temperature, is however left open and has not been reported so far. The nonlinear response in this intermediate regime is important since finite doping is usually present due to crystal imperfection and impurities, and the practical implementation of graphene-based device requires finite temperature information. Furthermore, nonlinear response usually occurs under strong external field. The strong-field-drive Dirac fermion (SDF) population redistribution due their externally perturbed dynamics and non-equilibrium carrier heating becomes inevitable in strong-field regime. Optical response of graphene with these strong-field effects considered has however not been reported. In this section, we fill in these gaps by constructing the full temperature spectrum of the nonlinear optical response of a finite-doped ( \mu \ne \mathrm{0} ) graphene single layer in both gapless and gapped cases under both weak-field and strongfield conditions. The dynamics of the quasiparticles when perturbed by a strong electric field are decomposed into linear and nonlinear components, and the optical nonlinearity of the graphene is investigated.

The effective Hamiltonian of single layer graphene around the K point is given as

(1) \widehat{\mathit{H}}={\mathit{v}}_F\left[\begin{array}{ll}\mathrm{0}& p_{-}\\ p_{+}& \mathrm{0}\end{array}\right],

where the Fermi velocity is v_F = 3ta/2ħ≈ 10⁶ m/s, t≈ 3 eV is the nearest neighbor hopping bandwidth, a≈ 0.142 nm is the carbon-carbon distance, and p_± = p_x±ip_y. The energy eigenvalue of Eq. (1) gives rise to the linear energy dispersion ϵ_s = sv_Fp, where s = ± 1. This energy dispersion results in a symmetric upper (s = + 1) and lower (s = - 1) Dirac cones, representing electrons and hole states respectively, and is analog to the charge conjugation symmetry in quantum electrodynamics. The velocity operator is given by \widehat{v}=\partial \widehat{\mathit{H}}/\partial p . Following Feynman [ 83], we write the expectation value of \widehat{v} as \langle {\widehat{v}}_s\rangle =\partial {ϵ}_s/\partial p .

This gives velocity eigenvector v_s = sv_Fp/p. We consider a time-dependent applied electric field in the form of

(2) \mathit{E}(\mathit{r},t)={\displaystyle \sum _{\mu }{\mathit{E}}_{\mu }\exp \{i({\mathit{q}}_{\mu }\cdot \mathit{r}-{\omega }_{\mu }t)\}},

where {\mathit{E}}_{\mu } , {\mathit{q}}_{\mu } and {\omega }_{\mu } are the amplitude, wavevector and frequency of the \mu -th wave of the electric field. Ignoring the weak magnetic component, the external field is minimally coupled to the quasiparticle by performing the substitution p→p-eA(r, t), where E(r,t) = -∂A(r, t)/∂t and e is the electric charge. The velocity becomes

(3) {\mathit{v}}_s=s{v}_F{\displaystyle \frac{\mathit{p}+e\mathit{A}}{\left|\mathit{p}+e\mathit{A}\right|}},

and for simplicity, we denote u= eA. We perform a Taylor expansion on v_s in terms of the externally applied electric field. Assuming that p»u where u = -eA(r, t), we obtain [ 68]

(4) {\mathit{v}}_s^{(\mathrm{0})}=s{v}_F\left({\displaystyle \frac{\mathit{p}}{p}}\right),

(5) {\mathit{v}}_s^{(\mathrm{1})}=s{v}_F\left[{\displaystyle \frac{\mathit{u}}{p}}-{\displaystyle \frac{\mathit{p}}{p}}\left({\displaystyle \frac{\mathit{p}\cdot \mathit{u}}{p^{\mathrm{2}}}}\right)\right],

(6) {\mathit{v}}_s^{(\mathrm{2})}=s{v}_F\left[{\displaystyle \frac{\mathit{p}}{\mathrm{2}p}}{\left({\displaystyle \frac{\mathit{u}}{p}}\right)}^{\mathrm{2}}+{\displaystyle \frac{\mathit{u}}{p}}\left({\displaystyle \frac{\mathit{p}\cdot \mathit{u}}{p^{\mathrm{2}}}}\right)\right],

(7) {\mathit{v}}_s^{(\mathrm{3})}=s{v}_F\left[-{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}{\displaystyle \frac{\mathit{u}}{p}}{\left({\displaystyle \frac{\mathit{u}}{p}}\right)}^{\mathrm{2}}+{\displaystyle \frac{\mathrm{2}}{}}{\displaystyle \frac{\mathit{p}}{p}}{\left({\displaystyle \frac{\mathit{u}}{p}}\right)}^{\mathrm{2}}\left({\displaystyle \frac{\mathit{p}\cdot \mathit{u}}{p^{\mathrm{2}}}}\right)+{\displaystyle \frac{\mathrm{2}}{}}{\displaystyle \frac{\mathit{u}}{p}}{\left({\displaystyle \frac{\mathit{p}\cdot \mathit{u}}{p^{\mathrm{2}}}}\right)}^{\mathrm{2}}-{\displaystyle \frac{\mathrm{2}}{}}{\displaystyle \frac{\mathit{p}}{p}}{\left({\displaystyle \frac{\mathit{p}\cdot \mathit{u}}{p^{\mathrm{2}}}}\right)}^{\mathrm{3}}\right],

where {\mathit{v}}_s^{(i)} represents i-th order velocity of graphene per spin and per valley degeneracy. The zero-order velocity is equal to the Fermi velocity v_F which is consistent with the unperturbed case. Note that the velocities only reverse their directions for between the two Dirac cones of s = ± 1. The magnitude remains unchanged due to the particle-hole symmetry of the energy band structure.

The i-th order current is given by [ 68, 84]

(8) {\mathit{J}}^{(i)}=e{\displaystyle \sum _s{\displaystyle {\int }_{\mathrm{0}}^{\mathrm{2}\mathrm{\pi }}{\displaystyle {\int }_{\mu -{ϵ}_{\mathrm{p}\mathrm{h}}-k_{B}T}^{\mathrm{\Lambda }}{\mathrm{d}}^{\mathrm{2}}p{\mathit{v}}_s^{(i)}f({ϵ}_s)}}},

where ϵ_ph is the energy of the incoming photons, k_B is the Boltzmann constant, T is the temperature, and f (ϵ_s) is the Fermi-Dirac distribution function. The integration cut-off Λ is equal to the Fermi level µ at T = 0 K, and is arbitrarily set to a large value of 0.5 eV for T > 0 K and µ > 0 numerical calculation. Up to room temperature, the Fermi-Dirac distribution terminates the momentum integration well before Λ and hence our choice of Λ is well-justified. For µ < 0, Λ cut off the momentum integration at µ + k_BT to avoid the low momentum regime where p»u fails. Deep charge carriers cannot respond to the external perturbation due to the unavailability of higher energy states. We qualitatively approximate this by choosing a lower momentum integration limit of µ-ϵ_s-k_BT.

Linear optical response

The linear current density for µ > ϵ_ph at T = 0 K, per spin and per valley, is given by

(9) {\mathit{J}}_{T=\mathrm{0}}^{(\mathrm{1})}=-{\displaystyle \frac{i{e}^{\mathrm{2}}}{\mathrm{4}\mathrm{\pi }\mathit{ћ}}}{\displaystyle \sum _{\mu }{\mathit{E}}_{\mu }{e}^{i({\mathit{q}}_{\mu }\cdot \mathit{r}-\omega t)}}.

Including spin and valley degeneracy, the zero temperature linear current density is given as

(10) {\mathit{J}}_{T=\mathrm{0}}^{(\mathrm{1})}=-i{\displaystyle \frac{e^{\mathrm{2}}}{\mathrm{\pi }\mathit{ћ}}}{\displaystyle \sum _{\mu }{\mathit{E}}_{\mu }\exp \{i({\mathit{q}}_{\mu }\cdot \mathit{r}-{\omega }_{\mu }t)\}}.

Equation (10) corresponds to a linear conductivity of {\sigma }_{T=\mathrm{0}}^{(\mathrm{1})}=e^{\mathrm{2}}/\mathrm{\pi }ћ . This is in agreement with the linear conductivity calculated using Kubo formula [ 27, 85]. For µ < 0, the current density reverses the direction since it is now contributed by s = - 1 carriers. For T > 0 K, we obtain

(11) {\mathit{J}}_T^{(\mathrm{1})}=-i{\displaystyle \frac{e^{\mathrm{2}}}{\mathrm{\pi }\mathit{ћ}}}{\displaystyle \frac{k_{B}T}{\mathit{ћ}\omega }}\ln \left[\mathrm{1}+\exp \left(\mathrm{1}+{\displaystyle \frac{\mathit{ћ}\omega }{k_{B}T}}\right)\right]\times {\displaystyle \sum _{\mu }{\mathit{E}}_{\mu }\exp \{i({\mathit{q}}_{\mu }\cdot \mathit{r}-{\omega }_{\mu }t)\}},

which reduces to Eq. (10) in the limit of T→ 0.

Third-order nonlinear response

Due to the inversion symmetry of graphene, it is straightforward to see that the second order velocity {\mathit{v}}_s^{(\mathrm{2})} does not generate any electric current. In fact, one can immediately see this by examining Eq. (6). All of the terms contained in {\mathit{v}}_s^{(\mathrm{2})} are proportional to cos \phi and hence will not survive in the angular-integration of Eq. (8).

The third-order nonlinear current at T = 0 K can be obtained as

(12) {\mathit{J}}_{T=\mathrm{0}}^{(\mathrm{3})}=-is{\displaystyle \frac{e^{\mathrm{4}}{v}_F^{\mathrm{2}}}{\mathrm{8}\mathrm{\pi }{\mathit{ћ}}^{\mathrm{2}}\mu }}{\displaystyle \sum _{\mu v\xi }\left({\displaystyle \frac{{ϵ}_{\mu v\xi }}{\mu -{ϵ}_{\mu v\xi }}}\right)}{\displaystyle \frac{{\mathit{E}}_{\mu }\cdot {\mathit{E}}_v{\mathit{E}}_{\xi }}{{\omega }_{\mu }{\omega }_v{\omega }_{\xi }}}\times \exp \{i[({\mathit{q}}_{\mu }+{\mathit{q}}_v+{\mathit{q}}_{\xi })\cdot \mathit{r}-({\omega }_{\mu }+{\omega }_v+{\omega }_{\xi })t]\},

ω/2π=where s = + 1( - 1) for µ > 0 (µ < 0) and \mu >{ϵ}_{\mu v\xi } where {ϵ}_{\mu v\xi }={ϵ}_{\mu }+{ϵ}_v+{ϵ}_{\xi } is the sum of three incoming photon energies. The magnitude of the zero temperature nonlinear current density is the same for electron filling (µ > 0) and hole filling (µ < 0) due to the up-down Dirac cones symmetry. At finite temperature, the nonlinear current density is obtained from

(13) {\mathit{J}}_T^{(\mathrm{3})}=-is{\displaystyle \frac{e^{\mathrm{4}}{v}_F^{\mathrm{2}}}{\mathrm{8}\mathrm{\pi }{\mathit{ћ}}^{\mathrm{2}}}}{\displaystyle \sum _{\mu v\xi }{\displaystyle \frac{{\mathit{E}}_{\mu }\cdot {\mathit{E}}_v{\mathit{E}}_{\xi }}{{\omega }_{\mu }{\omega }_v{\omega }_{\xi }}}}{\displaystyle \int {\displaystyle \frac{\mathrm{d}{ϵ}_p}{{ϵ}_p^{}}}{\displaystyle \frac{\mathrm{1}}{\mathrm{1}+\exp \left({\displaystyle \frac{{ϵ}_p-\mu }{k_{B}T}}\right)}}}\times \exp \{i[({\mathit{q}}_{\mu }+{\mathit{q}}_v+{\mathit{q}}_{\xi })\cdot \mathit{r}-({\omega }_{\mu }+{\omega }_v+{\omega }_{\xi })t]\},

where for simplicity we have suppressed the integration limit. We see that µ plays an important role in the finite temperature current density of graphene. As shown in Eq. (12), smaller µ generates stronger nonlinear current. However, the assumption of p»u in the derivation of the nonlinear velocities is no longer valid if µ is too small since this will involve charge carriers with momentum comparable to u. For THz waves at room temperature, the range of |µ| \underset{˜}{>} 0.05 eV will be adequate for p»u to hold and we choose µ = 60 meV as the smallest Fermi-level throughout this work. Experimentally, the Fermi level is continuously tunable up to ± (1~2) eV by an external gate voltage [ 86] and hence our choice of µ is practically achievable.

The numerical result of Eq. (13) is shown in Figs. 3 and 4. We observe three important and unusual features in the nonlinear optical response: the third-order nonlinear response is: (i) thermally enhanced up to room temperature; (ii) approximately inversely proportional to µ; and (iii) asymmetric between µ > 0 and µ < 0. Feature (i) is due to the thermal extension of the charge carrier lower-limit µ-ϵ_µvλ-k_BT at higher temperature. Thermally created vacancy at higher energy level allows more low-lying charge carriers to be excited and this amplifies the nonlinear current. However, it should be emphasized that the nonlinear current does not grow indefinitely with increasing temperature. At much higher temperature, the charge carriers in the opposite Dirac cone contribute to and opposite nonlinear current generation which eventually reduces in the net nonlinear current. This reduction is not observed in our case due to the largeness of µ we have chosen, i.e., the nonlinear current is always contributed by charge carriers in only one Dirac cone. For feature (ii), a small µ results in nonlinear current contributed by low-momentum charge carriers and this leads to the stronger current density. The combine effect of (i) and (ii) causes the superlinear growth of nonlinear current at µ = 60 meV and T > 150 K. Feature (iii) is explained by the finite temperature Dirac fermions population distribution in graphene. Consider the µ is in an arbitrary magnitude of µ = µ₀. Switching the Fermi-level from µ = + µ₀ to µ = -µ₀ is essentially equivalent to the mirror reflection of the upper Dirac cone across the zero energy point into the lower Dirac cone. However, the Fermi-Dirac is not reflected, but is shifted downwards by an amount of 2µ₀ and this breaks the overall up-down symmetry of the nonlinear currents at µ = ±µ₀. When µ = +µ₀, a larger amount of low-lying s = +1 electrons become excitable at finite temperature and this significantly enhances the nonlinear current, while in the case of µ = -µ₀, larger amount of deep s = - 1 electrons become excitable and the nonlienar current enhancement is relatively weaker.

Fig.3 Temperature dependence of third order nonlinear current density for µ < 0 at f = ω/2π= 1 THz (Ref. [63])

Full size|PPT slide

Fig.4 Temperature dependence of third order nonlinear current density for µ > 0 at f = ω/2π= 1 THz (Ref. [63])

Full size|PPT slide

The strong nonlinear response of single layer graphene is not surprising if we consider the quasiparticles dynamics in graphene. The massless Dirac fermions around K point is well-described by a ‘pseudospin’ Hamiltonian Eq. (1) and this ‘pseudospin’ nature mimics the ‘real-spin’ Rsahba spin-orbit interaction (RSOI) term in 2-dimensionally electron gas confined in a quantum well structure (2DEG) which has previously been shown to exhibit exceptionally strong nonlinear response [ 87]. In such system, the enhanced optical nonlinerity is caused by the highly non-parabolic band structure induced by RSOI [ 88]; while in graphene, the linear (and hence highly non-parabolic) Dirac conic band structure results in the same enhanced optical nonlinearity. The linear optical response is however much smaller in graphene (linear conductivity in the order of quantum conductance e²/h) and this gives rise to the relatively stronger optical nonlinearity in comparison to 2DEG with RSOI.

Two conclusions can be readily drawn from the above discussions. To achieve strong nonlinear optical effect in graphene: (i) small µ is preferred since a low-lying electron is strongly nonlinear; and (ii) electron filling µ > 0 is preferred due to the broken Dirac fermion population symmetry at finite temperature.

We remark that the total optical conductivity should include both intraband and interband contributions. It can however be seen that σ_inter is forbidden in few-THz regime due to the largeness of µ. By the virtue of momentum conservation, the requirement for vertical interband transition can be written as 3ϵ_photon > 2µ (where for simplicity, the three incoming photons are assumed to have the same energy ϵ_photon). For µ > 0.06 eV, each photon has to exceed 0.04 eV, or frequency higher than 10 THz, for vertical interband transition to become possible and this is well beyond the few-THz regime considered here. Therefore, it is reasonable to drop the σ_inter contribution and to consider σ_intra as the sole contributor to σ_total.

Critical electric field and photon-mixing effect

We now discuss the electric field strength required to create non-negligible photon-mixing effect in graphene. We define a critical field strength such that |J⁽³⁾|/|J⁽¹⁾| = 1. The physical importance of critical field is that it quantifies the optical nonlinearity of a system by comparing both of the linear and nonlinear response. A small critical electric field represents strong optical nonlinearity because the nonlinear response easily dominates over the linear response by only a small electric field.

By combining Eqs. (10) and (12), we obtain the T = 0 K critical field as

(14) E_c(\omega ,T=\mathrm{0}\mathrm{K})={\displaystyle \frac{\mathrm{2}\omega }{v_{\mathrm{F}}}}{\left[{\displaystyle \frac{\mathrm{2}\mu }{e^{\mathrm{2}}}}\left({\displaystyle \frac{\mu }{\mathrm{3}}}-\mathit{ћ}\omega \right)\right]}^{\mathrm{1}/\mathrm{2}},

where the two incident fields are assumed to have the same intensity and polarization. For ω = 1 THz and µ = 0.1 eV, the zero temperature critical field is approximately 10⁴ V/cm. This electric field strength is about one order of magnitude larger than the critical electric field of the 3-photon nonlinear interband conductivity in intrinsic graphene [ 77].

At T > 0, the critical field is E_c(ω, T) = βE_c(ω, T = 0 K) where the temperature dependence is embedded in the dimensionless parameter β:

(15) \beta ={\left\{{\displaystyle \frac{k_{\mathrm{B}}T}{\mathit{ћ}\omega }}{\displaystyle \frac{\ln \left[\mathrm{1}+\exp \left({\displaystyle \frac{\mathit{ћ}\omega }{k_{\mathrm{B}}T}}+\mathrm{1}\right)\right]}{\left|{\mathit{J}}_{T>\mathrm{0}}^{(\mathrm{3})}\right|/\left|{\mathit{J}}_{T=\mathrm{0}}^{(\mathrm{3})}\right|}}\right\}}^{\mathrm{1}/\mathrm{2}}.

β describes the temperature dependence of the optical nonlinearity in graphene. The temperature dependence of β is plotted in Fig. 5. β exhibits contrasting behavior at low and high temperature regime. At low temperature regime, β increases with increasing temperature due to the stronger linear current. At higher temperature, the rate of increase of J⁽³⁾ eventually exceeds J⁽¹⁾ and this leads to the peaking of β, and further increment of temperature results in the lowering of β. For µ = 60 meV, the β peaking is clearly observable at T≈ 150 K. The room temperature E_c is approximately 10% lower than E_c at T≈ 150 K. For µ = 0.1 eV and at room temperature, E_c is increased by about 60%, i.e., E_c≈ 2 × 10⁴ V/cm and this is consistent with the experimental electric field strength where gigahertz waves mixing occurred [ 54].

The nonlinear optical absorption in graphene creates an oscillating current density J⁽³⁾. This oscillation in turns induces an electromagnetic wave giving rise to the well-known four-wave mixing phenomenon. The strong nonlinear current density in graphene immediately suggests the occurrence of strong four-wave mixing effect. The strength of the electric field E⁽³⁾ induced by the nonlinear mixing of ω₃ = 2ω₁±ω₂ can be estimated by solving Maxwell’s inhomogeneous electromagnetic wave equation □E⁽³⁾ = (4π/c²)∂J⁽³⁾/∂t where □ is the d’Alembert operator. At distance far away from the graphene single layer, the solution is approximately given by ∂²E⁽³⁾/∂²t∝∂J⁽³⁾/∂t and the corresponding third-order polarizability is given as

(16) {\chi }^{(\mathrm{3})}={\displaystyle \frac{e^{\mathrm{4}}{v}_{\mathrm{F}}^{\mathrm{2}}}{\mathrm{8}\mathrm{\pi }{ћ}^{\mathrm{2}}\mu }}\left({\displaystyle \frac{ћ{\omega }_{\mathrm{3}}}{\mu -ћ{\omega }_{\mathrm{3}}}}\right){\displaystyle \frac{\mathrm{1}}{{\omega }_{\mathrm{2}}{\omega }_{\mathrm{3}}}}{\left({\displaystyle \frac{\mathrm{1}}{{\omega }_{\mathrm{1}}{ϵ}_{\mathrm{0}}}}\right)}^{\mathrm{2}}.

Fig.5 Temperature dependence of β at f = ω/2π= 1 THz. β exhibits contrasting behavior at low and high temperature regimes [68]

Full size|PPT slide

Optical response of graphene in a strong electric field

In previous sections, the optical response is derived by assuming that the Dirac fermion population is well described by f(ϵ₀) where {ϵ}_{\mathrm{0}}={\mathit{v}}_s^{(\mathrm{0})}\cdot \mathit{p} is the unperturbed linear energy spectrum. Under strong-field condition, the simple assumption of f(ϵ₀) is however no longer valid since the externally acquired dynamics \mathrm{\Delta }ϵ=({\mathit{v}}_s^{(\mathrm{1})}+{\mathit{v}}_s^{(\mathrm{2})}+{\mathit{v}}_s^{(\mathrm{3})})\cdot \mathit{p} is no longer negligible. This additional dynamics causes the Dirac fermions to redistribute themselves via a completely different distribution function of f({ϵ}_{\mathrm{0}})\to f({ϵ}_{\mathrm{0}}+\mathrm{\Delta }ϵ) . In this section, we study the optical response of strong-field driven Dirac fermions (SDF) in graphene with the strong-field induced carrier population redistribution taken into account.

The dynamics of the SDF in graphene can be expressed as

(17) {\mathit{v}}_s^{(\mathrm{0})}\cdot \mathit{p}\to ({\mathit{v}}_s^{()}+{\mathit{v}}_s^{()}+{\mathit{v}}_s^{()}+{\mathit{v}}_s^{()})\cdot \mathit{p}={ϵ}_{\mathrm{0}}+\mathrm{\Delta }ϵ,

where \left|{\mathit{v}}_s^{()}\right|=v_{\mathrm{F}}\approx {\mathrm{10}}^{\mathrm{6}}\mathrm{m}/\mathrm{s} is the Fermi velocity, ϵ₀ is the linear energy dispersion and \mathrm{\Delta }ϵ=({\mathit{v}}_s^{(\mathrm{1})}+{\mathit{v}}_s^{(\mathrm{2})}+{\mathit{v}}_s^{(\mathrm{3})})\cdot \mathit{p} represents the perturbed dynamics due to the strong field. The Fermi-Dirac distribution function can be expanded for small Δϵ. Up to third-order, the expansion yields

(18) f(ϵ)=f_{\mathrm{0}}+\mathrm{\Delta }{f}^{(\mathrm{1})}+\mathrm{\Delta }{f}^{(\mathrm{2})}+\mathrm{\Delta }{f}^{(\mathrm{3})},

where

(19) \mathrm{\Delta }{f}^{(\mathrm{1})}=\mathrm{0}, \mathrm{\Delta }{f}^{(\mathrm{2})}=({\mathit{v}}_s^{(\mathrm{2})}\cdot \mathit{p}){f^{\prime }}_{\mathrm{0}}, \mathrm{\Delta }{f}^{(\mathrm{3})}=({\mathit{v}}_s^{\mathrm{3}}\cdot \mathit{p}){f^{\prime }}_{\mathrm{0}},

where {f^{\prime }}_{\mathrm{0}} is the first derivative of the Fermi-Dirac distribution function with respect to ϵ₀. Equation (18) is the strong-field Fermi-Dirac distribution which includes nonlinear terms up to third-order in the external electric field. We can then calculate the total current of the SDF by the following equation:

(20) \mathit{J}=e{\displaystyle \sum _s{\displaystyle {\int }_{\mathrm{0}}^{\mathrm{2}\mathrm{\pi }}{\displaystyle {\int }_{\mu -{ϵ}_{ph}-k_{B}T}^{\mathrm{\Lambda }}{\mathrm{d}}^{\mathrm{2}}p({\mathit{v}}_s^{(\mathrm{0})}+{\mathit{v}}_s^{(\mathrm{1})}+{\mathit{v}}_s^{(\mathrm{2})}+{\mathit{v}}_s^{(\mathrm{3})})(f_{\mathrm{0}}+\mathrm{\Delta }{f}^{(\mathrm{1})}+\mathrm{\Delta }{f}^{(\mathrm{2})}+\mathrm{\Delta }{f}^{(\mathrm{3})})}}}.

Splitting the total current into linear and nonlinear components, we obtain [ 68]

(21) \begin{array}{l}{\mathit{J}}^{(\mathrm{1}S)}=e{\displaystyle \sum _s{\displaystyle {\int }_{\mathrm{0}}^{\mathrm{2}\mathrm{\pi }}{\displaystyle {\int }_{\mu -{ϵ}_{ph}-k_{B}T}^{\mathrm{\Lambda }}{\mathrm{d}}^{\mathrm{2}}p({\mathit{v}}_s^{(\mathrm{1})}{f}_{\mathrm{0}}+{\mathit{v}}_s^{(\mathrm{0})}\mathrm{\Delta }{f}^{(\mathrm{1})})}}},\hfill \\ {\mathit{J}}^{(\mathrm{2}S)}=e{\displaystyle \sum _s{\displaystyle {\int }_{\mathrm{0}}^{\mathrm{2}\mathrm{\pi }}{\displaystyle {\int }_{\mu -{ϵ}_{ph}-k_{B}T}^{\mathrm{\Lambda }}{\mathrm{d}}^{\mathrm{2}}p({\mathit{v}}_s^{(\mathrm{2})}{f}_{\mathrm{0}}+{\mathit{v}}_s^{(\mathrm{0})}\mathrm{\Delta }{f}^{(\mathrm{2})}+{\mathit{v}}_s^{(\mathrm{1})}\mathrm{\Delta }{f}^{(\mathrm{1})})}}},\hfill \\ {\mathit{J}}^{(3S)}=e{\displaystyle \sum _s{\displaystyle {\int }_{\mathrm{0}}^{\mathrm{2}\mathrm{\pi }}{\displaystyle {\int }_{\mu -{ϵ}_{ph}-k_{B}T}^{\mathrm{\Lambda }}{\mathrm{d}}^{\mathrm{2}}p({\mathit{v}}_s^{(\mathrm{3})}{f}_{\mathrm{0}}+{\mathit{v}}_s^{(\mathrm{0})}\mathrm{\Delta }{f}^{(\mathrm{3})}+{\mathit{v}}_s^{(\mathrm{1})}\mathrm{\Delta }{f}^{(\mathrm{2})}+{\mathit{v}}_s^{(\mathrm{2})}\mathrm{\Delta }{f}^{(\mathrm{1})})}}},\hfill \end{array}

where the superscript (S) emphasizes the optical response of SDF. The term {\mathit{J}}^{(iw)}=e{\displaystyle {\sum }_s\int {\mathrm{d}}^{\mathrm{2}}p{\mathit{v}}_s^{(i)}{f}_{\mathrm{0}}} is the weak-field Dirac fermions optical response. After some algebra, we obtain

(22) {\mathit{J}}^{(\mathrm{1}S)}={\mathit{J}}^{(\mathrm{1}w)},

(23) {\mathit{J}}^{(\mathrm{2}S)}=\mathrm{0},

(24) {\mathit{J}}^{(\mathrm{3}S)}={\mathit{J}}^{(\mathrm{3}w)}+{\mathit{J}}^{({\mathrm{3}}^{\prime })}.

The superscript (S) and (w) emphasize the optical response of SDF and weak-field Dirac fermions respectively. The consequences of Eqs. (22) and (23) are quite surprising: the linear and second-order nonlinear optical responses of graphene remain unchanged although the whole SDF population has redistributed themselves. This behavior can be understood by considering the nature of the strong-field induced population redistribution phenomena. Such process is a description of how strongly the Dirac fermions respond to an external perturbation and the degree of redistribution depends on the coupling between the externally acquired dynamics and the unperturbed dynamics of Dirac fermions, i.e., v_external ∙p. For first-order response, it can be seen that the externally acquired first-order dynamics is completely decoupled from the unperturbed dynamics, i.e., {\mathit{v}}_s^{(\mathrm{1})}\cdot \mathit{p}=0 . As a result, this orthogonality ensures that the linear response of graphene is always protected from the strong field effect. For second-order nonlinear response, the second-order coupling {\mathit{v}}_s^{(\mathrm{2})}\cdot \mathit{p} is finite and one would intuitively expect a finite second-order current to occur. This is however not the case as the additional second-order term vanishes after performing angular integration. In this case, although Dirac fermions are second-order-ly perturbed and redistributed, the crystal itself remains unaffected and retains its inversion symmetry. Therefore, second-order nonlinear response is still zero in the strong-field regime.

At T = 0 K, the third-order nonlinear optical response of SDF is

(25) {\mathit{J}}_{T=\mathrm{0}}^{({\mathrm{3}}^{\prime })}=-is{\displaystyle \sum _{\mu v\xi }{\displaystyle \frac{{\mathit{E}}_{\mu }\cdot {\mathit{E}}_v{\mathit{E}}_{\xi }}{{\omega }_{\mu }{\omega }_v{\omega }_{\xi }}}{\displaystyle \frac{e^{\mathrm{4}}{v}_F^{\mathrm{2}}}{\mathrm{8}\mathrm{\pi }{ћ}^{\mathrm{2}}\mu }}}\times \exp \{i[({\mathit{q}}_{\mu }+{\mathit{q}}_v+{\mathit{q}}_{\xi })\cdot \mathit{r}-({\omega }_{\mu }+{\omega }_v+{\omega }_{\xi })t]\}.

Finally, in the general case of T > 0 K, we have

(26) {\mathit{J}}^{({\mathrm{3}}^{\prime })}={\mathit{J}}_{T=\mathrm{0}}^{({\mathrm{3}}^{\prime })}{\displaystyle \frac{\mu }{k_{B}T}}{\displaystyle \int {\displaystyle \frac{\mathrm{d}p}{p}}{\displaystyle \frac{\exp \left({\displaystyle \frac{{ϵ}_{\mathrm{0}}-\mu }{k_{B}T}}\right)}{{\left(\exp \left({\displaystyle \frac{{ϵ}_{\mathrm{0}}-\mu }{k_{B}T}}\right)+\mathrm{1}\right)}^{\mathrm{2}}}}}.

Critical electric field in the strong-field regime

At T = 0 K, the critical electric field of strongly-driven massless Dirac fermion in graphene can be obtained by directly taking the ratio of the linear and nonlinear current densities derived in previous sections. This gives

(27) E_c^{(S)}(T=\mathrm{0})={\displaystyle \frac{\mathrm{2}\omega }{v_F}}\left[{\displaystyle \frac{\mathrm{2}\mathit{ћ}\omega }{e^{\mathrm{2}}}}(\mu -\mathrm{3}\mathit{ћ}\omega )\right].

Fig.6 Critical field of E_{\mathrm{c}}^{(\mathrm{S})} at f = ω/2π= 1 THz and µ = 0.1 eV. Weak-field critical field E_c is also shown [ 68]

Full size|PPT slide

For 1 THz and µ = 0.1 eV, E_{\mathrm{c}}^{(\mathrm{S})}(T=\mathrm{0})=\mathrm{3300}\mathrm{V}/\mathrm{c}\mathrm{m} and is 3 times smaller than that of the weak-field response (≈10⁴ V/cm). At finite temperature, we obtain E_{\mathrm{c}}^{(\mathrm{S})}(T)={\beta }^{(\mathrm{S})}(T)E_{\mathrm{c}}^{(\mathrm{S})}(T=\mathrm{0}) , where the dimensionless strong-field {\beta }^{(\mathrm{S})}(T) is given as

(28) {\beta }^{(\mathrm{S})}(T)={\left\{{\displaystyle \frac{k_{\mathrm{B}}T}{ћ\omega }}{\displaystyle \frac{\ln \left[\mathrm{1}+\exp \left({\displaystyle \frac{ћ\omega }{k_{\mathrm{B}}T}}+\mathrm{1}\right)\right]}{\left|{\mathit{J}}_{T>\mathrm{0}}^{(\mathrm{3}\mathrm{S})}\right|/\left|{\mathit{J}}_{T=\mathrm{0}}^{(\mathrm{3}\mathrm{S})}\right|}}\right\}}^{\mathrm{1}/\mathrm{2}}.

It can be seen in Fig. 6 that E_{\mathrm{c}}^{(\mathrm{S})} is significantly lower than weak-field E_c over a wide temperature regime from T = 0 K to T = 600 K. This indicates the stronger optical nonlinearity of SDF in comparison to the usual Dirac fermions. The stronger optical nonlinearity of SDF is due to the fact that the third-order nonlinear response is amplified by {\mathit{J}}^{({\mathrm{3}}^{\prime })} while the linear response remains unchanged.

We now discuss the optical response due to non-equilibrium hot Dirac fermions in graphene. The hot Dirac fermions in graphene are short-lived especially in the case of high lattice temperature where stronger electron-phonon coupling provides efficient pathway for the relaxation [ 89- 91]. Under weak-field condition, Dirac fermions rapidly thermalize themselves with the lattice, i.e., T = T_lattice. In strong-field regime, the non-equilibrium heating of SDF lifted the SDF temperature from lattice temperature and hence the temperature terms in Eqs. (22) and (24) has to be replaced by T→T_hot where T_hot is the hot SDF temperature and T_hot> T_lattice. For critical field varies between 10³ to 10⁴ V/cm, the hot SDF temperature reaches between T_hot = 350 K to T_hot = 600 K [ 92]. In contrast, equilibrium Dirac fermions are relatively ‘cold’ since the lattice temperature in most of the practical application is only up to T_lattice = 300 K. It can be seen from Fig. 7 that the nonlinear current of hot SDF in 350 K < T_hot < 600 K is generally stronger than that of the cold equilibrium Dirac fermions where T_lattice< 300 K.

Fig.7 Temperature dependence of strong-field third order nonlinear current density at f = ω/2π= 1 THz. Note that T = T_lattice if non-equilibrium heating is ignored and T = T_hot if non-equilibrium heating is considered. Since T_hot > T_lattice, the nonlinear optical response is significantly stronger if carrier heating is considered [68]

Full size|PPT slide

Nonlinear optical response of gapped graphene in a strong electric field

For completeness, we briefly discuss the nonlinear intraband conductivity of gapped graphene. The Hamiltonian of gapped graphene in low energy regime is given as

(29) \widehat{\mathit{H}}=\left[\begin{array}{ll}-\mathrm{\Delta }& v_{\mathrm{F}}{p}^{+}\\ v_{\mathrm{F}}{p}_{-}& \mathrm{\Delta }\end{array}\right],

where Δ is the onsite energy difference in the sublattice A and B respectively. The energy eigenvalue is given by {ϵ}_s=s\sqrt{v_{\mathrm{F}}^{\mathrm{2}}{p}^{\mathrm{2}}+{\mathrm{\Delta }}^{\mathrm{2}}} indicating a bandgap opening of 2Δ at the Dirac point. Similarly, up to third-order in external field, the strong-field perturbed carrier velocities can be written as

(30) \begin{array}{c}{v}_s^{\left(0\right)}=s{\displaystyle \frac{v_{F}p}{\sqrt{p^2+{\mathrm{\Delta }}^2}}},\\ v_s^{\left(1\right)}=s{\displaystyle \frac{v_F}{\sqrt{p^2+{\mathrm{\Delta }}^2}}}\left[u-p\left({\displaystyle \frac{p\cdot u}{p^2+{\mathrm{\Delta }}^2}}\right)\right],\\ v_s^{\left(2\right)}=s{\displaystyle \frac{v_F}{\sqrt{p^2+{\mathrm{\Delta }}^2}}}\left[-{\displaystyle \frac{p}{2}}{\displaystyle \frac{u^2}{p^2+{\mathrm{\Delta }}^2}}-u{\displaystyle \frac{p\cdot u}{p^2+{\mathrm{\Delta }}^2}}+{\displaystyle \frac{3p}{2}}{\left({\displaystyle \frac{p\cdot u}{p^2+{\mathrm{\Delta }}^2}}\right)}^2\right],\\ v_s^{\left(3\right)}=s{\displaystyle \frac{v_F}{\sqrt{p^2+{\mathrm{\Delta }}^2}}}\left[-{\displaystyle \frac{u}{2}}{\displaystyle \frac{u^2}{p^2+{\mathrm{\Delta }}^2}}+{\displaystyle \frac{3u^{2}p}{2}}{\displaystyle \frac{p\cdot u}{p^2+{\mathrm{\Delta }}^2}}+{\displaystyle \frac{3u}{2}}{\left({\displaystyle \frac{p\cdot u}{p^2+{\mathrm{\Delta }}^2}}\right)}^2-{\displaystyle \frac{5p}{2}}{\left({\displaystyle \frac{p\cdot u}{p^2+{\mathrm{\Delta }}^2}}\right)}^3\right].\end{array}

The carrier dynamics becomes ϵ = ϵ_s + Δϵ where {ϵ}_s=s\sqrt{v_{\mathrm{F}}^{\mathrm{2}}{p}^{\mathrm{2}}+{\mathrm{\Delta }}^{\mathrm{2}}} is the unperturbed energy spectrum and \mathrm{\Delta }ϵ=({\mathit{v}}_s^{(\mathrm{1})}+{\mathit{v}}_s^{(\mathrm{2})}+{\mathit{v}}_s^{(\mathrm{3})})\cdot \mathit{p} is the field-induced energy changed. In strong-field case, the carrier population redistribute themselves according to

(31) f({ϵ}_s+\mathrm{\Delta }ϵ)=f_{\mathrm{0}}+\mathrm{\Delta }{f}_{\mathrm{1}}+\mathrm{\Delta }{f}_{\mathrm{2}}+\mathrm{\Delta }{f}_{\mathrm{3}},

where

\begin{array}{c}\mathrm{\Delta }{f}_1=\left(v_s^{\left(1\right)}\cdot p\right){f^{\prime }}_0,\\ \mathrm{\Delta }{f}_2=\left(v_s^{\left(2\right)}\cdot p\right){f^{\prime }}_0+{\displaystyle \frac{{\left(v_s^{\left(1\right)}\cdot p\right)}^2}{2}}{f^{″}}_0,\\ \mathrm{\Delta }{f}_3=\left(v_s^{\left(1\right)}\cdot p\right){f^{\prime }}_0+\left(v_s^{\left(1\right)}\cdot p\right)\left(v_s^{\left(2\right)}\cdot p\right){f^{″}}_0+{\displaystyle \frac{{\left(v_s^{\left(1\right)}\cdot p\right)}^3}{6}}{f^{‴}}_0.\end{array}

We now write the current densities as

(32) \begin{array}{l}{J}^{(\mathrm{1})}=J^{(\mathrm{1w})}+J^{(\mathrm{1s})},\hfill \\ J^{(\mathrm{2})}=J^{(\mathrm{2w})}+J^{(\mathrm{2s})},\hfill \\ J^{(\mathrm{3})}=J^{(\mathrm{3w})}+J^{(\mathrm{3s})},\hfill \end{array}

where the superscript (w) and (s) denote weak-and strong-field term respectively. Using similar strategy as discussed in previous sections, we obtain the T = 0 K total responses

(33) {\mathit{J}}_{T=\mathrm{0}}^{(\mathrm{1})}=-s{\displaystyle \frac{e^{\mathrm{2}}\mathit{E}}{\mathrm{\pi }\mathit{ћ}}}[x+\mathrm{1}+{\mathrm{\Delta }}_{\omega }{\mathrm{\Delta }}_{\mu }(\mathrm{1}-{\mathrm{\Delta }}_{\mu }^{\mathrm{2}})],

(34) J_{T=\mathrm{0}}^{(\mathrm{2})}=\mathrm{0},

(35) J_{T=\mathrm{0}}^{(\mathrm{3})}={\displaystyle \frac{s{v}_{\mathrm{F}}^{\mathrm{2}}{u}^{}}{\mathrm{8}\mathrm{\pi }{\mathit{ћ}}^{\mathrm{2}}\mu }}(y+\mathrm{1}-{\mathrm{\Delta }}_{\mu }^{\mathrm{2}}+\mathrm{59}{\mathrm{\Delta }}_{\mu }^{\mathrm{4}}-\mathrm{177}{\mathrm{\Delta }}_{\mu }^{\mathrm{6}}+\mathrm{208}{\mathrm{\Delta }}_{\mu }^{\mathrm{8}}-\mathrm{90}{\mathrm{\Delta }}_{\mu }^{\mathrm{10}}),

where the weak-field term is represented by x and y, which can be explicitly written as

(36) x={\displaystyle \frac{{\displaystyle \frac{\mu }{\mathit{ћ}\omega }}}{{\left[\mathrm{1}+{\left({\displaystyle \frac{\mathrm{\Delta }}{\mu }}\right)}^{\mathrm{2}}\right]}^{\mathrm{1}/\mathrm{2}}}}-{\displaystyle \frac{{\displaystyle \frac{\mu }{\mathit{ћ}\omega }}-\mathrm{1}}{{\left[\mathrm{1}+{\left({\displaystyle \frac{\mathrm{\Delta }}{\mu -\mathit{ћ}\omega }}\right)}^{\mathrm{2}}\right]}^{\mathrm{1}/\mathrm{2}}}},

and

(37) y={\displaystyle \frac{\mu }{\mu -{ϵ}_{\mathrm{ph}}}}{\displaystyle \frac{\mathrm{1}+{\left({\displaystyle \frac{\mathrm{2}\mathrm{\Delta }}{\mu -{ϵ}_{\mathrm{p}\mathrm{h}}}}\right)}^{\mathrm{2}}}{{\left[\mathrm{1}+{\left({\displaystyle \frac{\mathrm{\Delta }}{\mu -{ϵ}_{\mathrm{p}\mathrm{h}}}}\right)}^{\mathrm{2}}\right]}^{\mathrm{5}/\mathrm{2}}}}-{\displaystyle \frac{\mathrm{1}+{\left({\displaystyle \frac{\mathrm{2}\mathrm{\Delta }}{\mu }}\right)}^{\mathrm{2}}}{{\left[\mathrm{1}+{\left({\displaystyle \frac{\mathrm{\Delta }}{\mu }}\right)}^{\mathrm{2}}\right]}^{\mathrm{5}/\mathrm{2}}}},

where ϵ_ph is the energy sum of the three photons. By solving the current integral in Eq. (8), one can show that the second order nonlinear current is zero since {\mathit{v}}_s^{(\mathrm{2})} is an odd function of \varphi (see Eq. (30)). One major difference between gapless and gapped graphene is that both of the linear and third-order nonlinear current density are altered by a strong electric field in gapped graphene whereas in gapless graphene, only the third-order nonlinear current density is altered. The third-order nonlinear conductivity at finite temperature is evaluated numerically. The temperature and bandgap dependence of the intraband third-order nonlinear current density at f = 1.5 THz and µ = 0.12 eV is shown in Fig. 8. Unlike the gapless graphene (Δ→ 0) which exhibits enhanced nonlinear third-order nonlinear optical response at elevated temperature, the third-order nonlinear response of gapped graphene is sensitively influenced by Δ and temperature. Due to the interplay between the bandgap opening and carrier thermal excitation, two distinct nonlinear optical response ‘hotspots’, in which an amplification factor of ≈ 3.5, are created at two regimes: (i) low temperature with large Δ; and (ii) high temperature with small Δ. These ‘hotspots’ are connected by a region in which the nonlinear optical response is two times higher than the linear optical response.

The low temperature hotspot (i) indicates that the bandgap opening in graphene effectively enhances the nonlinear response. The nonlinear response enhancement due to bandgap opening in graphene also occurs in the interband nonlinear optical response [ 70]. This suggests that both of the interband and the intraband optical absorption are universally enhanced by the bandgap opening in graphene. At very large bandgap, the nonlinear optical response however decreases. This is because the nonlinear velocity component is approximately \mathrm{\propto }{p}^{-\mathrm{3}} . A large bandgap destroys low momentum states and hence severely degrades the nonlinear optical response. The high temperature hotspot (ii) is a thermal effect. The thermal excitation vacates states beneath the Fermi level, allowing deeper charge carriers to become optically excitable. At very high temperature, the thermally excitable charge carrier population extends to the edge of the bandgap. Any further increment of the temperature does not increase the optically excitable charge carriers. On the other hand, the overall charge carrier momentum is elevated thermally and the p^-3 reduction of the nonlinear velocity takes place. The combination of these two aspects results the high temperature degradation of the nonlinear optical response.

Fig.8 ∆ and temperature dependence of the third-order nonlinear current density [68] at f = 1.5 THz and µ = 0.12eV

Full size|PPT slide

Bilayer graphene

In this section, we study the nonlinear interband optical response in BLG in the frequency regime of THz to FIR. The interband optical response is obtained by using a quantum mechanical treatment that couples the BLG quasiparticle to a time-dependent electric field [77]. We expressed the light-dressed electron wave function as an infinite sum in terms of the number of photons coupled to the massless Dirac fermion. This allows us to explicitly construct the nonlinear optical current density up to any arbitrary order in the external electric field.

We show that the optical response of BLG is significantly enhanced due to the third-order nonlinear process. The nonlinear effects are particularly strong in the low frequency regime, which covers the technologically important frequency band of THz to FIR. More importantly, the field intensity required for the onset of nonlinear response is rather low, indicating that BGL is an excellent material for nonlinear optics and photonics application. The third-order nonlinear optical response is composed of two terms: (i) single-frequency term which corresponds to the simultaneous absorption of two photons and the emission of one photon; and (ii) triple-frequency term which corresponds to the simultaneous absorption of three photons. Both of (i) and (ii) become comparable to the linear optical response at very moderate electric field of 10³ V/cm which is well within the experimental achievable range in laboratories. Furthermore, we investigate the temperature dependence of the nonlinear optical response. At room temperature, we found that the electric field required to produce nonlinear optical response comparable to the linear one is reduced to 10² V/cm. This thermally enhanced optical nonlinearity is not found in single layer graphene. This suggests that BLG is a preferred structure for developing graphene-based nonlinear photonics and optoelectronics device.

Recursion equations for n-photon-electron coupling

We consider the case where a time-dependent electric field E(t) = Ee^iωt is applied along the x-axis. We start with the low energy effective Hamiltonian of BLG in K valley [ 93- 95]:

(38) \mathit{H}=\alpha \left(\begin{array}{cc}0& {\left(p_{-}+eA\right)}^2\\ {\left(p_{+}+eA\right)}^2& 0\end{array}\right)-\beta \left(\begin{array}{cc}0& p_{+}+eA\\ p_{-}+eA& 0\end{array}\right),

where p_± = p_x ± ip_y, A={\displaystyle \frac{E}{i\omega }}{e}^{i\omega t} , α = 1/(2m^*), m^* = 0.033 m_e, and β = v_F≈ 10⁵ m/s [ 95]. Note that in Eq. (38), we have performed a Peierls substitution of p + eA. The Hamiltonian can be rearranged to the following form

(39) \mathit{H}=\left[\begin{array}{cc}0& Y_{-}+{\displaystyle \frac{eE}{i\omega }}{X}_{-}{e}^{i\omega t}-\alpha {\displaystyle \frac{e^2{E}^2}{{\omega }^2}}{e}^{i2\omega t}\\ Y_{+}+{\displaystyle \frac{eE}{i\omega }}{X}_{+}{e}^{i\omega t}-\alpha {\displaystyle \frac{e^2{E}^2}{{\omega }^2}}{e}^{i2\omega t}& 0\end{array}\right],

where for simplicity, we denote

(40) X_{\pm }=\mathrm{2}\alpha p_{\pm }-\beta ,

(41) Y_{\pm }=\alpha p_{\pm }^{\mathrm{2}}-\beta p_{\mp }.

The electron-photon coupled wave function can be written as an infinite sum in terms of the number of photons:

(42) \psi (\mathit{p},n)={\displaystyle \sum _{n=\mathrm{0}}^{\infty }\left(\begin{array}{l}{a}_n(\mathit{p})\\ b_n(\mathit{p})\end{array}\right)}{e}^{i(n\omega -ϵ/ћ)t},

where (a_n(p), b_n(p))^T are the spinor components representing n-photon coupling of the electron. The time derivatives of the wave function is

(43) {\displaystyle \frac{\partial \psi }{\partial t}}(\mathit{p},n)={\displaystyle \sum _{n=\mathrm{0}}^{\infty }i(n\omega -ϵ/\mathit{ћ})}\left(\begin{array}{l}{a}_n(\mathit{p})\\ b_n(\mathit{p})\end{array}\right)e^{i(n\omega -ϵ/\mathit{ћ})t}.

The spinor components can be obtained by solving the Schrodinger equation iћ{\displaystyle \frac{\partial \psi }{\partial t}}=H\psi . Combining Eqs. (38) and (43) with the Schrodinger equation, we have

(44) {\displaystyle \sum _{n=0}^{\infty }\left[\begin{array}{cc}0& Y_{-}+{\displaystyle \frac{eE}{i\omega }}{X}_{-}{e}^{i\omega t}-\alpha {\displaystyle \frac{e^2{E}^2}{{\omega }^2}}{e}^{i2\omega t}\\ Y_{+}+{\displaystyle \frac{eE}{i\omega }}{X}_{+}{e}^{i\omega t}-\alpha {\displaystyle \frac{e^2{E}^2}{{\omega }^2}}{e}^{i2\omega t}& 0\end{array}\right]}\left(\begin{array}{c}{a}_n\left(\mathit{p}\right)\\ b_n\left(\mathit{p}\right)\end{array}\right)e^{i\left(n\omega -ϵ/\mathit{ћ}\right)t}={\displaystyle \sum _{n=0}^{\infty }\left(ϵ-n\mathit{ћ}\omega \right)\left(\begin{array}{c}{a}_n\left(\mathit{p}\right)\\ b_n\left(\mathit{p}\right)\end{array}\right)e^{i\left(n\omega -ϵ/\mathit{ћ}\right)t}}.

The e^iωt and e^i2ωt terms cane be absorbed into the spinor components to obtain a_n-1, a_n-2 and b_n-1, b_n-2, respectively. Due to the off-diagonal nature of the Hamiltonian in Eq. (38), the upper and the lower spinor components a_n and b_n are coupled and two recursion equations can be obtained

(45) \begin{array}{l}(ϵ-n\mathit{ћ}\omega )a_n=Y_{-}{b}_n+{\displaystyle \frac{eE}{i\omega }}{X}_{-}{b}_{n-\mathrm{1}}-\alpha {\displaystyle \frac{e^{\mathrm{2}}{E}^{\mathrm{2}}}{{\omega }^{\mathrm{2}}}}{b}_{n-\mathrm{2}},\hfill \\ (ϵ-n\mathit{ћ}\omega )b_n=Y_{+}{a}_n+{\displaystyle \frac{eE}{i\omega }}{X}_{+}{a}_{n-\mathrm{1}}-\alpha {\displaystyle \frac{e^{\mathrm{2}}{E}^{\mathrm{2}}}{{\omega }^{\mathrm{2}}}}{a}_{n-\mathrm{2}}.\hfill \end{array}

The above equation contains information of all multiple photon processes in intrinsic graphene. The recursion relation couples the n photon processes to the n-1 photon processes. In general, we can write

(46) n\mathit{ћ}\omega (n\mathit{ћ}\omega -\mathrm{2}ϵ)a_n={\displaystyle \frac{eE}{i\omega }}[X_{-}(ϵ-n\mathit{ћ}\omega )b_{n-\mathrm{1}}+X_{+}{Y}_{-}{a}_{n-\mathrm{1}}]-\alpha {\displaystyle \frac{e^{\mathrm{2}}{E}^{\mathrm{2}}}{{\omega }^{\mathrm{2}}}}[Y_{-}{a}_{n-\mathrm{2}}+(ϵ-n\mathit{ћ}\omega )b_{n-\mathrm{2}}],

(47) n\mathit{ћ}\omega (n\mathit{ћ}\omega -\mathrm{2}ϵ)b_n={\displaystyle \frac{eE}{i\omega }}[X_{-}{Y}_{+}{b}_{n-\mathrm{1}}+X_{+}{a}_{n-\mathrm{1}}(ϵ-n\mathit{ћ}\omega )]-\alpha {\displaystyle \frac{e^{\mathrm{2}}{E}^{\mathrm{2}}}{{\omega }^{\mathrm{2}}}}[(ϵ-n\mathit{ћ}\omega )a_{n-\mathrm{2}}+Y_{+}{b}_{n-\mathrm{2}}].

For n = 0, there is no photon. The spinor components can be solved to obtain

(48) a_{\mathrm{0}}={\displaystyle \frac{Y_{-}}{ϵ\sqrt{\mathrm{2}}}},

(49) b_{\mathrm{0}}={\displaystyle \frac{\mathrm{1}}{\sqrt{\mathrm{2}}}},

where ϵ is the energy dispersion as given by

(50) ϵ=\sqrt{Y_{+}{Y}_{-}}=\pm \sqrt{{\alpha }^{\mathrm{2}}{p}^{\mathrm{4}}-\mathrm{2}\alpha \beta p^{\mathrm{3}}\cos \mathrm{3}\theta +{\beta }^{\mathrm{2}}{p}^{\mathrm{2}}}.

The n = 0 no-photon spinor components are in agreement with the single particle eigenfunction of BLG [ 95]. For n = 1, we obtain

(51) a_{\mathrm{1}}={\displaystyle \frac{eE}{i\sqrt{\mathrm{2}}\mathit{ћ}{\omega }^{\mathrm{2}}ϵ(\mathit{ћ}\omega -\mathrm{2}ϵ)}}[ϵ(ϵ-\mathit{ћ}\omega )X_{-}+X_{+}{Y}_{-}^{\mathrm{2}}],

(52) b_{\mathrm{1}}={\displaystyle \frac{eE}{i\sqrt{\mathrm{2}}\mathit{ћ}{\omega }^{\mathrm{2}}ϵ(\mathit{ћ}\omega -\mathrm{2}ϵ)}}[ϵX_{-}{Y}_{+}+(ϵ-\mathit{ћ}\omega )X_{+}{Y}_{-}].

This gives the spinor components of one-photon coupling. The two-photon terms can be recursively built by combining a₀, b₀ and a₁, b₁ into Eq. (45), and so on.

Optical current operator and density

We now construct the current density created by the external time-dependent electric field. The velocity operator in x-direction, v_x = ∂H/∂p_x, is given by

(53) {\widehat{\mathit{v}}}_x=2\alpha \left(\begin{array}{cc}0& \left(p_{-}+eA\right)\\ \left(p_{+}+eA\right)& 0\end{array}\right)-\beta \left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right)={\widehat{\mathit{v}}}_{\mathrm{A}}+{\widehat{\mathit{v}}}_{\mathrm{B}},

where

(54) {\widehat{\mathit{v}}}_{\mathrm{A}}=2\alpha \left(\begin{array}{cc}0& \left(p_{-}+eA\right)\\ (p_{+}+eA)& 0\end{array}\right),

(55) {\widehat{\mathit{v}}}_B=-\beta \left(\begin{array}{ll}\mathrm{0}& \mathrm{1}\\ \mathrm{1}& \mathrm{0}\end{array}\right),

(55)where {\widehat{\mathit{v}}}_{\mathrm{A}} originates from the p quadratic term of the Hamiltonian Eq. (38) and {\widehat{\mathit{v}}}_{\mathrm{B}} originates from the p linear term. In single layer graphene, the velocity operator is only contains {\widehat{\mathit{v}}}_{\mathrm{B}} . In BLG, the interlayer coupling creates an additional {\widehat{\mathit{v}}}_{\mathrm{B}} .

The total x directional optical current operator is given by

(56) \widehat{j}=-e{\psi }^{†}\widehat{\mathit{v}}\psi =-e\left({\psi }^{†}{\widehat{\mathit{v}}}_{\mathrm{A}}\psi +{\psi }^{†}{\widehat{\mathit{v}}}_{\mathrm{B}}\psi \right),

where the wave function is given by Eq. (43). The total current density can be obtained by integrating Eq. (56) in p-space, i.e.,

(57) J={\displaystyle \frac{\mathrm{1}}{{(\mathrm{2}\mathrm{\pi }\mathit{ћ})}^{\mathrm{2}}}}\Re {\displaystyle \int \mathrm{d}\mathit{p}\widehat{j}N(ϵ)},

where the temperature dependence of the total current density is encoded in N(ϵ) = f( -ϵ) -f(ϵ) = tanh(ϵ/2k_BT). Since \mathit{\psi } in Eq. (56) is a linear superposition of the spinor components n-th order, we can selectively construct J_n in any arbitrary order n in the external electric field. For example, the n = 1 linear optical current operator is

(58) j_{\mathrm{1}}={\widehat{j}}_{\mathrm{1}}^{(\mathrm{A})}+{\widehat{j}}_{\mathrm{1}}^{(\mathrm{B})}+{{\widehat{j}}^{\prime }}_{\mathrm{1}},

where

(59) {\widehat{j}}_{\mathrm{1}}^{(\mathrm{A})}=-\mathrm{2}e\alpha [(a_{\mathrm{1}}{b}_{\mathrm{0}}^{*}{p}_{+}+b_{\mathrm{1}}{a}_{\mathrm{0}}^{*}{p}_{-})+(a_{\mathrm{1}}^{*}{b}_{\mathrm{0}}{p}_{+}+b_{\mathrm{1}}^{*}{a}_{\mathrm{0}}{p}_{-})],

(60) {\widehat{j}}_{\mathrm{1}}^{(\mathrm{B})}=e\beta (a_{\mathrm{0}}^{*}{b}_{\mathrm{0}}+b_{\mathrm{0}}^{*}{a}_{\mathrm{0}}),

and

(61) {{\widehat{j}}^{\prime }}_{\mathrm{1}}={\displaystyle \frac{E}{i\omega }}(a_{\mathrm{0}}^{*}{b}_{\mathrm{1}}+b_{\mathrm{0}}^{*}{a}_{\mathrm{1}})+c.c.

The single-photon linear optical current density is hence given by

(62) J_{\mathrm{1}}={\displaystyle \frac{\mathrm{1}}{\mathrm{4}{\mathrm{\pi }}^{\mathrm{2}}{\mathit{ћ}}^{\mathrm{2}}}}\Re {\displaystyle \int {\widehat{j}}_{\mathrm{1}}p\mathrm{d}p\mathrm{d}\theta .}

It is obvious that the linear optical current is made up of the multiples of n = 1 and n = 0 spinor components (e.g., a_{\mathrm{0}}^{*}{b}_{\mathrm{1}} since it is a one-photon process. For the third-order nonlinear optical response, the three-photon process can be either composed of the multiples of n = 0 and n = 3 spinor components or the n = 1 and n = 2 spinor components. These respectively represent optical processes of simultaneous three-photon absorption and the simultaneous two-photon absorption followed by one-photon emission. J₃ is therefore composed of a single-frequency term ∝e^iωt and a triple-frequency term ∝e^i3ωt.

Fig.9 Frequency dependence of the linear optical conductivity at different interlayer coupling strength. Dotted curve: 0.1α, dashed curve: 0.5α, solid curve: α, dash-dotted curve: 1.5α, where α = 1/(2m*) and m* = 0.33 m_e [ 95]. The low frequency conductivity always approach 6σ₀ regardless the strength of the interlayer coupling [ 69]

Full size|PPT slide

Linear optical response

The frequency dependence of the linear optical conductivity is plotted in Fig. 9. There is a conductivity peak corresponding to the singularity in the density of states when the low energy Dirac pockets joint together [94]. Furthermore, it can be seen that in the limit of ω→ 0, the conductivity approaches 6σ₀ where σ₀ = πe²/2h = e²/4ħ. This is in agreement with the linear response result obtained from the Kubo formula [ 24]. Several values of the interlayer coupling strength α is chosen, and the low frequency conductivity is always 6σ₀. As suggested by Cserti et al., the universal minimum conductivity of 6σ₀ regardless the interlayer coupling strength is of topological origin [ 24].

Nonlinear optical response

The nonlinear optical response is numerically evaluated. In Fig. 10, we plot the nonlinear conductance versus frequency in unit of 6σ₀ for two different temperatures. The electric field is 1000 V/cm. All nonlinear terms decrease rapidly with frequency. This is expected as linear response dominates at high frequencies in almost all systems. For BLG, the nonlinear response at single frequency is about five times stronger than frequency tripled terms.

Figure 11 shows the temperature dependent nonlinear conductance at a field of 600 V/cm and at a frequency of 1 THz. At low temperature, the nonlinear conductance exceeds the linear conductance. The σ₃(ω) is greater than the linear conductance in the whole temperature regime. The all important σ₃(3ω) stays as the same as the linear conductance even at room temperatures.

Fig.10 Frequency dependence of the third-order nonlinear optical conductivities at zero and room temperatures [69]. The electric field strength is 1000 V/cm

Full size|PPT slide

Fig.11 Linear and nonlinaer conductances vs. temperature for frequency [69] of 1 THz. The electric field is 600 V/cm

Full size|PPT slide

There are two critical electric field strengths, E_c(ω) and E_c(3ω), at which the nonlinear response equals the linear response. Figure 12 shows the frequency dependence of the critical fields at zero and room temperature. Within the frequency range 0-3 THz, the critical fields are well within the field strength achievable in a laboratory. At f = ω/2π = 1 THz, E_c(ω) = 1100 V/cm at zero temperature and 800 V/cm at room temperature, and E_c(3ω) = 4700 V/cm at zero temperature and 3000 V/cm at room temperature. This is comparable to the nonlinear effect in single layer grapheme [ 77]. This result suggests that interlayer coupling and doubling the carrier numbers in BLG do not reduce the nonlinear effect. If this trend is maintained up to a few layers, the potential for developing graphene-based nonlinear devices can be significantly expanded. The small cusp at low frequency is due to a singularity in the density of states [ 94], which gives rise to a large value of linear current.

Fig.12 Frequency dependent critical fields at zero and room temperatures [69]

Full size|PPT slide

In Fig. 13, we present the temperature dependence of the critical field. The rapid decrease in the critical field at low temperature is mainly due to the decrease in linear current. The sole contribution to the linear current is frequencies, the contribution to the total nonlinear current from the central Dirac points and the three satellites Dirac points can be separated. We found that for both σ₃(ω) and σ₃(3ω), the contribution from the central Dirac point is less than 10% while each satellite Dirac point contributes around 30% of the total nonlinear current. This is a clear indication on the connection between the trigonal warping and nonlinear optical processes in BLG since the existence of the satellite Dirac point is a unique signature of the low energy trigonal warping effect.

Fig.13 Temperature dependent critical fields for frequency [69] of 1 THz

Full size|PPT slide

Semihydrogenated graphene

We now investigate the optical response in SHG in THz frequency regime. In general, for systems with a finite gap, the linear response, or one photon process for frequency below the bandgap Δ, is forbidden. However, multiphoton processes can still occur for frequencies below the gap. The strength of such nonlinear response is usually very weak. We found that the opening of a band gap at the Dirac point leads to a very strong nonlinear response below the gap. In fact, the low frequency nonlinear conductance can be as strong as the universal conductance in intrinsic graphene under a rather moderate electric field of the order of 10³ V/cm. This result is particularly useful for developing applications in nonlinear optics and nonlinear photonics since the linear process is fully suppressed in this frequency regime. Furthermore, we found that the nonlinear optical response at the onset frequency of the nonlinear subgap conductivity peak is universally enhanced by a factor of 31/13 ≈ 2.38 regardless the value of Δ. This suggests that this enhancement is related to the topological changes in the energy band structure of Dirac quasiparticles when a bandgap is created.

Recursion equations and linear optical current density

The interband optical conductivity is calculated by recursively solving the n-photon coupled spinor components in the presence of an external electric field. We first construct the recursion equations for the n-photon coupled spinor components. In the tight-binding approximation, the Hamiltonian for SHG under a time-dependent electric field along the x-axis Ee^-iωt can be written as

(63) H=\left(\begin{array}{cc}-\mathrm{\Delta }/2& v_0\left(p_{-}+eA\right)\\ v_0\left(p_{+}+eA\right)& \mathrm{\Delta }/2\end{array}\right),

where p_± = p_x ± ip_y, and A={\displaystyle \frac{E}{i\omega }}{e}^{i\omega t} . The on-site energies of the A-sublattice and B sublattice are -Δ/2 and Δ/2, respectively. This creates a bandgap opening of Δ at the Dirac point. The quasiparticle is equivalent to a massive Dirac fermion in this case.

We now write the two-component wave function in terms of two spinor components a_n(p), and b_n(p):

(64) \psi (\mathit{p})={\displaystyle \sum _{n=\mathrm{0}}^{\infty }\left(\begin{array}{l}{a}_n(\mathit{p})\\ b_n(\mathit{p})\end{array}\right)}{e}^{i(n\omega -ϵ/\mathit{ћ})t}.

By solving the Schrodinger’s equation, we obtain the following coupled recursion equations:

n\omega \mathit{ћ}(n\omega \mathit{ћ}-\mathrm{2}ϵ)a_n={\displaystyle \frac{eE{v}_{\mathrm{0}}}{i\omega }}\left[\left(ϵ-n\omega \mathit{ћ}-{\displaystyle \frac{\mathrm{\Delta }}{\mathrm{2}}}\right)b_{n-\mathrm{1}}+v_{\mathrm{0}}{p}_{-}{a}_{n-\mathrm{1}}\right],

(65) n\omega \mathit{ћ}(n\omega \mathit{ћ}-\mathrm{2}ϵ)b_n={\displaystyle \frac{eE{v}_{\mathrm{0}}}{i\omega }}\left[\left(ϵ-n\omega \mathit{ћ}+{\displaystyle \frac{\mathrm{\Delta }}{\mathrm{2}}}\right)a_{n-\mathrm{1}}+v_{\mathrm{0}}{p}_{+}{b}_{n-\mathrm{1}}\right].

Any higher order spinor components can hence be recursively constructed.

Linear and nonlinear optical responses

As discussed in previous sections, the optical current operator and the current density can be constructed using Eqs. (56) and (57). We found that the linear optical conductivity, σ₁, in SI unit is given by

(66) {\sigma }_{\mathrm{1}}(\omega )={\displaystyle \frac{e^{\mathrm{2}}}{\mathrm{4}\mathit{ћ}}}\left[\mathrm{1}+{\displaystyle \frac{{\mathrm{\Delta }}^{\mathrm{2}}}{{\omega }^{\mathrm{2}}{\mathit{ћ}}^{\mathrm{2}}}}\right]\tanh \left({\displaystyle \frac{\omega \mathit{ћ}}{\mathrm{4}{k}_{\mathrm{B}}T}}\right)\mathrm{\Theta }(\omega \mathit{ћ}-\mathrm{\Delta }).

Note that the step-function \mathrm{\Theta }(\omega \mathit{ћ}-\mathrm{\Delta }) function forbids any linear optical process to occur in the subgap regime [ 96]. The third-order nonlinear optical conductivities

(67) {\sigma }_{\mathrm{3}}(\omega )={\sigma }_{\mathrm{0}}{\displaystyle \frac{e^{\mathrm{2}}{E}^{\mathrm{2}}{v}_{\mathrm{0}}^{\mathrm{2}}}{{\omega }^{\mathrm{3}}\mathit{ћ}\left(\omega \mathit{ћ}+{\displaystyle \frac{\mathrm{\Delta }}{\mathrm{2}}}\right)}}\left[\mathrm{2}+{\displaystyle \frac{\mathrm{\Delta }}{\omega \mathit{ћ}}}+{\displaystyle \frac{{\mathrm{\Delta }}^{\mathrm{2}}}{{(\omega \mathit{ћ})}^{\mathrm{2}}}}+{\displaystyle \frac{{\mathrm{\Delta }}^{\mathrm{3}}}{\mathrm{2}{(\omega \mathit{ћ})}^{\mathrm{3}}}}-{\displaystyle \frac{\mathrm{3}{(\mathrm{\Delta })}^{\mathrm{4}}}{\mathrm{8}{(\omega \mathit{ћ})}^{\mathrm{4}}}}-{\displaystyle \frac{\mathrm{3}{\mathrm{\Delta }}^{\mathrm{5}}}{\mathrm{16}{(\omega \mathit{ћ})}^{\mathrm{5}}}}\right]\times \tanh \left({\displaystyle \frac{\omega \mathit{ћ}}{\mathrm{2}{k}_{\mathrm{B}}T}}\right)\mathrm{\Theta }(\omega \mathit{ћ}-\mathrm{\Delta }),

and

(68) {\sigma }_{\mathrm{3}}(\mathrm{3}\omega )={\sigma }_{\mathrm{0}}{\displaystyle \frac{e^{\mathrm{2}}{E}^{\mathrm{2}}{v}_{\mathrm{0}}^{\mathrm{2}}}{{\omega }^{\mathrm{4}}{\mathit{ћ}}^{\mathrm{2}}}}\left[X_{\mathrm{1}}\tanh \left({\displaystyle \frac{\omega \mathit{ћ}}{\mathrm{4}{k}_{\mathrm{B}}T}}\right)+X_{\mathrm{2}}\tanh \left({\displaystyle \frac{\omega \mathit{ћ}}{\mathrm{2}{k}_{\mathrm{B}}T}}\right)+X_{\mathrm{3}}\tanh \left({\displaystyle \frac{\mathrm{3}\omega \mathit{ћ}}{\mathrm{4}{k}_{\mathrm{B}}T}}\right)\right]\mathrm{\Theta }(\mathrm{3}\omega \mathit{ћ}-\mathrm{\Delta }),

where σ₀ = e²/4ħ and

(69) \begin{array}{l}{X}_{\mathrm{1}}\equiv -{\displaystyle \frac{\mathrm{1}}{\mathrm{48}}}\left[\mathrm{13}+{\displaystyle \frac{\mathrm{2}{\mathrm{\Delta }}^{\mathrm{2}}}{{(\omega \mathit{ћ})}^{\mathrm{2}}}}+{\displaystyle \frac{{\mathrm{\Delta }}^{\mathrm{4}}}{{(\omega \mathit{ћ})}^{\mathrm{4}}}}\right],\hfill \\ X_{\mathrm{2}}\equiv {\displaystyle \frac{\mathrm{1}}{\mathrm{3}}}\left[\mathrm{2}-{\displaystyle \frac{{\mathrm{\Delta }}^{\mathrm{2}}}{{(\omega \mathit{ћ})}^{\mathrm{2}}}}+{\displaystyle \frac{{\mathrm{\Delta }}^{\mathrm{4}}}{\mathrm{8}{(\omega \mathit{ћ})}^{\mathrm{4}}}}\right],\hfill \\ X_{\mathrm{3}}\equiv -{\displaystyle \frac{\mathrm{1}}{\mathrm{48}}}\left[\mathrm{45}-{\displaystyle \frac{\mathrm{14}{\mathrm{\Delta }}^{\mathrm{2}}}{{(\omega \mathit{ћ})}^{\mathrm{2}}}}+{\displaystyle \frac{{\mathrm{\Delta }}^{\mathrm{4}}}{{(\omega \mathit{ћ})}^{\mathrm{4}}}}\right].\hfill \end{array}

When Δ→ 0, the above equations reduces to the usual graphene conductivities [ 77], i.e.,

(70) {\sigma }_{\mathrm{3}}(\omega )\to {\sigma }_{\mathrm{0}}{\displaystyle \frac{e^{\mathrm{2}}{E}^{\mathrm{2}}{v}_{\mathrm{0}}^{\mathrm{2}}}{{\omega }^{\mathrm{4}}{\mathit{ћ}}^{\mathrm{2}}}}\tanh \left({\displaystyle \frac{\omega \mathit{ћ}}{\mathrm{2}{k}_{\mathrm{B}}T}}\right)\times \mathrm{2},

(71) {\sigma }_{\mathrm{3}}(\mathrm{3}\omega )={\sigma }_{\mathrm{0}}{\displaystyle \frac{e^{\mathrm{2}}{E}^{\mathrm{2}}{v}_{\mathrm{0}}^{\mathrm{2}}}{{\omega }^{\mathrm{4}}{\mathit{ћ}}^{\mathrm{2}}}}\left[-{\displaystyle \frac{\mathrm{13}}{\mathrm{48}}}\tanh \left({\displaystyle \frac{\omega \mathit{ћ}}{\mathrm{4}{k}_{\mathrm{B}}T}}\right)+{\displaystyle \frac{\mathrm{2}}{\mathrm{3}}}\tanh \left({\displaystyle \frac{\omega \mathit{ћ}}{\mathrm{2}{k}_{\mathrm{B}}T}}\right)-{\displaystyle \frac{\mathrm{45}}{\mathrm{48}}}\tanh \left({\displaystyle \frac{\mathrm{3}\omega \mathit{ћ}}{\mathrm{4}{k}_{\mathrm{B}}T}}\right)\right].

In Fig. 14, we plot the optical conductance versus frequency for a typical value of Δ = 0.03 eV . The on-site energy due to semihydrogenation removes the universal conductance. For ħω < Δ, the linear conductance is zero for any temperature by the virtue of energy conservation. The third-order current at single frequency, σ₃(3ω), is also zero for ħω < Δ. The triple-frequency third-order term σ₃(3ω)persists to a low frequency of ħω = Δ/3. The nonlinear effect in SHG is unique in that the response peak of the linear term and frequency tripled term is well separated by δħω = 2Δ/3. This provides a useful mechanism for two-color excitation and detection, one color is associated with the linear response and the other is associated with the nonlinear response. The relative intensities of the two colors can be tuned with the electric field. At a rather moderate electric field of 3600 V/cm, the magnitude of two peaks is roughly the same at 77 K. At room temperature, the peak in linear conductance disappears while the nonlinear conductance still exhibits a resonance.

Fig.14 Frequency dependent optical conductance in the low frequency regime for two temperatures. The electric field is 3600 V/cm. The absorption edge for the frequency tripled response is shifted to Δ/3. The inset is a schematic showing different optical processes [70]

Full size|PPT slide

We now discuss one interesting behavior of the nonlinear response peak. We compare the σ₃(3ω) at the optical response peak with that of the gapless graphene at the same frequency. The ratio at T = 0 K is given by σ₃(3ω)_Δ/σ₃(3ω) = 31/13 ≈ 2.38. This 2.38 times enhanced optical response at the onset frequency of σ₃(3ω) is universal for any value of Δ. This suggests that the 2.38 enhancement is related to the topological changes in the band structure of the Dirac fermion when a bandgap created.

At the frequencies close to the energy gap, the onset linear conductance is twice the universal conductance [ 96] σ_1c = 2σ₀. The onset triple-frequency nonlinear conductivity σ₃(3ω) = σ_1c at ħω = Δ/3 requires an applied field of E = 3600 V/cm. This is a rather weak field for typical experimental conditions. On the other hand, the onset single-frequency nonlinear conductivity σ₃(ω) = σ_1c at ħω = Δ requires an electric field of around three times greater. Therefore, the potential of using the frequency tripled nonlinear effect in the frequency below the gap is very significant. The electric field required for σ₃(ω) = σ_1c at the vertical absorption edge can be determined,

(72) E_{\mathrm{c}}(\mathit{ћ}\omega =\mathrm{\Delta }/\mathrm{3})={\displaystyle \frac{{\mathrm{\Delta }}^{\mathrm{2}}}{\mathrm{9}e\mathit{ћ}{v}_F}}{\left[{\displaystyle \frac{\mathrm{24}}{\mathrm{56}\tanh \left({\displaystyle \frac{\mathrm{\Delta }}{\mathrm{12}{k}_{\mathrm{B}}T}}\right)-\mathrm{25}\tanh \left({\displaystyle \frac{\mathrm{\Delta }}{\mathrm{6}{k}_{\mathrm{B}}T}}\right)}}\right]}^{\mathrm{1}/\mathrm{2}}.

Fig.15 Frequency dependence of the critical field E(3ω) for SHG and pure (gapless) graphene [70]. The inset shows the reduction of the critical field in SHG. Note that there exists a cut-off frequency f_c =ω_c/2π = Δ/3h≈ 2.4 THz since σ₃(3ω) = 0 at frequency smaller than f_c

Full size|PPT slide

In Fig. 15, we show the frequency dependence of the critical field at which σ₃(3ω)/σ_1c = 1. This field measures the nonlinearity of the system at a given frequency. In the entire low frequency regime, Δ/3 < ħω < Δ, we found that the critical field for SHG is smaller than that in pure graphene by around 10%-40%. This indicates that SHG is a strong nonlinear system at low frequencies and low temperatures. The reason for this is that the density of states near the band edge has a van Hove-like singularity, D(ϵ) ≈ϵ^1/2. This is qualitatively different from the case of normal two-dimensional semiconductors. In normal semiconductors, the energy dispersion near the band edge is parabolic and the density of states is constant. Here in SHG the large density of state near the band edge leads to a strong nonlinear effect.

In Fig. 16, we show the temperature dependence of the critical field E_c(3ω) at two different frequencies. At low temperature E_c(3ω) is nearly constant and is smaller than that required in pure graphene. At high temperature, E_c(3ω) in SHG is larger than that required in pure graphene. As temperature increases, the Van Hove singularity becomes weaker and the critical field increases. At high temperature, E_c(3ω) increases with temperature as E_c(3ω) ≈T^1/2. It should be pointed out that a high critical field in SHG at room temperature will not remove the key property of two-color optical response in SHG. In pure graphene, the response maximum of the linear term and frequency tripled term is not resolved.

Fig.16 Temperature dependence of the critical field [70] at two different frequencies of 2.4 and 5 THz

Full size|PPT slide

The nonlinear effect reported is more general than that in SHG. Many effects can lead to a finite gap in the Dirac point in graphene. For example, the spin-orbit coupling can result in a gap of the size of Δ≈ 0.2 meV. This is a very small gap but will produce qualitatively the same nonlinear effect as in SHG. Impurity scattering included gap which is also in the form of Eq. (72). The critical field mentioned earlier for σ₃(3ω) is proportional to Δ². Therefore, in general, the smaller the gap, the weaker the critical field at the onset frequency of the nonlinear response peak. If the gap can be controlled by external means, then the distance between the two peaks also becomes tunable. However, smaller gaps will result in a smaller distance between the peaks of linear response and frequency tripled response and the linear and nonlinear response peaks becomes less resolved.

Graphene superlattice

We now study the optical response of a Kronig-Penney type graphene superlattice. In this structure, the band structure of the massless Dirac fermion is no longer symetrical in k-space. The effect of anisotropy on the optical response in THz frequency regime is investigated. It is found that the optical absorption, both linear and nonlinear response, are universally enhanced by the anisotropy when the external field aligns with the superlattice periodicity. Since both linear and nonlinear response are enhanced by the same amount, the optical nonlinearity (i.e., the relative magnitude between linear and nonlinear responses) is unexpectedly preserved regardless how strong the band structure anisotropy is. The enhanced optical absorption and the preserved optical nonlinearity reveal that anisotropy has transformed graphene superlattice into a stronger nonlinear material which produces larger nonlinear optical current than isotropic case under the same critical electric field strength. Such enhanced optical absorption and well-preserved optical nonlinearity also occurs in gapped graphene in which the quasiparticle is in the form of massive Dirac fermions. The anisotropic massive Dirac femrion is a bizzare quasiparticle not only with non-uniform ‘light speed’ but also non-uniform mass dependent of the propagation direction. The results suggest that the enhanced electron-single-photon and electron-multiple-photon couplings is a universal feature of relativistic Dirac fermions of both massless or massive types, and the band structure isotropy is not a pre-requisite for the strong optical nonlinearity in graphene.

Recursion equations of anisotropic massless Dirac fermion

In a graphene superlattice created by applying a Kronig-Penney potential [ 48, 49], the K-point electrons no longer travels with uniform v_F in all direction. Instead, the group velocity in the direction perpendicular to the periodicity of the Kronig-Penney potential is reduced by a factor of λ dependent on the strength and periodicity of the potential. The band structure of the superlattice is no longer circular, but is deformed to an elliptic cone. Such quasiparticle nature is analog to a massless Dirac fermion traveling in anisotropic spacetime [ 48, 49, 97].

Fig.17 Band structure of graphene superlattice (inset). In the p_x-p_y plan, the Dirac cone is elongated elliptically in the y-direction. L, w and U are the superlattice periodicity, potential width and potential height, respectively, of the Kronig-Penney type graphene superlattice [71]

Full size|PPT slide

In topological insulator (TI), the quasiparticle residing in its surface state is also in the massless Dirac form with Fermi velocity approximately half of the graphene [ 63]. Interestingly, the anisotropic massless Dirac fermion can also be found in the (2, 2, 1) side-surface state of Bi₂Se₃ TI with a rather strong anisotropy of v_x = 3.1 × 10⁵ m/s and v_y = 1.4 × 10⁵ m/s [ 98]. In a Bi square net of SrMnBi₂ TI, highly anisotropic Fermi velocity differs by a factor of 8 was experimentally observed [ 99].

The anisotropy can be modeled by defining an anisotropy parameter, λ, which modifies the y-direction group velocity by v_y = λv_F where v_F = 10⁶ m/s is the Fermi velocity, and the anisotropy parameter is continuously tunable, 0 ≤λ≤ 1, by varying the superlattice periodicity L, potential width w and potential height U (Fig. 17) [ 48, 49]. The graphene superlattice Hamiltonian is written as \widehat{\mathscr{H}}={\sigma }_x{p}_x+\lambda {\sigma }_y{p}_y , where the λ term has created the desired anisotropy in y-direction. The energy dispersion is given as {ϵ}_s(p,\theta )=s{v}_{F}p\sqrt{{\cos }^{\mathrm{2}}\theta +{\lambda }^{\mathrm{2}}{\sin }^{\mathrm{2}}\theta } where s = ± 1 denotes electron and hole state. The group velocity in θ-direction is given as v=v_F\sqrt{{\cos }^{\mathrm{2}}\theta +{\lambda }^{\mathrm{2}}{\sin }^{\mathrm{2}}\theta } , which gives the expected x- and y-components of v(\theta =\mathrm{0})=v_F=v_x and v(\theta =\mathrm{\pi }/\mathrm{2})=\lambda v_F=v_y . The eigenfunction is given as {\psi }_{\mathrm{0}}(s,p)={\displaystyle \frac{\mathrm{1}}{\sqrt{\mathrm{2}}}}{(\mathrm{1},v_F(p_x+i\lambda p_y)/{ϵ}_s)}^{\mathrm{T}} , where T denotes transpose. The band structure is plotted in Fig. 17. It can be seen that due to the reduced group velocity in y-direction, the conic Dirac cone is elongated in y-direction, forming an anisotropic elliptic Dirac cone. When an external field E = E₀e^iωt is applied along the x-direction, the quasiparticle is minimally coupled to the photon according to p_x→p_x-eA, where A = -∂E/∂t. The Hamiltonian is then given by

(73) \widehat{H}=v_F\left[\begin{array}{ll}\mathrm{0}& p_x-\lambda i{p}_y-eA\\ p_x+\lambda i{p}_y-eA& \mathrm{0}\end{array}\right].

The energy dispersion is given as {ϵ}_s(p,\varphi )=s{v}_{F}p\sqrt{{\cos }^{\mathrm{2}}\theta +{\lambda }^{\mathrm{2}}{\sin }^{\mathrm{2}}\theta } , where s = ± 1. The single electron eigenstate is given as

(74) {\psi }_{\mathrm{0}}(s,p)={\displaystyle \frac{\mathrm{1}}{\sqrt{\mathrm{2}}}}\left[\begin{array}{l}\mathrm{1}\\ s{v}_F{\displaystyle \frac{p_x+i\lambda p_y}{ϵ}}\end{array}\right],

where ϵ=\left|{ϵ}_s\right| . The wave function in the presents of an external electric field is written as

(75) {\psi }_n(p)={\displaystyle \sum _n\left[\begin{array}{l}{a}_n\\ b_n\end{array}\right]}{e}^{i\left(\frac{ϵs}{\mathit{ћ}}-n\omega \right)t}.

Similarly, the Schrodinger’s equation H\psi =i\mathit{ћ}\partial \psi /\partial t can then be solved to obtain a pair of recursive equations for the spinor components a_n and b_n:

(76) (ϵ-n\mathit{ћ}\omega )a^n=v_F{\tilde{p}}_{-}{b}_n+{\displaystyle \frac{eE{v}_F}{i\omega }}{b}_{n-\mathrm{1}},

(77) (ϵ-n\mathit{ћ}\omega )b^n=v_F{\tilde{p}}_{+}{a}_n+{\displaystyle \frac{eE{v}_F}{i\omega }}{a}_{n-\mathrm{1}},

where {\tilde{p}}_{\pm }=p_x\pm i\lambda p_y .

Linear and nonlinear optical responses

The linear optical response is found to be σ₁(ω) = (σ₀/λ) tanh (ħω/2), where the spin and valley degeneracies have been included (a factor of 4). The second order and third order spinors can be constructed similarly using the recursion equation, Eq. (76). Following exactly the same procedures, we found that

(78) {\sigma }_{\mathrm{3}}(\omega )=\mathrm{2}{\sigma }_{\mathrm{0}}{\displaystyle \frac{E^{\mathrm{2}}{v}_F^{\mathrm{2}}{e}^{\mathrm{2}}}{{\mathit{ћ}}^{\mathrm{2}}{\omega }^{\mathrm{4}}}}{\displaystyle \frac{\mathrm{1}}{\mathrm{2}\mathrm{\pi }}}{\displaystyle \int {\displaystyle \frac{\mathrm{d}\varphi }{R_{+}}}\left(\mathrm{1}-{\displaystyle \frac{R_{-}}{R_{+}}}\right)}\tanh \left({\displaystyle \frac{\mathit{ћ}\omega }{k_{\mathrm{B}}T}}\right),

(79) {\sigma }_{\mathrm{3}}(\omega )=\mathrm{2}{\displaystyle \frac{{\sigma }_{\mathrm{0}}}{\lambda }}{\displaystyle \frac{E^{\mathrm{2}}{v}_F^{\mathrm{2}}{e}^{\mathrm{2}}}{{\mathit{ћ}}^{\mathrm{2}}{\omega }^{\mathrm{4}}}}\tanh \left({\displaystyle \frac{\mathit{ћ}\omega }{k_{\mathrm{B}}T}}\right).

For the third-order triple-frequency (TF) conductivity, we obtain

(80) {\sigma }_{\mathrm{3}}(\mathrm{3}\omega )={\displaystyle \frac{{\sigma }_{\mathrm{0}}}{\lambda }}{\displaystyle \frac{E^{\mathrm{2}}{v}_F^{\mathrm{2}}{e}^{\mathrm{2}}}{{\mathit{ћ}}^{\mathrm{2}}{\omega }^{\mathrm{4}}}}\left({\displaystyle \frac{\mathrm{13}}{\mathrm{48}}}\tanh \left({\displaystyle \frac{\mathit{ћ}\omega }{\mathrm{2}{k}_{\mathrm{B}}T}}\right)-{\displaystyle \frac{\mathrm{2}}{\mathrm{3}}}\tanh \left({\displaystyle \frac{\mathit{ћ}\omega }{k_{\mathrm{B}}T}}\right)+{\displaystyle \frac{\mathrm{15}}{\mathrm{16}}}\tanh \left({\displaystyle \frac{\mathrm{3}\mathit{ћ}\omega }{\mathrm{2}{k}_{\mathrm{B}}T}}\right)\right).

(80)

In Fig. 18, the nonlinear optical conductivities at different band anisotropy λ is shown. We see that σ₁(ω), σ₃(ω) and σ₃(3ω) are all universally enhanced by a factor of 1/λ, in comparison with that of the isotropic case [ 77]. For λ = 0.1, which can be achieved by applying spatial period of L≈ 20 nm, potential width of w = 10 nm and potential height of U = 0.3 eV, the total optical absorption is enhanced by 10 times. In the extremely anisotropic case of λ = 0.01, which can be achieved by L≈ 25 nm, w = 10 nm and U = 0.3 eV [ 48], 100 times amplification is achieved. The 1/λ enhanced optical absorption is quite a surprising result. Intuitively, one might expect a reduced optical response in the anisotropic case since the y-component of the group velocity v_y = λv_F is reduced by a factor of λ and the resulting ‘slower’ charge carrier should degrade the optical current. This is however not the complete picture since E is directed along x-direction and the x-directional optical response is only minimally affected by the reduced y-directional group velocity v_y = λv_F. On the other hand, when λ < 1, the p_y components in a equi-energy slice actually becomes larger in comparison to the isotropic Fermi velocity case because of the smaller slope (or equivalently the reduced v_y) in y-direction (see Fig. 1). The overall larger momentum of the charge carrier across an equi-energy surface is the underlying reason of the anisotropy-induced enhancement of the interband optical absorption in the Kronig-Penney type graphene superlattice. As the anisotropy increases, i.e., λ→ 0, the band structure becomes more y-directionally elongated across an equienergy surface and this generates the 1/λ dependence.

The critical field strength E_c remains the same regardless the level of anisotropy since both linear and nonlinear response is enhanced by the same factor of 1/λ. Therefore, just like normal graphene, graphene superlattice is also an exceptionally strong nonlinear material with E_c≈ 10³ V/cm for up to room temperature at f = 1 THz. The strong optical nonlinearity observed in normal graphene and graphene superlattice is a general feature of the relativistic behavior of the quasiparticle. The band structure isotropy is not necessarily required to achieve the strong optical nonlinearity. As long as the quasiparticle energy dispersion maintains its linear form, the strong optical nonlinearity is always guaranteed and is well protected from any band structure anisotropy. The total integrated optical absorption is given as \mathrm{\Sigma }(\lambda )=\int \sigma (\omega ,\lambda )\mathrm{d}\omega and it can be immediately seen that the total nonlinear absorption is increased by a factor 1/λ for all THz frequency regime as shown in Fig. 18. Although graphene superlattice is equally advantageous as normal graphene in terms of the smallness of E_c, the 1/λ increased total response indicates that the nonlinear optical current output of graphene superlattice is still larger than that of the normal graphene at a given electric field strength. This suggests the improved THz photon detection and THz frequency up-conversion in graphene superlattice which are potentially useful in the development of graphene-based THz optical device. Finally, we briefly discuss the optical response of gapped graphene with anisotropic band structure. We found that the linear and nonlinear optical conductivity is in the same form as Eqs. (66), (70) and (70) multiplied by a factor of 1/λ. The band anisotropy enhanced subgap triple-frequency conductivity is plotted in Fig. (19). We conclude that the 1/λ enhancement is universal in both gapped and gapless cases in the presence of band anisotropy.

Fig.18 Frequency dependence [71] of σ₃(3ω) at E = 1000 V/cm and T = 300 K

Full size|PPT slide

Discussion

Terahertz photon mixing effect of gapless and gapped single layer graphene

In graphene, the nonlinear effect is approximately inversely proportional to the Fermi-level and grows rapidly with temperature up to room temperature. The critical electric field required to generate nonlinear effect comparable to linear effect is in a rather moderate value of 10⁴ V/cm even in room temperature. The optical response of single layer graphene under strong-field condition exhibits the following interesting behavior: (1) the linear and second-order nonlinear responses are well protected from external field due to the unique Dirac fermions dynamics and the preservation of crystal inversion symmetry; (2) the third-order nonlinear optical response is enhanced by three distinct mechanisms: (i) third-order response is intrinsically proportional to E³; (ii) strong-field induces Dirac fermion population redistribution creates an additional contribution to third-order response; and (iii) the nonequilibrium heating raises the carrier temperature to T_hot> T_lattice and further enhances the nonlinear current. The strong and temperature-robust nonlinear optical nonlinearity suggests that graphene can potentially be an excellent candidate in nonlinear photon-mixing applications. In gapped graphene, the nonlinear optical response is strongly influenced by the bandgap value and the temperature. To maximize the nonlinear optical response of a gapped graphene-based photon mixer, the bandgap value and the operating temperature has to be carefully chosen.

Fig.19 Anisotropic gapped graphene frequency-tripling conductivity [71] at T = 300 K and E = 3400 V/cm and ∆ = 0.03 eV

Full size|PPT slide

We point out several experiments which can potentially be used to verify our theoretical calculations. Several experimental works emphasizing the visible and near-infrared nonlinear optical response of graphene has been reported recently [ 73, 100, 101]. Multiple-photon absorption/transmission experiments [ 100, 101] can be repeated in the THz regime to qualitatively estimate the optical nonlinearity of graphene. The nonlinear wave-mixing effect can be more accurately quantified by irradiating an graphene sample with two THz waves of frequencies ω₁ and ω₂, and selectively filtering the outgoing waves to determine the strength of the mixed wave (2ω₁±ω₂) [ 73]. The temperature dependence of the wave-mixing effect can be probed by performing these experiments under controlled temperature condition.

Bilayer graphene

We found that BLG is a rather strong nonlinear material. This nonlinear effect is robust from low to room temperatures. The frequency tripling nonlinear term is comparable to the linear term in the THz frequency regime. This suggests that BLG has a potential in a THz emitter/detector at frequencies, which are traditionally difficult to obtain by using an existing emitter at one-third the frequency. We now briefly present on the role of phonon excitation in BLG. In the temperature range of up to room temperature, the dominant electron-phonon coupling is via longitudinal acoustic (LA) phonons since either the couplings to other graphene lattice phonon modes are too weak or the energy scales of these optical phonon modes are far too high [ 102]. The velocity of the LA phonon is around 2 × 10⁴ m/s [ 102]. Under an electric field around 1000 V/cm with a frequency of 1 THz, the energy of the photoexcited electron is around 1 THz. These electrons are located very close to the Dirac point, and the electron velocity is around 0.6 × 10⁶ m/s. In the absence of other disorders and due to the energy conversation, the probability of single phonon emission is negligible. The multiple phonon excitations are possible but the probability is also very low due to the high order electron-phonon coupling. Therefore, in the absence of disorders, we do not expect that phonon excitation will play a significant role in altering the nonlinear electrical current in this energy regime. In the presence of impurities, electron-phonon coupling in single layer graphene can be enhanced by disorder-assisted supercollision process [ 103- 105].

In conclusion, we have shown that BLG exhibits a strong nonlinear effect in the THz to FIR regime under an electric field of around 10³ V/cm. In particular, a moderate field can induce the frequency tripling term at room temperature. This suggests a potential for developing graphene-based optics and photonics applications.

Semihydrogenated graphene

It is found that SHG with a bandgap in its Dirac point exhibits strong nonlinear optical response at frequency range of Δ/3 < ωħ < Δ. In this frequency range, the optical response is solely contributed by three-photon nonlinear process and hence has a zero critical field. The nonlinear response peak and the linear response peaks are well-separated giving rise to a two-color characteristic. Furthermore, the triple-frequency nonlinear optical response is universally enhanced by a factor of 31/13 ≈ 2.38, suggesting a topological origin due to the bandgap opening in the Dirac fermion energy spectrum.

Graphene superlattice

The anisotropic Dirac fermion in the graphene superlattice tunes up the total optical conductance while maintaining the same critical electric field. This also occurs in anisotropic graphene with a gap. Furthermore, the optical nonlinearity is perfectly protected from band anisotropy while the total optical responses, including both linear and nonlinear processes, are enhanced by a factor of 1/λ. Since λ is dependent on the superlattice parameters, a graphene superlattice can potentially be used as a tunable THz source/detector. It can be noticed that the fifth order nonlinear term can potentially play a role in the optical nonlinearity of a superlattice structure. The nth-order conductance is proportional to a dimensionless parameter Z = (eEv_F /ħω²)^n-1 and an overlap integral of eigenstates of different orders ＜ {\varphi }_{n-m}|{\varphi }_m ＞ . Because the overlap integral decreases very rapidly with n, the third-order nonlinear effect persists for Z > 1 while the fifth-order term is negligible. At frequency around 1 THz, the critical field (the field at which the third-order current equals the linear current) is around 2000 V/cm [ 77]. For v_F = 10⁶ m/s, ω = 1 THz, and E = 2000 V/cm, the resulting Z = 50. At this value of Z, the third-order current equals approximately the linear current, but the fifth-order current is about 10^-5 of the linear current, totally negligible [ 106]. Finally, as a weak sinusoidal term can be added to a graphene via holographic illumination [ 107] or by patterning the substrate, possible experimental verification of the results could be performed with direct measurement of the optical conductivity of such a system.

Conclusion

In conclusion, we review and discuss the nonlinear optical response of graphene and its related sister-structures in the THz and FIR frequency regimes. It is found that not only single layer graphene exhibits strong optical nonlinearity, stacking up graphene layers, bandgap opening at the Dirac points (such as SHG), and the construction of graphene superlattice via electrostatic gating also render the material with strong optical nonlinearity. Finally, we propose future experiments on graphene structures to be performed in order to verify our theoretical results.

References

1 Novoselov K S, Geim A K, Morozov S V, Jiang D, Katsnelson M I, Grigorieva I V, Dubonos S V, Firsov A A. Two-dimensional gas of massless Dirac fermions in graphene. Nature, 2005, 438(7065): 197-200

PMID

2 Novoselov K S, Geim A K, Morozov S V, Jiang D, Zhang Y, Dubonos S V, Grigorieva I V, Firsov A A. Electric field effect in atomically thin carbon films. Science, 2004, 306(5696): 666-669

PMID

3 Wallace P R. The band theory of graphite. Physical Review, 1947, 71(9): 622-634

4 Katsnelson M I, Novoselov K S, Geim A K. Chiral tunnelling and the Klein paradox in graphene. Nature Physics, 2006, 2(9): 620-625

5 Klein O. Die reflexion von elektronen an einem potentialsprung nach der relativistischen dynamik von Dirac. Zeitschrift fur Physik, 1929, 53(3-4): 157-165

6 Young A F, Kim P. Quantum interference and Klein tunnelling in graphene heterojunctions. Nature Physics, 2009, 5(3): 222-226

7 Stander N, Huard B, Goldhaber-Gordon D. Evidence for Klein tunneling in graphene p-n junctions. Physical Review Letters, 2009, 102(2): 026807

PMID

8 Wright A R, Cao J C, Zhang C. Enhanced optical conductivity of bilayer graphene nanoribbons in the terahertz regime. Physical Review Letters, 2009, 103(20): 207401

PMID

9 Wang X L, Dou S X, Zhang C. Zero-gap materials for future spintronics, electronics and optics. NPG Asia Materials, 2010, 2(1): 31-38

10 Liu J, Ma Z, Wright A R, Zhang C. Orbital magnetization of graphene and graphene nanoribbons. Journal of Applied Physics, 2008, 103(10): 103711

11 YuD C, Lupton E M, Gao H J, Zhang C, Liu F. A unified geometric rule for designing nanomagnetism in graphene. Nano Research, 2008, 1(6): 497-501

12 Cai J Z, Lu L, Kong W J, Zhu H W, Zhang C, Wei B Q, Wu D H, Liu F. Pressure-induced transition in magnetoresistance of single-walled carbon nanotubes. Physical Review Letters, 2006, 97(2): 026402

PMID

13 Bolotin K I, Sikes K J, Jiang Z, Klima M, Fudenberg G, Hone J, Kim P, Stormer H L. Ultrahigh electron mobility in suspended graphene. Solid State Communications, 2008, 146(9-10): 351-355

14 Chen J H, Jang C, Xiao S, Ishigami M, Fuhrer M S. Intrinsic and extrinsic performance limits of graphene devices on SiO₂. Nature Nanotechnology, 2008, 3(4): 206-209

PMID

15 Geim A K, Novoselov K S. The rise of graphene. Nature Materials, 2007, 6(3): 183-191

PMID

16 Xia F, Farmer D B, Lin Y M, Avouris P. Graphene field-effect transistors with high on/off current ratio and large transport band gap at room temperature. Nano Letters, 2010, 10(2): 715-718

PMID

17 Schwierz F. Graphene transistors. Nature Nanotechnology, 2010, 5(7): 487-496

PMID

18 Zhang Y, Tan Y W, Stormer H L, Kim P. Experimental observation of the quantum Hall effect and Berry’s phase in graphene. Nature, 2005, 438(7065): 201-204

PMID

19 Gusynin V P, Sharapov S G. Unconventional integer quantum Hall effect in graphene. Physical Review Letters, 2005, 95(14): 146801

PMID

20 Novoselov K S, Jiang Z, Zhang Y, Morozov S V, Stormer H L, Zeitler U, Maan J C, Boebinger G S, Kim P, Geim A K. Room-temperature quantum Hall effect in graphene. Science, 2007, 315(5817): 1379

PMID

21 Ziegler K. Minimal conductivity of graphene: nonuniversal values from the Kubo formula. Physical Review B: Condensed Matter and Materials Physics, 2007, 75(23): 233407

22 Herbut I F, Juricic V, Vafek O. Coulomb interaction, ripples, and the minimal conductivity of graphene. Physical Review Letters, 2008, 100(4): 046403

23 Peres N M R, Guinea F, Castro Neto A H. Electronic properties of disordered two-dimensional carbon. Physical Review B: Condensed Matter and Materials Physics, 2006, 73(12): 125411

24 Cserti J, Csordás A, Dávid G. Role of the trigonal warping on the minimal conductivity of bilayer graphene. Physical Review Letters, 2007, 99(6): 066802

PMID

25 Martin J, Akerman N, Ulbricht G, Lohmann T, Smet J H, von Klitzing K, Yacoby A. Observation of electron-hole puddles in graphene using a scanning single-electron transistor. Nature Physics, 2008, 4(2): 144-148

26 Falkovsky L A, Pershoguba S S. Optical far-infrared properties of a graphene monolayer and multilayer. Physical Review B: Condensed Matter and Materials Physics, 2007, 76(15): 153410

27 Zhang C, Chen L, Ma Z. Orientation dependence of the optical spectra in graphene at high frequencies. Physical Review B, 2008, 77(24): 241402

28 Gusynin V P, Sharapov S G, Carbotte J P. Unusual microwave response of dirac quasiparticles in graphene. Physical Review Letters, 2006, 96(25): 256802

PMID

29 Nair R R, Blake P, Grigorenko A N, Novoselov K S, Booth T J, Stauber T, Peres N M R, Geim A K. Fine structure constant defines visual transparency of graphene. Science, 2008, 320(5881): 1308

PMID

30 Li Z Q, Henriksen E A, Jiang Z, Hao Z, Martin M C, Kim P, Stormer H L, Basov D N. Dirac charge dynamics in graphene by infrared spectroscopy. Nature Physics, 2008, 4(7): 532-535

31 Kuzmenko A B, van Heumen E, Carbone F, van der Marel D. Universal optical conductance of graphite. Physical Review Letters, 2008, 100(11): 117401

PMID

32 Rycerz A, Tworzydło J, Beenakker C W J. Valley filter and valley valve in graphene. Nature Physics, 2007, 3(3): 172-175

33 Gunlycke D, White C T. Graphene valley filter using a line defect. Physical Review Letters, 2011, 106(13): 136806

PMID

34 Garcia-Pomar J L, Cortijo A, Nieto-Vesperinas M. Fully valley-polarized electron beams in graphene. Physical Review Letters, 2008, 100(23): 236801

PMID

35 Pereira J M Jr, Peeters F M, Costa Filho R N, Farias G A. Valley polarization due to trigonal warping on tunneling electrons in graphene. Journal of Physics Condensed Matter, 2009, 21(4): 045301

36 Chaves A, Covaci L, Rakhimov K Y, Farias G A, Peeters F M. Wave-packet dynamics and valley filter in strained graphene. Physical Review B: Condensed Matter and Materials Physics, 2010, 82(20): 205430

37 Moldovan D, Masir M R, Covaci L, Peeters F M. Resonant valley filtering of massive Dirac electrons. Physical Review B: Condensed Matter and Materials Physics, 2012, 86(11): 115431

38 Zhai F, Chang K. Valley filtering in graphene with a Dirac gap. Physical Review B: Condensed Matter and Materials Physics, 2012, 85(15): 155415

39 Péterfalvi C G, Oroszlány L, Lambert C J, Cserti J. Intraband electron focusing in bilayer graphene. New Journal of Physics, 2012, 14(6): 063028

40 Majidi L, Zareyan M. Pseudospin polarized quantum transport in monolayer graphene. Physical Review B: Condensed Matter and Materials Physics, 2011, 83(11): 115422

41 San-Jose P, Prada E, McCann E, Schomerus H. Pseudospin valve in bilayer graphene: towards graphene-based pseudospintronics. Physical Review Letters, 2009, 102(24): 247204

PMID

42 Trushin M, Schliemann J. Pseudospin in optical and transport properties of graphene. Physical Review Letters, 2011, 107(15): 156801

PMID

43 Min H, Borghi G, Polini M, MacDonald A H. Pseudospin magnetism in graphene. Physical Review B, 2008, 77(4): 041407

44 Majidi L, Zareyan M. Enhanced Andreev reflection in gapped graphene. Physical Review B: Condensed Matter and Materials Physics, 2012, 86(7): 075443

45 Brey L, Fertig H A. Electronic states of graphene nanoribbons studied with the Dirac equation. Physical Review B: Condensed Matter and Materials Physics, 2006, 73(23): 235411

46 Han M Y, Ozyilmaz B, Zhang Y, Kim P. Energy band-gap engineering of graphene nanoribbons. Physical Review Letters, 2007, 98(20): 206805

PMID

47 Ezawa M. Peculiar width dependence of the electronic properties of carbon nanoribbons. Physical Review B: Condensed Matter and Materials Physics, 2006, 73(4): 045432

48 Park C H, Yang L, Son Y W, Cohen M L, Louie S G. Anisotropic behaviours of massless Dirac fermions in graphene under periodic potentials. Nature Physics, 2008, 4(3): 213-217

49 Park C H, Yang L, Son Y W, Cohen M L, Louie S G. New generation of massless Dirac fermions in graphene under external periodic potentials. Physical Review Letters, 2008, 101(12): 126804

PMID

50 Park C H, Son Y W, Yang L, Cohen M L, Louie S G. Electron beam supercollimation in graphene superlattices. Nano Letters, 2008, 8(9): 2920-2924

PMID

51 Morozov S V, Novoselov K S, Katsnelson M I, Schedin F, Ponomarenko L A, Jiang D, Geim A K. Strong suppression of weak localization in graphene. Physical Review Letters, 2006, 97(1): 016801

PMID

52 Suzuura H, Ando T. Crossover from symplectic to orthogonal class in a two-dimensional honeycomb lattice. Physical Review Letters, 2002, 89(26): 266603

PMID

53 Khveschenko D V. Electron localization properties in graphene. Physical Review Letters, 2006, 97: 036802

PMID

54 Dragoman D, Dragoman M. Giant thermoelectric effect in graphene. Applied Physics Letters, 2007, 91(20): 203116

55 Wei P, Bao W, Pu Y, Lau C N, Shi J. Anomalous thermoelectric transport of Dirac particles in graphene. Physical Review Letters, 2009, 102(16): 166808

PMID

56 Balandin A A, Ghosh S, Bao W, Calizo I, Teweldebrhan D, Miao F, Lau C N. Superior thermal conductivity of single-layer graphene. Nano Letters, 2008, 8(3): 902-907

PMID

57 Kane C L, Mele E J. Quantum spin Hall effect in graphene. Physical Review Letters, 2005, 95(22): 226801

PMID

58 Nandkishore R, Levitov L S, Chubukov A V. Chiral superconductivity from repulsive interactions in doped graphene. Nature Physics, 2012, 8(2): 158-163

59 Sarma S D, Adam S, Hwang E H. Electronic transport in two-dimensional graphene. Reviews of Modern Physics, 2011, 83(2): 407-439

60 Bonaccorso F, Sun Z, Hasan T, Ferrari A C. Graphene photonics and optoelectronics. Nature Photonics, 2010, 4(9): 611-622

61 Castro Neto A H, Guinea F, Peres N M R, Novoselov K S, Geim A K. The electronic properties of graphene. Reviews of Modern Physics, 2009, 81: 109-162

62 Beenakker C W J. Colloquium: Andreev reflection and Klein tunneling in graphene. Reviews of Modern Physics, 2008, 80(4): 1337-1354

63 Hasan M Z, Kane C L. Colloquium: topological insulators. Reviews of Modern Physics, 2010, 82(4): 3045-3067

64 Fleurence A, Friedlein R, Ozaki T, Kawai H, Wang Y, Yamada-Takamura Y. Experimental evidence for epitaxial silicene on diboride thin films. Physical Review Letters, 2012, 108(24): 245501

PMID

65 Vogt P, De Padova P, Quaresima C, Avila J, Frantzeskakis E, Asensio M C, Resta A, Ealet B, Le Lay G. Silicene: compelling experimental evidence for graphenelike two-dimensional silicon. Physical Review Letters, 2012, 108(15): 155501

PMID

66 Bianco E, Butler S, Jiang S, Restrepo O D, Windl W, Goldberger J E. Stability and exfoliation of germanane: a germanium graphane analogue. ACS Nano, 2013, 7(5): 4414-4421

67 Xu Y, Yan B, Zhang H J, Wang J, Xu G, Tang P, Duan W, Zhang S C. Large-gap quantum spin Hall insulators in tin films. Physical Review Letters, 2013, 111(13): 136804

PMID

68 Shareef S, Ang Y S, Zhang C. Room-temperature strong terahertz photon mixing in graphene. Journal of the Optical Society of America. B, Optical Physics, 2012, 29(3): 274-279

69 Ang Y S, Sultan S, Zhang C. Nonlinear optical spectrum of bilayer graphene in the terahertz regime. Applied Physics Letters, 2010, 97(24): 243110

70 Ang Y S, Zhang C. Subgap optical conductivity in semihydrogenated graphene. Applied Physics Letters, 2011, 98(4): 042107

71 Ang Y S, Zhang C. Enhanced optical conductance in graphene superlattice due to anisotropic band dispersion. Journal of Physics. D, Applied Physics, 2012, 45(39): 395303

72 Siegel P H. Terahertz technology. IEEE Transactions on Microwave Theory and Techniques, 2002, 50(3): 910-928

73 Hendry E, Hale P J, Moger J, Savchenko A K, Mikhailov S A. Coherent nonlinear optical response of graphene. Physical Review Letters, 2010, 105(9): 097401

PMID

74 Mikhailov S A. Non-linear electromagnetic response of graphene. Europhysics Letters, 2007, 79(2): 27002

75 Mikhailov S A, Ziegler K. Nonlinear electromagnetic response of graphene: frequency multiplication and the self-consistent-field effects. Journal of Physics Condensed Matter, 2008, 20(38): 384204

76 Dragoman M, Neculoiu D, Deligeorgis G, Konstantinidis G, Dragoman D, Cismaru A, Muller A A, Plana R. Millimeter-wave generation via frequency multiplication in graphene. Applied Physics Letters, 2010, 97(9): 093101

77 Wright A R, Xu X G, Cao J C, Zhang C. Strong nonlinear optical response of graphene in the terahertz regime. Applied Physics Letters, 2009, 95(7): 072101

78 Lim G K, Chen Z L, Clark J, Goh R G S, Ng W H, Tan H W, Friend R H, Ho P K H, Chua L L. Giant broadband nonlinear optical absorption response in dispersed graphene single sheets. Nature Photonics, 2011, 5(9): 554-560

79 Wang J, Hernandez Y, Lotya M, Coleman J N, Blau W J. Broadband nonlinear optical response of graphene dispersions. Advanced Materials, 2009, 21(23): 2430-2435

80 Hong S Y, Dadap J I, Petrone N, Yeh P C, Hone J, Osgood R M Jr. Optical third-harmonic generation in graphene. Physical Review X, 2013, 3(2): 021014

81 Wu S, Mao L, Jones A M, Yao W, Zhang C, Xu X. Quantum-enhanced tunable second-order optical nonlinearity in bilayer graphene. Nano Letters, 2012, 12(4): 2032-2036

PMID

82 Ishikawa K L.Nonlinear optical response of graphene in time domain. Physical Review B, 2010, 82(20): 201402

83 Feynman R P. Forces in molecules. Physical Review, 1939, 56(4): 340-343

84 Zhang C. Frequency-dependent electrical transport under intense terahertz radiation. Physical Review B: Condensed Matter and Materials Physics, 2002, 66(8): 081105

85 Ludwig A W W, Fisher M P A, Shankar R, Grinstein G. Integer quantum Hall transition: an alternative approach and exact results. Physical Review B: Condensed Matter and Materials Physics, 1994, 50(11): 7526-7552

86 Chen C F, Park C H, Boudouris B W, Horng J, Geng B, Girit C, Zettl A, Crommie M F, Segalman R A, Louie S G, Wang F. Controlling inelastic light scattering quantum pathways in graphene. Nature, 2011, 471(7340): 617-620

PMID

87 Gao F, Wang G, Zhang C. Strong photon-mixing of terahertz waves in semiconductor quantum wells induced by Rashba spin-orbit coupling. Nanotechnology, 2008, 19(46): 465401

PMID

88 Wolff P A, Pearson G A. Theory of optical mixing by mobile carriers in semiconductors. Physical Review Letters, 1966, 17(19): 1015-1017

89 Dong H M, Xu W, Tan R B. Temperature relaxation and energy loss of hot carriers in graphene. Solid State Communications, 2010, 150(37-38): 1770-1773

90 Sun D, Wu Z K, Divin C, Li X, Berger C, de Heer W, First P, Norris T. Ultrafast relaxation of excited Dirac fermions in epitaxial graphene using optical differential transmission spectroscopy. Physical Review Letters, 2008, 101(15): 157402

PMID

91 Butscher S, Milde F, Hirtschulz M, Malic E, Knorr A. Hot electron relaxation and phonon dynamics in graphene. Applied Physics Letters, 2007, 91(20): 203103

92 Bao W S, Liu S Y, Lei X L. Hot-electron transport in graphene driven by intense terahertz fields. Physics Letters. [Part A], 2010, 374(10): 1266-1269

93 McCann E, Fal’ko V I. Landau-level degeneracy and quantum Hall effect in a graphite bilayer. Physical Review Letters, 2006, 96(8): 086805

PMID

94 Koshino M, Ando T. Transport in bilayer graphene: calculations within a self-consistent Born approximation. Physical Review B: Condensed Matter and Materials Physics, 2006, 73(24): 245403

95 McCann C, Abergel D S L, Fal’ko V I. Electrons in bilayer graphene. Solid State Communications, 2007, 143(-2): 110-115

96 Chen L, Ma Z, Zhang C. Vertical absorption edge and temperature dependent resistivity in semihydrogenated graphene. Applied Physics Letters, 2010, 96(2): 023107

97 Edwards W F. Special relativity in anisotropic space. American Journal of Physics, 1963, 31(7): 482-489

98 Moon C Y, Han J, Lee H, Choi H J. Low-velocity anisotropic Dirac fermions on the side surface of topological insulators. Physical Review B: Condensed Matter and Materials Physics, 2011, 84(19): 195425

99 Park J, Lee G, Wolff-Fabris F, Koh Y Y, Eom M J, Kim Y K, Farhan M A, Jo Y J, Kim C, Shim J H, Kim J S. Anisotropic Dirac fermions in a Bi square net of SrMnBi₂. Physical Review Letters, 2011, 107(12): 126402

PMID

100 Wang J, Hernandez Y, Lotya M, Coleman J N, Blau W J. Broadband nonlinear optical response of graphene dispersions. Advanced Materials, 2009, 21(23): 2430-2435

101 Lim G K, Chen Z L, Clark J, Goh R G S, Ng W H, Tan H W, Friend R H, Ho P K H, Chua L L. Giant broadband nonlinear optical absorption response in dispersed graphene single sheets. Nature photonics, 2011, 5(9): 554-560

102 Hwang E H, Das Sarma S. Acoustic phonon scattering limited carrier mobility in two-dimensional extrinsic graphene. Physical Review B: Condensed Matter and Materials Physics, 2008, 77(11): 115449

103 Song J C, Reizer M Y, Levitov L S. Disorder-assisted electron-phonon scattering and cooling pathways in graphene. Physical Review Letters, 2012, 109(10): 106602

PMID

104 Betz A C, Jhang S H, Pallecchi E, Ferreira R, Feve G, Berroir J M, Placais B. Supercollision cooling in undoped graphene. Nature Physics, 2012, 9(2): 109-112

105 Graham M W, Shi S F, Ralph D C, Park J, McEuen P L. Photocurrent measurements of supercollision cooling in graphene. Nature Physics, 2012, 9(2): 103-108

106 Xu X G, Cao J C. Nonlinear response induced strong absorptance of graphene in the terahertz regime. Modern Physics Letters B, 2010, 24(21): 2243-2249

107 Weiss D, Zhang C, Gerhardts R R, Klitzing K, Weimann G. Density of states in a two-dimensional electron gas in the presence of a one-dimensional superlattice potential. Physical Review B: Condensed Matter and Materials Physics, 1989, 39(17): 13020-13023