Role of electronic excitations in magneto-Raman spectra of graphene

Experimental observation of electronic excitations in graphene is mostly associated with angle-resolved photoemission spectroscopy [1, 2] and optical absorption in monolayers [3–9] and bilayers [5, 6, 10–16] where Raman spectroscopy has been used to study graphene phonons [17], and became the method of choice for determining the number of atomic layers in graphitic flakes [17–23]. In particular, single- and multiple phonon-emission lines in the Raman spectrum of graphene and the influence of the electron–phonon coupling on the phonon spectrum have been investigated in great detail in [24–32].

This paper presents a theory of inelastic light scattering in the visible range of photon energies accompanied by electronic excitations in mono- and bilayer graphene both with and without an external magnetic field. We classify the relevant modes according to their symmetry, analyse peculiar selection rules for the Raman-active excitations of electrons between Landau levels (LLs) in graphene at strong magnetic fields, review the fine structure of the phonon-related Raman under conditions of magneto-phonon resonance (when the inter-LL transition is in resonance with a phonon) and model the dependence of the spectra on the value and homogeneity of carrier density and strain in the sample.

Graphene is a gapless semiconductor [33–36], with an almost linear Dirac-type spectrum ε = αvp in the monolayer case and parabolic ε = αp²/2m for bilayer, in the conduction (α = +) and valence (α = −) bands, which touch each other in the corners of the hexagonal Brillouin zone (BZ), usually called valleys. In bilayer graphene, there are two additional parabolic bands with the gap 2γ₁. The conduction–valence band degeneracy in the BZ corners of the band structure of graphene is prescribed by the hexagonal symmetry C_6v of its honeycomb lattice, and it is natural to relate Raman-active modes to the irreducible representations of the point group C_6v.

We argue that the dominant electronic modes generated by inelastic scattering of photons with energy Ω less than the bandwidth of graphene are superpositions of the interband electron–hole (e–h) pairs which have symmetry of the representation A₂ of the group C_6v and are odd with respect to the inversion of time. Their excitation process may go in two ways both of which consist of the same steps, but in reversed order. In the first process, the absorption of a photon with energy Ω transfers an electron from an occupied state in the valence band into a virtual state in the conduction band, followed by emission of the second photon with energy $\tilde {\Omega }=\Omega -\omega$. In the second process the reversed order of events leads to the virtual electron state with large deficit of energy as opposed to the first process with large excess. Net amplitude is determined by the sum of partial amplitudes of the opposite processes. The dominance of such a process over the one involving the contact interaction [37–39] of an electron with two photons is a peculiarity of the Dirac-type electrons in graphene.

A strong magnetic field leads to the quantization of the electronic spectrum into LLs [40]. In graphene, LLs appear symmetrically in the conduction bands ^α_n (α = ±) with peculiar zero-energy LLs in monolayer (n = 0) and bilayer (n = 0.1). For the inelastic light scattering in graphene, we find peculiar selection rules, n⁻ → n⁺ of the dominant Raman-active transitions (solid line in figure 3), in contrast to the Δn = ± 1 transitions between LLs, which are dominant in the absorption of left- and right-handed circularly polarized infrared photons [7]. The weaker Raman mode with a similar inter-LL selection rule also differs from the absorption process, having different symmetry properties which allow the interaction with Raman scattering of light on optical phonons creating anticrossings at the G-peak line; see figure 4. Raman spectroscopy, therefore, provides data supplementary to that obtained in optical absorption.

The following theory is based on the effective low-energy Hamiltonian derived from the tight-binding model of electron states in graphene. It describes electrons in the conduction and valence bands around the BZ corners K (ξ = +) and K' (ξ = −) [41, 42]. For monolayer graphene the effective Hamiltonian reads

$$\begin{equation} H_{\mathrm{m}}=\xi v\boldsymbol{\sigma }\cdot \mathbf{p} - \frac{v^{2}}{6\gamma_{0}}\bigl(\sigma ^{x}(p_{x}^{2}-p_{y}^{2})-2\sigma^{y}p_{x}p_{y}\bigr). \end{equation}\noindent \tag{ 1 }$$

Here, σ = (σ^x,σ^y) are Pauli matrices acting on the sublattice components of electronic states on A and B sublattices of graphene honeycomb lattice, and p being the in-plane momentum counted from the BZ corner. The first term in equation (1) determines Dirac spectrum electrons. The monolayer hopping parameter γ₀ ≈ 3 eV determines the bandwidth, ∼6γ₀, and Dirac velocity of electrons v ≈ 10⁸ cm s⁻¹ [10, 12, 19, 43]. The second term takes into account weak trigonal warping [41] breaking the continuous rotational symmetry of the Dirac Hamiltonian. The basis is constructed of wavefunctions corresponding to atomic sites A,B in the valley K and B,A in K'. Since four-spinors {φ_AK,φ_BK,φ_BK',φ_AK'} realize four-dimensional (4D) irreducible representation of the full symmetry group of the crystal, matrix operators can be combined into irreducible representations [26, 27, 44] of the group C_6v; see table 1.

Table 1. Classification of the Hamiltonian components by C_6v irreducible representations.

C6v rep.	A1	A2	B1	B2	E1	E2
Matrix	1	σz	ξ	ξσz	ξσ	ez×σ
t→−t	+	−	−	+	−	+
	$(\mathbf {l}\cdot {\tilde{\bf l}}^{*})$	$(\mathbf {l}\times {\tilde{\bf l}}^{*})_{z}$			$\mathbf {l}, {\tilde{\bf l}}^{*}$	M,d
	p2				p	(2pxpy,p2x−p2y)
						$\mathbf {u}=\frac {1}{\sqrt 2}(\mathbf {u}_{A}-\mathbf {u}_{B})$

Bilayer graphene, figure 1(b), consists of two coupled sheets of monolayer with AB (Bernal) stacking (as in bulk graphite [45]), and has symmetry D_3h which differs from C_6v by the accompaniment of a z → − z reflection to the rotations by π/3, comparing to monolayer graphene, but, otherwise, D_3h and C_6v are isomorphic. Its unit cell contains four inequivalent atoms A1, B1, A2 and B2 where numbers denote the layer. The conventional Hamiltonian in the vicinity of the valleys is

$$\begin{equation} H_{b}= \left(\begin{array}{@{}cc@{}} \frac{v_{3}}{v} \sigma^{x} H_{m} \sigma^{x} & H_{m} \\ H_{m} & \gamma_{1}\sigma^{x} \end{array} \right), \end{equation}\noindent \tag{ 2 }$$

where v₃ ∼ 0.1v is related to the weak direct A1 → B2 interlayer hops [19, 36, 43], and γ₁ ∼ 0.4 eV is the direct interlayer coupling [1, 5, 6, 10, 12, 19, 31, 43, 46]. The basis is constructed using components corresponding to atomic sites A1,B2,A2,B1 in the valley K and B2,A1,B1,A2 in K'.

Zoom In Zoom Out Reset image size

The electronic dispersions around the K point for mono- and bilayer are compared in figure 1. Two low-energy bands touch each other at the centre of the valley, which is also the position of the Fermi energy in the neutral structure (a neutrality point), and are described by the effective two-band low-energy Hamiltonian written in the basis of orbitals on sites A1 and B2 [36],

$$\begin{equation} {\skew3\hat{H}}_{b\,{\rm eff}} = -\frac{v^{2}}{\gamma_{1}}\bigl(\sigma^{x}( p_{x}^{2} - p_{y}^{2}) + 2 \sigma^{y} p_{x}p_{y} \bigr). \end{equation} \tag{ 3 }$$

Two others, referred to as high-energy bands, are split by the interlayer coupling γ₁ from the neutrality point.

To include the interaction of the electrons with photons in the effective Hamiltonians (1) and (2), we construct the canonical momentum $\mathbf {p}-\frac {e}{c}\mathbf {A}$, where the vector potential of light

$$\begin{equation} \mathbf{A}= \sum_{\mathbf{l},\mathbf{q},q_{z}} \frac{\hbar c}{\sqrt{2\Omega}} \left(\mathbf{l}\,{\rm e}^{{\rm i}(\mathbf{q}\mathbf{r}-\Omega t)/\hbar }b_{\mathbf{q},q_{z},\mathbf{l}}+\mathrm{h.c.}\right) \end{equation} \tag{ 4 }$$

includes an annihilation operator b_q,qz,l for a photon with in-plane momentum q, energy Ω (which determines its out-of-plane momentum component $q_{z}=\sqrt {\Omega ^{2}/c^{2}-\mathbf {q}^{2}}$) and polarization l. Expanding the Hamiltonians up to the second order in the vector potential, one obtains the interaction part

$$\begin{equation} H_{{\rm int}} = -\frac{ev}{c}\mathbf{J}\cdot\mathbf{A} + \frac{e^{2}}{2c^{2}}\sum_{i,j}(\partial^{2}_{p_{i} p_{j}}H) \, A_{i}A_{j}, \end{equation} \tag{ 5 }$$

where $(ev/c) J_{i}=(e/c)\partial _{p_{i}}\hat {H}$ is the current vertex and $\frac {e^{2}}{2c^{2}} \partial ^{2}_{p_{i} p_{j}}H$ is the two-photon contact-interaction tensor.

The full amplitude R = R_D + R_w of inelastic Raman scattering of the photon on electrons changing its energy Ω to $\tilde {\Omega } = \Omega - \omega$ (ω = _f − _i is the Raman shift) and momentum and polarization from q,l to $\tilde {\mathbf {q}}, {\tilde{\bf l}}$ is described by the Feynman diagrams shown in figure 2.

Zoom In Zoom Out Reset image size

The building blocks of the diagrams include Green's functions for the electrons and the electron–photon interaction vertices:

Thus, the inelastic light scattering amplitude is contributed by a one-step process R_w (contact interaction [38]) and a two-step process R_D involving an intermediate virtual state. The two-step process, such as shown in figure 2, consists of an absorption (emission) of a photon with energy Ω ( $\tilde {\Omega }$) transferring an electron with momentum p from an occupied state in the valence band into a virtual intermediate state, followed by another electron emission (absorption) of the second photon with energy $\tilde {\Omega }$ (Ω), which moves the electron to the final state with momentum $\mathbf {p+q-}\tilde{\bf q}$,

$$\begin{equation} R_{{\rm D}}=\frac{e^{2}\hbar^{2}v^{2}}{2\sqrt{\Omega\tilde\Omega}}\Bigl( (\mathbf{J}_{\mathbf{q}}\cdot \mathbf{l}) G^{A}_{\Omega+\epsilon_{\rm i}} (\mathbf{J}_{-\tilde{\mathbf{q}}}\cdot {\tilde{\bf l}}^{*}) + (\mathbf{J}_{-\tilde{\mathbf{q}}}\cdot {\tilde{\bf l}}^{*}) G^{A}_{-\Omega+\epsilon_{\rm f}} (\mathbf{J}_{\mathbf{q}}\cdot \mathbf{l}) \Bigr). \end{equation} \tag{ 7 }$$

In graphene, the two-step process via the virtual state and the contact interaction process have different properties and polarization selection rules, which reflects their intricate relation with the irreducible representations of the symmetry group of the crystal. Generally, the components of the scattering amplitude R^ij realize a representation of the symmetry group C_6v. Since vectors l, $\tilde {\mathbf {l}}^{*}$ belong to E₁ representation, the R representation can be expanded into irreducible ones

$$\begin{equation} E_{1}\otimes E_{1}=A_{1}\oplus A_{2}\oplus E_{2}. \end{equation} \tag{ 8 }$$

Forming the corresponding combinations from the components of the polarization vector, see table 1,

$$\begin{eqnarray} &&\Xi_{s}=\bigl|\mathbf{l}\times {\tilde{\bf l}}^{*}\bigr|^{2}, \quad \Xi_{s2}=\bigl|\mathbf{l}\cdot {\tilde{\bf l}}^{*}\bigr|^{2}, \quad \Xi_{\mathrm{o}}=|\mathbf{d}|^{2}=1+(\mathbf{l}\times \mathbf{l}^{*})({\tilde{\bf l}}\times {\tilde{\bf l}}^{*}), \end{eqnarray} \noindent \tag{ 9 }$$

$$\begin{eqnarray} &&\mathrm{where}\quad \mathbf{d}=(l_{x}\tilde{l}_{y}^{*}+l_{y}\tilde{l}_{x}^{*},l_{x}\tilde{l}_{x}^{*}-l_{y}\tilde{l}_{y}^{*}) \quad \mathrm{and}\quad \Xi_{s}+\Xi_{s2}+\Xi_{\mathrm{o}}=2, \end{eqnarray} \noindent \tag{ 10 }$$

we can write the scattering probability in the most general form

$$\begin{equation} w=w_{s}\Xi_{s}+w_{s2}\Xi_{s2}+w_{o}\Xi_{\mathrm{o}}. \end{equation} \tag{ 11 }$$

The first two terms with polarization factors Ξ_s and Ξ_s2, corresponding to A₂ and A₁ representations, respectively, describe the contribution of photons scattered with the same circular polarization as the incoming beam. The third term, with polarization factor Ξ_o, corresponding to E₂ representation, represents the scattered photons with circular polarization opposite to the incoming beam.

Microscopically, the probability for a photon to undergo inelastic scattering from the state (q,q_z) with energy Ω into a state ${(\tilde {\bf q}},\tilde {q}_{z})$ with energy $\tilde {\Omega }$, by exciting an e–h pair in graphene, is

$$\begin{equation} w=\frac{2}{\pi\hbar^{3}}\sum_{\alpha_{\mathrm{i}}\,\alpha_{\mathrm{f}}}\int \mathrm{d}^{2}\mathbf{p} |\langle f|R|i\rangle|^{2} \times f_{\rm i}\left( 1-f_{\rm f} \right)\delta\left( \epsilon_{\mathrm{i}} + \omega -\epsilon_{\mathrm{f}} \right), \end{equation} \tag{ 12 }$$

where f_i and f_f are filling factors of the initial and final electronic states, respectively, the spin and valley degeneracies have already been taken into account, and the initial and final states are defined by the momentum p and (sub)band index α_i/f: $|i(f)\rangle =|\mathbf {p}(\mathbf {p+q}-\tilde {\bf q},\alpha _{\mathrm {i/f}})\rangle$. The angle-integrated spectral density of Raman scattering g(ω) is

$$\begin{equation} g(\omega) = \frac{1}{c}\int\!\!\!\int \frac{\mathrm{d}^{2}\tilde{\mathbf{q}}\mathrm{\,d}\tilde{q}_{z}}{(2\pi\hbar)^3}\, w\,\, \delta \left( \tilde{\Omega} - c\sqrt{\tilde{\mathbf{q}}^{2}+{\tilde{q}_{z}}^{2}} \right), \end{equation} \tag{ 13 }$$

and quantum efficiency $I=\int\!g \,\mathrm {d}\omega$ expresses the total probability for a single incoming photon to scatter inelastically on an electron and excite an e–h pair in the low-energy part of the spectrum.

2.1. Raman scattering in monolayer graphene

Calculation of the Raman spectrum and intensity of the signal in each of the distinct polarization components, w_s, w_s2 and w_o, requires an explicit expression for the Raman scattering amplitude R. Raman measurements are usually carried out at $\Omega ,\tilde {\Omega }\sim 1$–3 eV, which is much smaller than the bandwidth [17]. Raman shift, which measures the electronic excitation energy, is ω ≲ 0.2 eV , allowing us to expand the Green's functions over the large Ω and to perform summation over the intermediate virtual states of the process.

Both the electron and the hole of the e–h excitation have almost the same momentum ($\mathbf {p}+\mathbf {q}-\tilde {\mathbf {q}}$ and p, respectively), since $\mathbf {q},\tilde {\mathbf {q}}\ll \mathbf {p}$ and the momentum transfer from light is negligible (v/c ∼ 3 × 10⁻³), so that below we use J_q = ξσδ_p,p+q. For vp ≪ Ω the dominant contributions to the Raman scattering amplitude are

$$\begin{eqnarray} R_{\mathrm{{\rm D}}} &\approx \frac{(e\hbar v)^{2}}{\Omega^{2}}\left(-\mathrm{i} \sigma^{z} (\mathbf{l}\times {\tilde{\bf l}}^{*})_{z}+ \frac{\mathbf{M}\cdot\mathbf{d}}{\Omega}\right), \end{eqnarray} \tag{ 14 }$$

$$\begin{eqnarray} \mathbf{M}&=\xi v (\sigma_{x}p_{y}+\sigma_{y}p_{x},\sigma_{x}p_{x}-\sigma_{y}p_{y}),\nonumber\\ R_{w} &=\frac{e^{2}v^{2}\hbar ^{2}}{6\Omega \gamma_{0}}(\mathbf{e}_{z}\times \boldsymbol{\sigma })\cdot \mathbf{d}, \end{eqnarray} \tag{ 15 }$$

where the main order of R_D has the matrix form of representation A₂ of C_6v. Contact interaction R_w is responsible for the creation of the excitations with the symmetry of the representation E₂; see table 1. For the scattering of photons with Ω < γ₀, R_w ≪  R_D [47].

In undoped graphene the inter-band e–h pairs with _f ≈ −_i ≈ ω/2 are the only allowed electronic excitations. The probability of the creating an electronic excitation by photons with Ω < γ₀, $\mathbf {q}-\tilde {\mathbf {q}}\ll \mathbf {p}$ ($\omega \gg \frac {v}{c}\Omega$) in the main order of the expansion in ω/Ω ≪ 1 is

$$\begin{equation} w\approx \hbar e^{4}v^{2}\frac{\omega}{\Omega^{4}}\left\{ \Xi_{s}+\frac{\Omega^{2}}{2(6\gamma_{0})^{2}}\Xi_{\mathrm{o}} \right\}. \end{equation} \tag{ 16 }$$

It is dominated by the contribution R_D of the first two diagrams in figure 2, describing two-step processes involving a virtual intermediate state.

Doping in graphene, with chemical potential μ ≫ Ωv/c, blocks inter-band electronic excitations with energy transfer ω < 2μ. After integrating over all directions of the propagation of scattered photons, we find the spectral density of the angle-integrated Raman signal,

$$\begin{equation} g(\omega )\approx \frac{1}{4}\Xi_{s}\left( \frac{e^{2}}{\pi \hbar c}\,\frac{v}{c}\right)^{2} \frac{\omega}{\Omega^{2}} \theta(\omega-2\mu). \end{equation} \tag{ 17 }$$

The total contribution of the electronic scattering to the total Raman efficiency is

$$\begin{equation} I_{0}=\int_{0}^{\Omega/2}g(\omega )\mathrm{d}\omega \sim \left( \frac{e^{2}}{2hc}\,\frac{v}{c}\right)^{2}\sim 10^{-10}. \end{equation} \tag{ 18 }$$

In doped graphene, one may also expect to see some manifestation of the intra-band e–h excitations in the vicinity of the Fermi level and electron energies _i ≈ _f ≈ μ. Their analysis requires taking into account all terms in the expression for R_D, equations (15). Then, we find that, for ω ⩽ (v/c)Ω ≪ Ω < γ₀, the yield of this low-energy feature is $\delta I=\int \delta g(\omega )\mathrm {d}\omega \sim 10^{-15}$ for Ω ∼ 1 eV .

2.2. Raman scattering in bilayer graphene

For bilayer graphene with the interlayer coupling less than the base photon energy, $\gamma _{1}\ll \Omega ,\tilde {\Omega }$, we also study the low-energy excitations in the final states, in particular, those with ω ≪ γ₁. The reason for restricting the following analysis to the low-energy excitations is that the higher-energy excitations decay fast due to electron–electron and electron–phonon processes, leaving very little chance to observe any structure in the Raman signal. As a result, we employ the same expansion ω/Ω ≪ 1 as in section 2.1. Keeping terms up to the Ω⁻³ order (the latter appears when the virtual state is taken to be in the high-energy bands) in the expansion of the Green's functions and performing summation over the intermediate states of the process, we obtain the amplitude R in the form of a 4 × 4 matrix

$$\begin{eqnarray*} R &\approx \frac{(e\hbar v)^{2}}{\Omega^{2}} \left( -\mathrm{i} \left( \begin{array}{@{}cc@{}} \sigma_{z} & 0 \\ 0 & \sigma_{z} \end{array} \right) ( \mathbf{l}\times\tilde{\mathbf{l}}^{*} )_{z} + \frac{\mathcal{M}\cdot\mathbf{d}}{\Omega} \right),\\ \mathcal{M} &= \left( \begin{array}{@{}cc@{}} \gamma_{1}\sigma_{x}(\mathbf{e}_{z}\times\boldsymbol{\sigma})\sigma_{x} & \mathbf{M} \\ \mathbf{M} & 0 \end{array} \right) . \end{eqnarray*}$$

Since the electronic modes with the excitation energies ω < γ₁/2 result from the e–h pairs in the low-energy bands described by the low-energy Hamiltonian H_eff in equation (3), we take only the part of R that acts in that 2D Hilbert space, leaving the terms in the lowest relevant order in vp/γ₁ ≪ 1 and γ₁/Ω ≪ 1. This projection results in the effective amplitude of the inelastic two-phonon process,

$$\begin{equation} R_{{\rm eff}}\approx \frac{e^{2}\hbar^{2}v^{2}}{\Omega^{2}}\left\{ -\mathrm{i}\sigma_{z} ( \mathbf{l}\times\tilde{\mathbf{l}}^{*})_{z} + \frac{\gamma_{1}}{\Omega} \left[ \sigma_{x}d_{y} + \sigma_{y}d_{x} \right] \right\}. \end{equation} \tag{ 19 }$$

It is important to stress that for Ω ≫ γ₁, R_eff cannot be correctly obtained within an effective two-band theory, equation (3), without the above-described intermediate steps, although for Ω < γ₁ one could formally use H_eff and equations (4)–(7) to describe inelastic scattering of far-infrared (FIR) light. Within these approximations the scattering probability

$$\begin{equation} w \approx e^{4}\hbar v^{2} \frac{\gamma_{1}}{\Omega^{4}}\left\{\Xi_{s} + \frac{\gamma_{1}^{2}}{2\Omega^{2}} \Xi_{\mathrm{o}} \right\} \end{equation} \tag{ 20 }$$

retains the characteristics of monolayer graphene—the crossed polarization of in/out photons in the dominant mode [48]. The angle-integrated spectral density of Raman scattering

$$\begin{equation} g(\omega) = \frac{1}{4}\left(\frac{e^{2}}{\pi\hbar c}\frac{v}{c}\right)^{2}\frac{\gamma_{1}}{\Omega^{2}} \left\{\Xi_{s} + \frac{\gamma_{1}^{2}}{2\Omega^{2}} \Xi_{\mathrm{o}} \right\}\theta(\omega-2\mu) \end{equation} \tag{ 21 }$$

reflects the constant density of states characteristic of the parabolic spectrum of the bilayers. This is different from monolayer graphene, where g(ω)∝ω, reflecting the linear energy dependence of the density of states [47].

Upon quantization of electron states into LLs, Raman scattering inevitably causes inter-LL transitions, which, having fixed energies, manifest themselves in the electronic contribution to the Raman plot as a pronounced structure that can be used to detect their contribution experimentally [49, 50]

3.1. Monolayer graphene

The electronic spectrum of monolayer graphene in a strong magnetic field can be described as a sequence n^α of LLs, corresponding to the states [7, 40]

$$\begin{eqnarray} &&|n^{\alpha }\rangle =\frac{1}{\sqrt{2}} \left({\Phi_{n} \atop \mathrm{i}\alpha \Phi_{n-1}}\!\right)\quad \mbox{for}\; n\geqslant 1, \quad |0\rangle =\left(\!{\Phi _{0} \atop 0}\!\right),\\ &&\mbox{with energy}\;\alpha \varepsilon_{n},\quad \varepsilon_{n}=\sqrt{2n}\,\hbar v/\lambda_{B},\nonumber \end{eqnarray} \tag{ 22 }$$

where $\lambda _{B}=\sqrt {\hbar c/eB}$ is the magnetic length, n enumerates the LLs, α = +  denotes the conduction band and α = − the valence band and Φ_n are the normalized envelopes of LL wavefunctions. Then, the corresponding elements of Feynman diagrams in figure 2 are

Here $\mathbf {e}_{\pm }=\frac {1}{\sqrt {2}} ( \mathbf {e}_{x}\pm \mathrm {i}\mathbf {e}_{y})$ is used to stress that a circularly polarized photon carries angular momentum m = ± 1. The matrix structure of the interaction terms allows the following selection rules for optically active inter-LL excitations in monolayer graphene:

$$\begin{equation} \begin{array}{@{}llll@{}} ({\rm i}) & n^{-}\rightarrow n^{+}; & {\rm (iii)} & n^{-}\rightarrow (n+1)^{+}, \\ ({\rm ii}) & (n\mp 1)^{-}\rightarrow (n\pm 1)^{+}; & & (n+1)^{-}\rightarrow n^{+}. \end{array} \end{equation} \tag{ 23 }$$

The excitation of the e–h pairs by Raman scattering in graphene at strong magnetic fields produces the electronic transition (i) between LLs, with angular momentum transfer Δm = 0 and excitation energy ω = 2ε_n, see figure 3, and transitions (ii), with Δm = ± 2 and ω = ε_n−1 + ε_n+1. The amplitudes of these two processes,

$$\begin{eqnarray*} &&\fl R_{n^{-}\rightarrow n^{+}}=\frac{1}{4}\frac{(ev\hbar )^{2}}{c^{2}\Omega} \sum_{\alpha =\pm }\left[\! \frac{(\mathbf{le}_{+})({\tilde{\bf l}}^{*}\mathbf{e}_{-})}{\Omega\!-\!\varepsilon_{n}\!-\!\alpha\varepsilon_{n+1}} \!-\! \frac{(\mathbf{le}_{+})({\tilde{\bf l}}^{*}\mathbf{e}_{-})}{\varepsilon_{n}\!-\!\Omega\!-\!\alpha\varepsilon_{n-1}}- \frac{(\mathbf{le}_{-})({\tilde{\bf l}}^{*}\mathbf{e}_{+})}{\Omega-\varepsilon_{n}-\alpha\varepsilon_{n+1}} + \frac{(\mathbf{le}_{-})({\tilde{\bf l}}^{*}\mathbf{e}_{+})}{\varepsilon_{n}-\Omega-\alpha\varepsilon_{n-1}}\!\right], \end{eqnarray*}$$

$$\begin{eqnarray*} &&\fl R_{(n\mp 1)^{-}\rightarrow (n\pm 1)^{+}}=\mp \frac{1}{4}\frac{(ev\hbar )^{2}}{c^{2}\Omega } (\mathbf{le}_{\pm })({\tilde{\bf l}}^{*}\mathbf{e}_{\pm }) \sum_{\alpha =\pm }\left[ \frac{\alpha}{\Omega-\varepsilon_{n+1}-\alpha\varepsilon_{n}} + \frac{\alpha }{\varepsilon_{n-1}-\Omega-\alpha\varepsilon_{n}}\right], \end{eqnarray*} \noindent$$

are such that R_n⁻→n⁺ ≫ R_{(n∓1)⁻→(n±1)⁺} for ω ≪ Ω. The latter relation is determined by a partial cancellation of the two diagrams constituting R_D. Note that the inter-LL modes n⁻ → n⁺ have the symmetry of the representation A₂ in table 2 and the same circular polarization of 'in' and 'out' photons involved in its excitation. Finally, the contact term R_w in figure 2 allows for a weak transition (iii) n⁻ → (n ± 1)⁺, with the amplitude R_w ≪ R_n⁻→n⁺. This transition, together with selection rules Δm = ± 1, resembles the inter-LL transition involved in the FIR absorption [3, 7]. However, the FIR-active excitation is 'valley-symmetric' [7, 51] and belongs to the representation E₁, whereas the Raman-active n⁻ → (n ± 1)⁺ mode belongs to E₂, allowing the latter to couple to the Γ-point optical phonon and, thus, leading to the magneto-phonon resonance feature in the Raman spectrum [29]. Also, R_w originates from the trigonal warping term in $\mathcal {H}$, which violates the rotational symmetry of the Dirac Hamiltonian by transferring angular momentum ±3 from electrons to the lattice, allowing thus a change of angular momentum by ±1.

Zoom In Zoom Out Reset image size

Table 2. Raman-active inter-LL excitations in graphene.

C6v rep.	Transition	Intensity	Polarization
E2	$\displaystyle n^{-}\rightarrow (n+1)^{+}$ $\displaystyle (n+1)^{-}\rightarrow n^{+}$	Weak in Raman, strong in magnetophonon  resonance	σ±→σ∓
A1	${\displaystyle (n-1)^{-}\rightarrow (n+1)^{+} \atop \displaystyle (n+1)^{-}\rightarrow (n-1)^{+}}$	Weak in Raman	σ±→σ∓
A2	n−→n+	Dominant in Raman	σ±→σ±

The dominant inter-LL transitions n⁻ → n⁺ determine the spectral density of light scattered from electronic excitations in graphene at high magnetic fields:

$$\begin{equation} g_{n^{-}\rightarrow n^{+}}(\omega )\approx \Xi_{s}\left( \frac{v^{2}}{c^{2}}\frac{e^{2}/\lambda _{B}}{\pi \Omega }\right) ^{2} \sum_{n\geqslant 1}\gamma_{n}(\omega-2\varepsilon_{n}). \end{equation} \tag{ 24 }$$

Here, we use Lorentzian γ_n(x) = π⁻¹Γ_n/[x² + Γ²_n], and inelastic LL broadening Γ_n, which increases with the LL number, and the factor Ξ_s in equation (24) indicates that 'in' and 'out' photons have the same circular polarization.

3.2. Bilayer graphene

The same as in the non-magnetic case, e–h excitations in bilayer graphene are created in low-energy LLs, as at high energies the LL broadening due to, for example, electron–phonon interaction will smear out the LL spectrum. In strong magnetic fields, low-energy LLs in bilayer graphene can be described [52] by

$$\begin{eqnarray} &&|n^{\alpha}\rangle = \frac{1}{\sqrt{2}}\left({\Phi_{n} \atop \alpha\Phi_{n-2}}\!\right)\quad \mbox{for}\ \ n\geqslant 2, \quad |0/1\rangle = \left(\!{\Phi_{0/1} \atop 0}\!\right),\\ &&\mbox{with energy } \epsilon_{n^{\alpha}}=2\alpha\frac{\hbar^{2}v^{2}}{\gamma_{1}\lambda_{B}^{2}}\sqrt{n(n-1)}.\nonumber \end{eqnarray} \tag{ 25 }$$

In a neutral bilayer, all LLs have additional fourfold degeneracy (two due to the electron spin and two due to the valley). The peculiar n = 0 and n = 1 LLs are degenerate at  = 0, producing an eightfold degenerate level. Projecting the effective transition amplitude R_eff onto the eigenstates |n^α〉, we find the electronic Raman spectrum. This procedure leads to the same selection rules as in monolayer graphene; see equations (23). Among these, (i) are the dominant transitions. For an undoped bilayer, the angle-integrated spectral density g(ω) of Raman scattering in the strong magnetic field is equal to:

$$\begin{eqnarray} &&\fl g(\omega) \approx \left( \frac{v^{2}}{c^{2}}\frac{e^{2}/\lambda _{B}}{\pi \Omega }\right) ^{2} \Biggl\{\Xi_{s} \sum_{n\geqslant 2}\gamma(\omega - 2\epsilon_{n^{+}}) +\Xi_{\mathrm{o}}\left( \frac{\gamma_{1}}{\Omega} \right)^{2} \nonumber\\ &&\times \Biggl(\sum_{n=1,2}\gamma(\omega-\epsilon_{(n+1)^{+}}) + \frac{1}{2}\sum_{n\geqslant 3}\gamma(\omega-\epsilon_{(n+1)^{+}} - \epsilon_{(n-1)^{+}})\Biggr)\Biggl\}. \end{eqnarray} \tag{ 26 }$$

The first term describes the dominant contribution due to the n⁻ → n⁺ transitions, whereas the last two terms describe the spectral density of the (n ∓ 1)⁻ → (n ± 1)⁺ transitions.

Low-energy electronic contribution to the Raman plot is shown by a solid line in figure 3(b). The dominant peaks are due to the n⁻ → n⁺ transitions. The quantum efficiency of a single n⁻ → n⁺ peak in figure 3(b) is similar to the value I₁ in monolayer graphene [47].

A weaker feature in figure 3(b) corresponds to both 2⁻ → 0 and 0 → 2⁺ transitions and is visible, as it is positioned to the left of the first 2⁻ → 2⁺ peak. The quantum efficiencies of the (n ± 1)⁻ → (n ∓ 1)⁺ transitions are smaller by the factor $\left (\frac {\gamma _{1}}{\Omega }\right )^{2}$ in comparison with the n⁻ → n⁺ transitions. This is different from the monolayer graphene case, where the corresponding ratio between quantum efficiencies of (n ± 1)⁻ → (n ∓ 1)⁺ and n⁻ → n⁺ transitions is $\left (\frac {\omega }{\Omega }\right )^{2}$, much smaller than for the bilayer.

Single and multiple phonon-emission lines in the Raman spectrum of graphene and the influence of electron–phonon coupling on the phonon spectrum have been investigated in great detail in [24–32, 53]. The coupling of electromagnetic field and phonon excitations is realized via the excitation of the virtual e–h pair, with the same symmetry properties as the Γ-point optical phonon that is created upon recombination of the pair. Analysing all possible combinations of polarization vectors, one may see from table 1 that the amplitude of Raman scattering with excitation of the phonon mode is proportional to

$$\begin{equation} R = R_{G} \sum_{\nu}\mathbf{u}_{\nu}\cdot\mathbf{d} \end{equation} \tag{ 27 }$$

that describes quantized sublattice displacement

$$\begin{equation} \mathcal{U}(\mathbf{r},t)=(\mathcal{U}_{A}-\mathcal{U}_{B})/\sqrt{2} = \sum_{\mathbf{k},\nu}\frac{\hbar}{\sqrt{2M\omega_{0}}} \left( \mathbf{u}_{\nu} c_{\nu,\mathbf{k}} \,{\rm e}^{{\rm i}(\mathbf{k}\mathbf{r}-\omega_{0}t)/\hbar} + \mathrm{h.c.} \right), \end{equation} \tag{ 28 }$$

where M is the mass of a single carbon atom, u is the polarization of the phonon, and the expression for R_G can be taken from Basko [53]. Raman efficiency then is proportional to

$$\begin{equation} g \propto \sum_{\upsilon\nu} d_{\upsilon}d_{\nu}^{*} \, \mathrm{Im}\,\tilde{D}^{A}_{\upsilon\nu}, \end{equation} \tag{ 29 }$$

where $\tilde D$ is the phonon propagator renormalized by electron–phonon interaction. The latter strongly affects the energy of the phonon under the conditions of its resonance mixing with the inter-LL electronic transitions, an effect known as magnetophonon resonance.

We use both the phenomenological approach and the tight-binding model to study the interaction of electrons with Γ-point optical phonons in the case of a strong magnetic field. The phenomenology is suitable for describing the low-energy excitations of electrons in the vicinity of Dirac points (K-points), while the tight-binding model is needed to perform the integration over the whole BZ [53] when estimating the intensity of the G-phonon line in Raman. To discuss the electron–phonon interaction it is instructive to trace the electron Hamiltonian

$$\begin{equation} H=H_{0}-\mathbf{Q}\cdot\mathcal{U} \end{equation} \tag{ 30 }$$

back to the tight-binding model

$$\begin{eqnarray} H_{0} &= \gamma_{0} \sum_{s} \left(\!\begin{array}{@{}cc@{}} 0 & {\rm e}^{{\rm i} \mathbf{k}\mathbf{r}_{s}} \\ \mathrm{e}^{-\mathrm{i} \mathbf{k}\mathbf{r}_{s}} & 0 \end{array}\!\right) \stackrel{\mathbf{k\rightarrow \xi K+p}}{\longrightarrow} \xi v \, \boldsymbol{\sigma}\mathbf{p}, \end{eqnarray} \tag{ 31 }$$

$$\begin{eqnarray} \mathbf{Q} &= \sqrt{6}\frac{\mathrm{d} \gamma_{0}}{\mathrm{d} r} \sum_{s} \left(\begin{array}{@{}cc@{}} 0 & {\rm e}^{{\rm i} \mathbf{k}\mathbf{r}_{s}} \\ \mathrm{e}^{-\mathrm{i} \mathbf{k}\mathbf{r}_{s}} & 0 \end{array}\right) \frac{\mathbf{r}_{s}}{a} \stackrel{\mathbf{k\rightarrow \xi K+p}}{\longrightarrow} F \, \boldsymbol{\sigma}\times\mathbf{e}_{z}, \end{eqnarray} \tag{ 32 }$$

where k is a momentum within the BZ, $\mathbf {K}= (-\frac {4\pi }{3a},0)$, $v=\frac {\sqrt 3}{2} \gamma _{0} a /\hbar$, the vectors $\mathbf {r}_1=(-\frac {a}{2},\frac {a}{2\sqrt 3})$, $\mathbf {r}_2=(\frac {a}{2},\frac {a}{2\sqrt 3})$ and $\mathbf {r}_3=(0,-\frac {a}{\sqrt 3})$ are the positions of the nearest 'B' atoms with respect to 'A' atom, and $F=\frac {3}{\sqrt 2}\,\frac {\partial \gamma _{0}}{\partial r}$, where $\frac {\partial \gamma _{0}}{\partial r}$ is a response of the hopping amplitude to the variation of the bond length.

The renormalization of the phonon energy by its interaction with electrons is taken into account in the random phase approximation. Defining electrons' bare Green's functions phonon bare propagator, and electron–phonon interaction vertex,

we use the Dyson equation to find the renormalized phonon propagator:

Here $\lambda =\frac {\sqrt 3}{M\omega _0}\left (\frac {Fa}{2v}\right )^2= \frac {3\sqrt 3}{2M\omega _0}\left (\frac {\hbar }{\gamma _0}\frac {\mathrm {d} \gamma _0}{\mathrm {d} r}\right )^{2}\approx 0.018$ is dimensionless constant of the electron–phonon coupling [26]. Trace is taken in AB sites space and a summation over spins is already performed. As the momentum k of a phonon created by the Raman process is small by parameter v/c, in the following it can be neglected.

The integral for the polarization loop in equation (33) is divergent in the frame of the linear Dirac model, so it can be split into two parts:

$$\begin{equation} \Pi_{\upsilon\nu}^A = \Pi_{0} \delta_{\upsilon\nu} + \Delta\Pi_{\upsilon\nu}, \quad \Pi_{0}\delta_{\upsilon\nu}=\Pi_{\upsilon\nu}^A|_{\omega,\mu,B\rightarrow 0}\approx - 0.64\gamma_{0}\delta_{\upsilon\nu}. \end{equation} \tag{ 35 }$$

Here Π₀ is calculated numerically within the tight-binding model, and ΔΠ_υν is determined within the effective Dirac model, after subtracting the divergent part from the integrand. For B = 0 the answer is ΔΠ^A_υν∝δ_υν and it repeats the derivation of Raman line shape found in  [25, 28–30, 32, 54, 55]:

$$\begin{equation*} g = \Xi_{\mathrm{o}} \, I_{G} \frac{\frac{1}{\pi}\gamma_{G}}{(\omega-\omega_{G})^{2}+\gamma_{G}^{2}}, \end{equation*}$$

where the width of the G-peak and its position are determined by

$$\begin{eqnarray} &&\gamma_{G}=\gamma_{\rm ph} - \frac{\lambda}{4} \tilde{\omega}_{0} \, \theta(\tilde{\omega}_{0}-2\mu), \end{eqnarray} \tag{ 36 }$$

$$\begin{eqnarray} &&\omega_{G}=\tilde{\omega}_{0} - \lambda\frac{\mu}{\pi} + \lambda\frac{\tilde{\omega}_{0}}{4\pi} \left( \log\left|\frac{2\mu+ \tilde{\omega}_{0}}{2\mu-\tilde{\omega}_{0}}\right| \right) \end{eqnarray} \tag{ 37 }$$

where $\tilde {\omega }_{0}=\omega _0 + 0.64\lambda \gamma _{0}$ (we implied that ω₀ ≫ 0.64λγ₀ ≫ λε_F), and the expression for I_G can be taken from Basko [53].

4.1. Magnetophonon resonance in a homogeneously doped monolayer

In the presence of a magnetic field the continuous electronic spectrum splits into the set of discrete LLs. Working in the basis of LL equations (22)–(34) and choosing the circularly polarized basis $\mathbf {e}_{\circlearrowright /\circlearrowleft }= \frac {1}{\sqrt 2}(\mathbf {e}_{x} \pm \mathrm {i} \mathbf {e}_{y})$ for incoming and scattered light, we find that

$$\begin{equation*} \Pi_{\upsilon\nu} = (\mathbf{e}_{\circlearrowright})_{\upsilon}(\mathbf{e}_{\circlearrowright})_{\nu}^{*} \Pi_{\circlearrowright} + (\mathbf{e}_{\circlearrowleft})_{\upsilon}(\mathbf{e}_{\circlearrowleft})_{\nu}^{*} \Pi_{\circlearrowleft}, \mbox{ where } \Pi_{\circlearrowright/\circlearrowleft}=\Pi_{0}+\Delta\Pi_{\circlearrowright/\circlearrowleft}. \end{equation*}$$

For the chemical potential pinned at the zero LL,

$$\begin{equation} \Delta\Pi_{\circlearrowright/\circlearrowleft} = \frac{\varepsilon_{1}^{2}}{\pi} \Biggl[ \frac{-f_{0}}{\varepsilon_{1}\pm \omega} + \frac{1-f_{0}}{-\varepsilon_{1}\pm \omega} + \frac{1}{\varepsilon_{1}} -\omega^{2} \sum_{n\geqslant 2}^{\infty} \frac{1}{\Omega_{n,+}(\Omega_{n,+}^{2}-\omega^{2})} \Biggr], \end{equation} \tag{ 38 }$$

where f₀ is a density-dependent filling factor of the zero level and Ω_n,± = ε_n ± ε_n−1 are the frequencies of e–h excitations. If μ is large enough, so that the mth level (m ⩾ 1) is partially filled, then

$$\begin{eqnarray} &&\fl \Delta\Pi_{\circlearrowright/\circlearrowleft} = \frac{\varepsilon_{1}^{2}}{2\pi} \Biggl[ \frac{f_{m^+}-1}{\Omega_{m,-}\pm \omega} + \frac{-f_{m^+}}{\Omega_{m+1,-}\pm \omega} + \frac{f_{m^+}-1}{\Omega_{m,+}\pm \omega}\nonumber\\ &&+ \frac{1-f_{m^+}}{-\Omega_{m+1,+}\pm \omega} + \frac{-1}{\Omega_{m,+}\pm \omega} + 2\frac{\varepsilon_{m+1}}{\varepsilon_{1}^{2}} -2\omega^{2} \sum_{n\geqslant m+2}^{\infty} \frac{1}{\Omega_{n,+}(\Omega_{n,+}^{2}-\omega^{2})} \Biggr].\qquad \end{eqnarray} \tag{ 39 }$$

Then, the net Raman efficiency of the Γ-point phonon mode in graphene in a strong magnetic field is

$$\begin{equation} g = \sum_{\mathcal{A}= \circlearrowleft/\circlearrowright}|\mathbf{d}\mathbf{e}_{\mathcal{A}}|^{2} \frac{I_{G}\frac{1}{\pi}\gamma_{\rm ph}}{(\omega-\tilde{\omega}_{0}-\lambda\,\Delta\Pi_{\mathcal{A}}(\omega))^{2}+\gamma_{\rm ph}^{2}}, \end{equation} \tag{ 40 }$$

where the two terms in the sum correspond to the opposite circular polarizations of the optical phonon, and it is plotted in figure 4 for several doping densities.

Zoom In Zoom Out Reset image size

The line described by equation (40) displays two types of typical behaviour depending on the range of parameters. Far from the resonance only a slight shift of the electron and phonon excitations exists,

$$\begin{equation} \omega_{\circlearrowleft/\circlearrowright} \approx \tilde{\omega}_{0} + \lambda \, \Delta\Pi_{\circlearrowleft/\circlearrowright}(\tilde{\omega}_{0}). \end{equation} \tag{ 41 }$$

Since the line shift may differ for the opposite circular polarizations of the phonon, this leads to a small splitting of the G-line in Raman. In the vicinity of the intersection of the phonon line with the e–h excitation line for several doping densities, the polarization operator becomes large and causes significant change in the spectra, forming an anticrossing, in addition to the line splitting. The form of the line is determined by four solutions $\omega_{\circlearrowleft/\circlearrowright}^{\pm}$

$$\begin{equation} \omega - \tilde{\omega}_{0} - \lambda\,\Delta\Pi_{\circlearrowleft/\circlearrowright}(\omega)=0, \end{equation} \tag{ 42 }$$

two for each polarization with energies above (+) or below (−) the bare phonon line. The form of such anticrossings in undoped graphene is shown in figure 4. The resonance takes a more sophisticated form in the case of the doped graphene. With increasing the electron density (and correspondingly the filling factor), the transition 1⁻ → 0⁺ invoked by anticlockwise polarized light is weakened by Pauli blocking (see the first anticrossing in figure 5(a)), while the 0⁻ → 1⁺ transition remains. When the n = 0 LL is fully occupied (n_e = 1/πλ²_B ≈ 1.4 × 10¹² cm⁻¹²), in anticlockwise polarization the anticrossing is completely blocked (see figure 5(b)), which coincides with the experimentally observed features [49].

Zoom In Zoom Out Reset image size

4.2. Influence of inhomogeneous doping and strain on the magnetophonon resonance line shape

4.3. Inhomogeneous density

4.4. Inhomogeneously strained graphene

Inhomogeneous strain in graphene is often associated with magnetic field, B_strain for electrons in K valley and −B_strain in K', which is estimated as $B_{\mathrm {strain}}\sim \frac {\hbar c}{e a} \frac {h^{2}}{l^{3}}$, where a is the bond length and l and h are the linear scale and height of a ripple in graphene [34]. In the external magnetic field B, the total effective field acting on electrons in the valleys K and K' is different,

$$\begin{equation*} B_{\xi}=B+\xi B_{\mathrm{strain}}, \end{equation*}$$

where ξ = +  for K point and ξ = − for K'. The difference in magnetic field leads to different Landau spectra in each valley, thus splitting the values of the external field B, for which, e.g., the resonance of 0⁻_K → 1⁺_K and 0⁻_K' → 1⁺_K' with the phonon is realized. Formally, this can be taken into account by a substitution $\Delta \Pi (B)\rightarrow \frac {1}{2}\sum _{\xi =\pm }\Delta \Pi (B_{\xi })$ in equation (40).

For the local Raman spectrum related to an area with a fixed B_strain, this results in a fine structure of the anticrossing shown in figure 6(a), where an additional line appears seemingly crossing the phonon mode without splitting. Such a behaviour occurs due to the fact that instead of an anticrossing/mixing of two modes, now three modes are coupled, forming one phonon–electron mixed excitation with energy pinned in the middle of two strongly split side bands. Here, the central phonon line corresponds to the circular polarization of the phonon for which mixing with electronic modes is suppressed by Pauli blocking, due to the high density of electrons used in the calculation. Averaging out a normal distribution of B_strain with zero average expectation and variance δB_strain = 3 T, the anticrossing branches become asymmetrically broadened, as shown in figure 6(b).

Zoom In Zoom Out Reset image size

We have presented a theory of the inelastic scattering of photons in mono- and bilayer graphene accompanied by an excitation of e–h pairs. The dominant scattering processes are characterized by crossed polarization of in/out photons, while the Raman efficiency is proportional to the density of states (a linear increase with Raman shift for monolayer and constant for bilayer) and in the case of doping a gap of width 2μ opens up.

In a strong magnetic field, Raman spectra acquire a pronounced structure with peaks corresponding to the inter-LL transitions. The selection rules for the strongest mode preserve LL number n: n⁻ → n⁺, and the quantum efficiency of the single peak in figure 3 is $I_{1}\sim ( \frac {v^{2}}{c^{2}}\,\frac {e^{2}/\lambda _{B}}{\pi \Omega } )^{2}\sim 10^{-12}$ per incoming photon for Ω ∼ 1 eV photons in the magnetic field B ∼ 10 T. The theoretical prediction for the intensity of the phonon-induced G peak [53], which is a well-known Raman feature in carbon materials [54], estimates I_G ∼ 10⁻¹¹. This result is close in order of magnitude to the intensity of a single n⁻ → n⁺ peak. Hence, the spectral features of inter-LL transitions in graphene can be observed experimentally [49]. The weaker mode with the same selection rules, similar to the absorption process, n⁻ → (n ± 1)⁺, but belonging to different symmetry group representations (E₂, while the absorption has E₁), is coupled to the phonon mode, creating the anticrossing and the G-peak line.

We thank A Ferrari, D Smirnov, M Mucha-Kruczynski, K Kechedzhi, M Goerbig, F Guinea, T Ando and D Basko for useful discussions. This reviewed work was supported by the EPSRC, the Royal Society, the EC ERC, the Royal Society's Wolfson Research Merit Award and the European Research Council's Advance Investigator grant 'Graphene and Beyond'.