Giant Faraday rotation induced by Berry phase in bilayer graphene under strong terahertz fields

High-order terahertz (THz) sideband generation (HSG) in semiconductors is a phenomenon with physics similar to high-order harmonic generation but in a much lower frequency regime. It was found that the electron-hole pairs excited by a weak optical laser can accumulate Berry phases along a cyclic path under the driving of a strong THz field. The Berry phases appear as the Faraday rotation angles of the emission signal under short-pulse excitation in monolayer MoS$_2$. In this paper, the theory of Berry phase in THz extreme nonlinear optics is applied to biased bilayer graphene with Bernal stacking, which has similar Bloch band features and optical properties to the monolayer MoS$_2$, such as time-reversal related valleys and valley contrasting optical selection rules. The bilayer graphene has much larger Berry curvature than monolayer MoS$_2$, which leads to a giant Faraday rotation of the optical emission ($\sim$ 1 rad for a THz field with frequency 1 THz and strength 8 kV/cm). This provides opportunities to use bilayer graphene and low-power THz lasers for ultrafast electro-optical devices.


I. INTRODUCTION
In semiconductors irradiated by a strong terahertz (THz) laser of frequency ω, the electron-hole pairs excited by a weak optical laser of frequency Ω are driven into large amplitude oscillations. They subsequently acquire kinetic energies much greater than the THz photon energy during the oscillations, which leads to optical emission at frequencies Ω±2Nω, where N is an integer. 1,2 This high-order THz-sideband generation (HSG) in semiconductors is a generalization of the high-order harmonic generation (HHG) [3][4][5][6] to the THz frequency regime. The sidebands in the frequency domain indicate THz modulation in the time domain, which has potential applications in Tbit/s optical communications, trillion-FPS optical imaging and wideband optical multiplexers.
The core physics of HHG and HSG is captured by the quantum trajectory theory. 1,5 When the electrons (or electronhole pairs) are driven by the strong laser field, their oscillation amplitudes are much larger than the wavepacket diffusion range. Therefore the quantum evolution of the electrons (or electron-hole pairs) is well described by quantum trajectories that satisfy the stationary phase condition (i.e., the least action condition in classical mechanics) in the formalism of Dirac-Feynmann path integrals. 1,5 A fundamental difference between HHG and HSG is that the electron-hole pair in semiconductors can have nontrivial Bloch states. 7 As a result, when the eletron (or hole) is driven by the THz field F (t) in the conduction (or valence) bands, not only does the quasi-momentum change according tok = −eF (t), 8 but also does its Bloch wavefunction evolve with k. Since the THz frequency is much smaller than the band gaps of semiconductors, the evolution is adiabatic (i.e., no interband tunneling is induced by the THz field). This adiabatic evolution along the quantum trajectory in semiconductors can accumulate a geometrical phase, in particular a Berry phase for a cyclic evolution. 9 A rich structures of the Bloch states in condensed matter systems (such as in topological in-sulators 10,11 ) thus leads to a variety of phase effects in extreme nonlinear optics.
An additional important difference between the HHG and HSG physics is that the electron-hole pairs are elementary optical excitations in solids. Therefore, the quantum trajectories in semiconductors, unlike those in atoms, can be triggered by lasers on demand at designed frequencies 12 or times (relative to the THz field oscillation). This excitation at will provides a great deal of controllability and flexibility for studying the quantum trajectories in extreme nonlinear optics.
In Ref. 7, the Berry phase dependent quantum trajectory theory was applied to monolayer MoS 2 with a band gap in the visible wavelength regime. The optical emission delayed by integer multiples of the THz period after the pulse excitation acquires a Faraday rotation. The rotation angle was shown to be exactly the Berry phase of the quantum trajectory (∼ 0.01 rad for a THz field with frequency 1 THz and strength 8 kV/cm). That result provides new opportunities to utilize THz extreme nonlinear optics of thin film materials for ultrafast electro-optical devices.
In this paper, we consider the transient optical response of Bernal stacked bilayer graphene excited by a weak optical pulse and a strong THz field. The bilayer graphene with an interlayer gate bias [13][14][15][16] has similar Bloch band features and optical properties as the monolayer MoS 2 but with a smaller band gap tunable up to ∼ 250 meV. Its conduction and valence-band edges are located at the corners of the 2D hexagonal Brillouin zone, and similar to MoS 2 , the optical interband transitions at the two time-reversal (TR) related valleys have nearly perfect but opposite polarization selection rules. 17,18 These features make bilayer graphene an excellent candidate for the observation of the Berry phases of quantum trajectories. In bilayer graphene without interlayer bias, the conduction (valence) state acquires a ±2π Berry phase (pseudospin winding number being ±2) along any closed path around the Dirac points, 14,19 and the Berry curvature is completely concentrated on the singular monopole at the Dirac point. In bi-ased bilayer graphene that has an energy gap opened at the Dirac points, the Berry curvature becomes non-singular while its distribution is still concentrated in a small region in the vicinity of the Dirac points (since the Berry curvature of the region far from the Dirac point is zero as in the unbiased case) (see Fig. 1). Thus the Berry curvature in bilayer graphene is much (about 100 times) larger than in monolayer MoS 2 . This leads to a much larger Berry phase of a quantum trajectory and in turn giant Faraday rotation. In other words, bilayer graphene requires a much weaker THz laser to induce an observable Faraday rotation, and serves as a better material for ultrafast electro-optical devices than MoS 2 . Besides, the band gap at the Dirac points is tunable by an applied gate bias between zero and midinfrared energies, 15,16 which offers outstanding controllability.
This paper is organized as follows. In Sec. II, we present a general theory of optical response in semiconductors under a strong THz field, with the Berry phase effects included. The optical response under a pulsed optical excitation is studied using the quantum trajectory theory. In materials with both TR symmetry and optical selection rules, the Faraday rotation angle of the optical emission delayed by multiples of the THz period is shown to be exactly equal to the Berry phase of the stationary quantum trajectory. In Sec. III, we show the valley contrasted optical selection rules of the biased bilayer graphene, which makes it a good system for the observation of the Berry phases of quantum trajectories. Furthermore, calculation shows that the bilayer graphene has a Berry curvature much lager than that in monolayer MoS 2 , which leads to a much larger Berry phase of the stationary trajectory and in turn much larger Faraday rotation. The giant Faraday rotation is verified by numerical simulations. Sec. IV concludes this paper.

II. MODEL AND FORMALISM
We consider a semiconductor under an elliptically polarized THz field with ω much smaller than the energy gap of the material (so that the THz field does not induce interband tunneling). The Hamiltonian in the Bloch state representation H (k) evolves adiabatically in the k-space H (k) → H k (t) under the driving of this field, with where A (t) is the electromagnetism vector potential with F = −∂A/∂t, and k 0 = eF/ω. Then we study the interaction of this system with a weak optical laser that creates electron-hole pairs at the band edge. The interaction Hamiltonian isĤ I = −P · E I e −iΩt + h.c.. The interband polarization operatorP in the Bloch state representation iŝ whereê andĥ are electron and hole operators, and the interband dipole moment d µν,k is given by 20 Here + and − denote the conduction and valence bands, respectively, µ and ν are the spin or pseudo-spin indices introduced to label degenerate bands, 21 and |±, µ, k are the Bloch states of the conduction (valence) bands with eigenenergies E ± k . We assume that the initial state is the vacuum state |G with empty conduction bands and filled valence bands. After excitation by an optical laser, the electron (hole) is driven into adiabatic evolution in the conduction (valence) band by the THz field, and thus obtains a geometric phase being the instantaneous eigenstates of H k (t) with instantaneous eigenenergies E ± k . In general, the Berry connection can be non-Abelian (i.e. A ± k µν = i ±, µ,k ∇ k ±, ν,k 0 for µ ν) and the geometric phase factor becomes a unitary ma- , whereT is the time-ordering operator. The electron and hole recombine at a later time, leading to optical emission modified by the geometric phase.
The linear optical response of THz field driven semiconductors has been formulated in Ref. 7. If the THz field does not mix different spin or pseudo-spin states of the degenerate energy bands, the Berry connections are Abelian. The optical response is simplified to k is the energy of the electron-hole pair, and Ak = A + k − A − k is the combined Berry connection of the electron-hole pair. Now we consider a system that has TR symmetry and nontrivial Berry curvatures in its energy bands. If one state of the Kramers pair is denoted as the pseudospin state ⇑ and the other as ⇓, we have the relations according to the TR symmetry. We apply a short optical laser pulse to the system at time t = 0, with the duration much shorter than the THz period T = 2π/ω. The pulse can be approximated by a δ-pulse E I (t ) ≈ Eδ (t ). After an integer multiple of the THz period t n = nT , the electron (hole) undergoes a cyclic evolution in the conduction (valence) band and the geometric phase becomes the gauge invariant Berry phase. 9 Using equation (6), we get i.e. the two TR related states have opposite Berry phases.
Then the response at t n is simplified as We see that P (t n ) is given by the interference between two responses with the same dynamical phase but opposite Berry phases, which, as shown below, results in a Faraday rotation of the optical emission.
The main consequence of equation (8) can be understood using the stationary phase formalism (or the quantum trajectory theory). 1,5,7 In the path integral, the electron-hole pair moves along all possible trajectories when driven by the THz field, with the phase given by the action S (k) = t n 0 εk (τ) dτ. The Berry phase is dropped from this action that determines the electron-hole motion, since it is generally much smaller than the dynamical phase (see Fig. 2). As the THz field is strong, the motion amplitude of the electron-hole pair is much larger than the wavepacket diffusion range. Thus the response is dominated by the stationary phase points of the action (i.e. trajectories that satisfy the least action condition) where vk = ∇ k εk is the semiclassical velocity of the electronhole pair. The stationary phase condition in Eq. (9) means the electron under the acceleration by the THz field returns to the hole after nT for recombination. Then the Berry phase is determined mainly by the stationary trajectory: where d −k(0) = dk (0) near the band edge has been used and k s is the solution to equation (9). The Berry phase φ (n) B (k s ) is accumulated along the stationary trajectory. When multiple solutions exist, the response is given by the interference of all possible stationary trajectories. The optical selection rule near the Dirac points in bilayer graphene is such that the interband dipole moment d k ≈ d cv,k e x − ie y / √ 2, i.e. the optical transitions at valley ⇑ (⇓) is coupled to the σ+ (σ−) polarized light (see Sec. III). The linear susceptibilities for the σ ± -polarized lights are respectively If the excitation light is linearly polarized in the x-y plane, the Faraday rotation of the optical emission at t n is exactly given by the Berry phase φ FR (t n ) = φ (n) B (k s ). The Faraday rotation caused by the elliptically polarized THz field can be intuitively understood as following. 7 The linearly polarized optical laser is a superposition of two opposite circular polarizations, which cause respective transitions From the preceding discussion it is clear that the theory can be applied to any material with TR symmetry, optical selection rules and nontrivial Berry phases, such as the topological insulators, 10,11 monolayer MoS 2 and other group-VI dichalcogenides 17,22-25 and bilayer graphene. [13][14][15][16] In Ref. 7, the monolayer MoS 2 was studied and the result -Faraday rotation equals the Berry phase of the stationary trajectory -was verified by numerical calculation. In Sec. III, we will investigate the bilayer graphene. It has about 100 times larger Berry curvature at the band edges than the monolayer MoS 2 does, which leads to a giant Faraday rotation of the optical emission.

III. BERRY PHASE INDUCED FARADAY ROTATION IN BILAYER GRAPHENE
In bilayer graphene with Bernal (A-B) stacking, the band structure is well described by the tight-binding Hamiltonian with an intralayer nearest-neighbor hopping t ≈ 3eV, an inter-layer nearest-neighbor hopping γ ≈ 0.4eV, and an interlayer bias 2∆, [13][14][15] Here , 0 with a ≈ 2.46Å being the lattice constant. Expanding the momentum near the two Dirac 2 at (≈ 10 6 m/s). The states at K ± valleys are related by TR operation and the corresponding Hamiltonians are given by where k ± = k x ±ik y and the K ± valley is denoted by the pseduospin notation ⇑ / ⇓. In this paper, it is assumed that the valley coherence time in bilayer graphene is longer than the THz period and thus the intervalley scattering can be neglected. The calculated band structure (as shown in Fig. 1 (a)) reproduces the well-konwn Mexican hat structure.
The bilayer graphene has two conduction bands (positive energy bands) and two valence bands (negative energy bands). We assume that the Fermi level is kept in the energy gap while tuning the gap value. 15 Hence the system is initially in the vacuum state |G with empty conduction bands and filled valence bands. Since the optical laser is near resonant with the energy gap between the bottom conduction band and the top valence band, only the optical transitions between these two bands are considered, with |±, µ, k (µ =⇑, ⇓) now denoting the corresponding Bloch states. Furthermore, it was shown in Ref. 18 that there is a nearly perfect optical selection rule for the interband transitions between the two bands near the Dirac points, where the valley K ± favors the σ ± polarized transition, respectively. This can be seen from Fig. 1 (b), which plots the degree of circular polarization with d ± µµ,k = d µµ,k · e x ± ie y / √ 2. Therefore, bilayer graphene serves as an ideal system for studying the Faraday rotation and Berry phases of quantum trajectories. Now we calculate the combined Berry curvature of the bottom conduction band and the top valence band, defined as where with j = ± and |n, µ, k denoting the four eigenstates of H µ (k) with eigenenergy E n k . Note that the Berry curvature is Abelian. The states at the two valleys have opposite Berry curvatures and thus give opposite Berry phases as pointed out in Sec. II. If there is no interlayer bias, the electron (hole) acquires a ±2π Berry phase (or pseudospin winding number ±2) after traversing a path around the Dirac points K ± . 14,19 When a bias 2∆ = 0.3 eV is applied, the distribution of Ω xy is shown Fig. 1  (c). The Berry curvature is concentrated in a small region enclosed by the band edge of the Mexican hatlike energy dispersion, being much larger than in monolayer MoS 2 (≈ 22.3 (π/a) −2 ). 7,17 Since the Berry phase equals the Berry curvature flux through the area enclosed by the path in the k-space, the electron-hole pair in bilayer graphene acquires a much larger Berry phase during the cyclic evolution, which leads to a giant Faraday rotation of the optical emission.
Then we study the stationary trajectories, i.e. the stationary phase points of the dynamical phase φ (n) D (k) = t n 0 εk (τ) dτ, un-der the driving of the THz field. Figure 2 plots the distribution of the dynamical phase φ (1) D (k) over one THz period and the corresponding Berry phase φ (1) B (k) at the K + valley for different filed strengths F. The dynamical phase is indeed much larger than the Berry phase. This justifies the approximation in Eq. (9) for determining the quantum trajectories. When the THz field is weak (i.e. k 0 , the amplitude of the quantum trajectory in k-space, is small as compared with the radius of the "Mexican hat" ring (≈ 0.0233π/a) around the Dirac point), we can expand ε˜k (τ) around k and get Hence the stationary phase points are exactly the extreme points of the energy band, i.e. the Dirac points and the ringshaped bandedges (Fig. 2 (a)). Note that the trajectory around the Dirac point has a Berry phase close to that of the trajectory around the bandedge point ( Fig. 2 (a)). Thus the Faraday  (6) and (7). The interlayer bias is 2∆ = 0.3 eV, the frequency of the THz field is ω = 4 meV and the THz field is circularly polarized (i.e. θ = π/4). In (a), the THz field strength F = 2, 4 and 8 kV/cm (i.e., k 0 = 0.0039, 0.0078 and 0.0157 π/a), and the dynamical phase has two kinds of stationary phase points (indicated by the arrows) at the Dirac point and the bandedge. In (b), the THz field strength F = 12, 16 and 20 kV/cm (i.e., k 0 = 0.0235, 0.0313 and 0.0392 π/a), and only the Dirac point is the stationary phase point (indicated by the arrow). rotation is well approximated by the Berry phase of the stationary trajectory around the Dirac point. When the THz field is strong enough (i.e. k 0 is comparable to the "Mexican hat" ring radius (≈ 0.0233π/a)), only the Dirac points are the stationary phase points ( Fig. 2 (b)). In this case, the Faraday rotation is also given by the Berry phase of the stationary trajectory centered at the Dirac point.
To verify the prediction of the quantum trajectory method (i.e. Faraday rotation equals the Berry phase of the stationary trajectory), we calculate the Faraday rotation directly through numerical integration of Eq. (5) for bilayer graphene and compare it with the Berry phase φ (n) B of the quantum trajectories. In the numerical calculation, the frequency of the THz field is ω = 4 meV and the optical pulse has the gaussian form Ee x exp −iΩt − t 2 /δt 2 , with Ω = 2∆ and ωδt = 0.2 ( 2π) (the optical laser spectrum width is ∼ 2/δt = 40 meV, much smaller than the valence/conduction band splitting). The re-sults calculated for various THz field ellipticity θ and THz field strength F are shown in Fig. 3, which confirms that the Faraday rotation is well approximated by the Berry phase of the trajectory around the Dirac point. The Faraday rotation angle is about 100 times larger than in monolayer MoS 2 , 7 consistent with the Berry curvature difference between the two materials. The Berry phase is given by the Berry curvature flux through the area element πk 2 0 sin θ cos θ enclosed by the closed path in k-space (equation (2)), which explains the approximate sinusoidal and parabolic dependence of the rotation angle on θ and F, respectively. Also as expected from the quantum trajectory analysis, the Berry phases around the bandedge points deviate appreciably from the numerical results when the THz field is strong enough (see Fig. 3 (b)), which means the bandedge points fail to be the stationary phase points.
In the linear response regime, the giant Faraday rotation proposed is not affected by the recombination of the electronhole pair. However, it requires the valley coherence time to be longer than the period of the THz field. Although the valley coherence time in bilayer graphene has not been determined to the best of our knowledge, the experiments for an analogous material (monolayer WSe 2 ) shows that the valley coherence time can be even longer than the recombination time of electron-hole pairs, 26 which indicates the observation of the predicted giant Faraday rotation is promising. On the other hand, the Farday rotation can also be employed to study the intervalley coherence.  3: (color online). Faraday rotation angle φ FR of the optical emission at one THz period after pulse excitation, as a function of (a) the polarization ellipticity θ and (b) the strength F of the THz field. The symbols show the numerical integration (NI) results and the lines are the Berry phases of the quantum trajectories centered at the stationary points (0 for the Dirac point, and k 2 for the bandedge point with k y = 0). In (a), the THz field strength F = 8 kV/cm (i.e., k 0 = 0.0157π/a). In (b), the THz field is circularly polarized (i.e., θ = π/4). The other parameters are the same as in Fig. 2.

IV. CONCLUSIONS
In summary, the Berry phase dependent theory of THz extreme nonlinear optics developed in Ref. 7 is applied to the biased bilayer graphene. The Faraday rotation angle of the emission light delayed by integer multiples of the THz period is the Berry phase of the quantum trajectory. As the bilayer graphene has large Berry curvatures near the Dirac points, the Berry phase of the stationary trajectory is large, which leads to a giant Faraday rotation. The predictions by the quantum trajectory analysis are confirmed by numerical simulations. This result provides an opportunity to use bilayer graphene and relatively low-power THz lasers for ultrafast electro-optical devices.