Ab initio calculation of the $G$ peak intensity of graphene: Combined study of the laser and Fermi energy dependence and importance of quantum interference effects

TL;DR

This study uses ab initio calculations to analyze the $G$ peak intensity in graphene, highlighting the roles of quantum interference and non-resonant states across various laser and Fermi energies.

Contribution

It introduces a general, ab initio approach that surpasses previous analytical models, emphasizing quantum interference effects in the $G$ peak of graphene.

Findings

Quantum interference significantly affects the $G$ peak intensity.

Non-resonant electronic states contribute notably to the $G$ peak.

Methodology is applicable to other materials.

Abstract

We present the results of a diagrammatic, fully ab initio calculation of the $G$ peak intensity of graphene. The flexibility and generality of our approach enables us to go beyond the previous analytical calculations in the low-energy regime. We study the laser and Fermi energy dependence of the $G$ peak intensity and analyze the contributions from resonant and non-resonant electronic transitions. In particular, we explicitly demonstrate the importance of quantum interference and non-resonant states for the $G$ peak process. Our method of analysis and computational concept is completely general and can easily be applied to study other materials as well.

Figures

Equations

$$\frac{\mathrm{d}{I}_G}{\mathrm{d}\Omega }\propto \frac{{\omega }_{\mathrm{out}}^2}{(2\pi {)}^2{c}^4}|\mathcal{M}{|}^2\times \delta ({\omega }_{\mathrm{in}}-{\omega }_{\mathrm{out}}-{\omega }_{\mathrm{ph}}),$$

$$i\mathcal{M}2\pi \delta ({\omega }_{\mathrm{in}}-{\omega }_{\mathrm{out}}-{\omega }_{\mathrm{ph}})\equiv ⟨{\omega }_{\mathrm{out}},\nu ;{\omega }_{\mathrm{ph}},\lambda |i\hat{T}|{\omega }_{\mathrm{in}},\mu ⟩$$

$$i\mathcal{M}\times 2\pi \delta ({\omega }_{\mathrm{in}}-{\omega }_{\mathrm{out}}-{\omega }_{\mathrm{ph}})={\left[\underset{t\to (\infty -i\epsilon )}{\lim }⟨{\omega }_{\mathrm{out}},\nu ;{\omega }_{\mathrm{ph}},\lambda |\mathcal{T}\exp \left(-i{\int }_{-t}^{+t}\mathrm{d}{t}^{\prime }{\hat{H}}_{\mathrm{int}}(t^{\prime })\right)|{\omega }_{\mathrm{in}},\mu ⟩\right]}_{\begin{array}{c}\text{connected,}\\ \text{amputated}\\ \text{diagrams only}\end{array}},$$

$$\begin{array}{rl}i{\mathcal{M}}_1& =\sum _{\text{spin}}\sum _{\mathbf{k}}\int \frac{\mathrm{d}\omega }{2\pi i}{\mathrm{e}}^{i\omega 0^{+}}\mathrm{tr}\left[G_{\mathbf{k}}(\omega ){\gamma }_{\mathbf{k}}^{\mu }{G}_{\mathbf{k}}(\omega -{\omega }_{\mathrm{in}}){\left(g_{\mathbf{k}}^{\lambda }\right)}^{†}{G}_{\mathbf{k}}(\omega -{\omega }_{\mathrm{out}}){\left({\gamma }_{\mathbf{k}}^{\nu }\right)}^{†}\right]\\ i{\mathcal{M}}_2& =\sum _{\text{spin}}\sum _{\mathbf{k}}\int \frac{\mathrm{d}\omega }{2\pi i}{\mathrm{e}}^{i\omega 0^{+}}\mathrm{tr}\left[G_{\mathbf{k}}(\omega ){\left({\gamma }_{\mathbf{k}}^{\nu }\right)}^{†}{G}_{\mathbf{k}}(\omega +{\omega }_{\mathrm{out}}){\left(g_{\mathbf{k}}^{\lambda }\right)}^{†}{G}_{\mathbf{k}}(\omega +{\omega }_{\mathrm{in}}){\gamma }_{\mathbf{k}}^{\mu }\right].\end{array}$$

$$i{\mathcal{M}}_{\mathbf{k}}=i{\mathcal{M}}_{\mathbf{k}}^{\mathrm{aDR}}+i{\mathcal{M}}_{\mathbf{k}}^{\mathrm{SR}}+i{\mathcal{M}}_{\mathbf{k}}^{\mathrm{NR}},$$

$$i\mathcal{M}_{\mathbf{k}}^{\mathrm{aDR}}=\frac{\Big{(}\left(\gamma^{\nu}_{\mathbf{k}}\right)^{\dagger}\Big{)}_{\pi,\pi^{*}}\Big{(}\left(g^{\lambda}_{\mathbf{k}}\right)^{\dagger}\Big{)}_{\pi^{*},\pi^{*}}\Big{(}\gamma^{\mu}_{\mathbf{k}}\Big{)}_{\pi^{*},\pi}}{\left[\omega_{\mathrm{in}}-\Delta\varepsilon_{\mathbf{k}}^{\pi^{*},\pi}+i\bar{\gamma}_{\mathbf{k}}^{\pi^{*},\pi}\right]\left[\omega_{\mathrm{out}}-\Delta\varepsilon_{\mathbf{k}}^{\pi^{*},\pi}+i\bar{\gamma}_{\mathbf{k}}^{\pi^{*},\pi}\right]}$$

$$i\mathcal{M}_{\mathbf{k}}^{\mathrm{SR}}=\sum_{s=\pi,\pi^{*}}\left\{\frac{\Big{(}\left(\gamma^{\nu}_{\mathbf{k}}\right)^{\dagger}\Big{)}_{\pi,\pi^{*}}\Big{(}\gamma^{\mu}_{\mathbf{k}}\Big{)}_{\pi^{*},s}\Big{(}\left(g^{\lambda}_{\mathbf{k}}\right)^{\dagger}\Big{)}_{s,\pi}}{\left[\omega_{\mathrm{out}}-\Delta\varepsilon_{\mathbf{k}}^{\pi^{*},\pi}+i\bar{\gamma}_{\mathbf{k}}^{\pi^{*},\pi}\right]\left[-\omega_{\mathrm{ph}}-\Delta\varepsilon_{\mathbf{k}}^{s,\pi}+i\bar{\gamma}_{\mathbf{k}}^{s,\pi}\right]}+\frac{\Big{(}\left(g^{\lambda}_{\mathbf{k}}\right)^{\dagger}\Big{)}_{\pi,s}\Big{(}\left(\gamma^{\nu}_{\mathbf{k}}\right)^{\dagger}\Big{)}_{s,\pi^{*}}\Big{(}\gamma^{\mu}_{\mathbf{k}}\Big{)}_{\pi^{*},\pi}}{\left[\omega_{\mathrm{ph}}-\Delta\varepsilon_{\mathbf{k}}^{s,\pi}+i\bar{\gamma}_{\mathbf{k}}^{s,\pi}\right]\left[\omega_{\mathrm{in}}-\Delta\varepsilon_{\mathbf{k}}^{\pi^{*},\pi}+i\bar{\gamma}_{\mathbf{k}}^{\pi^{*},\pi}\right]}\right\}$$

$$\mathcal{R}(x)=\left(\begin{array}{cc}0& c\\ c& 0\end{array}\right),\mathcal{R}(y)=\left(\begin{array}{cc}c& 0\\ 0& -c\end{array}\right),$$

$$\left|\Delta {\epsilon }_{\mathbf{k}}-\frac{{\omega }_{\mathrm{in}}+{\omega }_{\mathrm{out}}}{2}\right|<{\epsilon }_{\mathrm{cut}},$$

Full text

I Introduction

Raman spectroscopy of graphene has been a subject of considerable interest as it is a fast and non-destructive way for sample characterization. The typical Raman spectrum of pristine graphene displays only two prominent peaks, the so-called $G$ and $2D$ peaks, which are the result of a scattering process involving one and two phonons, respectively. Despite this apparent simplicity, the Raman spectrum yields a large amount of information on, amongst others, doping, strain, inter-layer interaction, and the underlying substrate Ferrari (2007); Malard et al. (2009); Ferrari and Basko (2013); Reichardt and Wirtz (2017). To understand the influence of these quantities on the Raman spectrum of graphene, considerable effort has been made to understand the shape Basko (2008); Venezuela et al. (2011); Berciaud et al. (2013); Hasdeo et al. (2014); Neumann et al. (2016), width Graf et al. (2007); Pisana et al. (2007); Herziger et al. (2014); Neumann et al. (2015), height Cançado et al. (2007); Basko (2009); Basko et al. (2009), and position Maultzsch et al. (2004a); Ferrari et al. (2006); Berciaud et al. (2008); Pisana et al. (2007); Mohiuddin et al. (2009); Yoon et al. (2011) of the $G$ and $2D$ peaks.

While a clear picture for the $2D$ peak has been established Thomsen and Reich (2000); Maultzsch et al. (2004b); Venezuela et al. (2011); Herziger et al. (2014), a corresponding simple picture for the $G$ peak is still missing. The latest theoretical effort by Basko Basko (2009) focused on an analytical tight-binding approach in the low-energy regime. In particular, the question of which electronic transitions contribute to the $G$ peak has still not been satisfactorily addressed. This question is also intimately tied to the role of quantum interference effects in the $G$ peak process, which has been probed experimentally by the tuning of destructive interference effects via variation of the Fermi level Chen et al. (2011).

To address these open questions, we present results of a fully ab initio calculation of the one-phonon Raman intensity of graphene. Our general and flexible approach allows us to study the laser and Fermi energy dependence of the $G$ peak intensity and enables us to analyze and understand the contributions from resonant and non-resonant electronic transitions. In particular, we explicitly demonstrate the importance of quantum interference and non-resonant states for the $G$ peak process in the low-energy regime.

This paper is organized as follows:

In section II, we describe our approach to the calculation of one-phonon Raman intensities from first principles. Subsection II.1 contains the computational details of our ab initio calculations, while we give details on the diagrammatic approach to the calculation of the Raman matrix element in subsection II.2. This subsection also introduces the two main concepts needed for the subsequent discussion: quantum interference effects and resonant contributions to the Raman matrix element.

The results of our calculations are presented in section III. In particular we discuss: (A) the laser energy dependence of the $G$ peak intensity, (B) its Fermi energy dependence for different laser energies, and (C) which electronic transitions play a non-negligible role for the $G$ peak process.

II Theory

In order to compute one-phonon Raman intensities from first principles, we use a perturbative approach organized via Feynman diagrams to calculate the quantum mechanical amplitude, i.e., the scattering matrix element, for a one-phonon Raman process. The quantities entering the expression for the matrix element, such as the electronic band structure, the electron-light, and the electron-phonon coupling were obtained from ab initio calculations on the level of density functional theory (DFT), as detailed below.

II.1 Computational details

The electronic density in the ground state was calculated with the PWscf code as included in the Quantum ESPRESSO suite Giannozzi et al. (2009) using an ultrasoft pseudopotential. The exchange-correlation functional was approximated on the level of the generalized gradient approximation (GGA) in the parametrization of Perdew, Burke, and Ernzerhof (PBE) Perdew et al. (1996). Integrals over the first Brillouin zone (BZ) were carried out on a regular mesh of 60$\times$60$\times$1 $\mathbf{k}$-points, while the plane-wave cutoff was set to 80 Ry. These values lead to converged results for a vacuum spacing of 14 Å, which separates periodic copies of the graphene sheet in the $z$-direction. As graphene has semi-metallic character, occupations were smeared out using a 0.002 Ry Fermi-Dirac smearing. For the lattice constant, we use a value of 2.46 Å, obtained from structure relaxation.

The frequency of the doubly-degenerate optical phonon at $\Gamma$ was fixed at the experimentally Lee et al. (2012) obtained value of 1581.6 cm*-1* for pristine, freestanding graphene. Since we are interested in the intensity of the one-phonon Raman peak, we neglect the dependence of the phonon frequency on the Fermi level as this shift is only on the order of a few cm*-1* (see Ref. Pisana et al., 2007), corresponding to a change of 1-2 meV. This change is negligible compared to the dominating energy scale appearing in the expression for the Raman intensity set by the laser energy, which is on the order of 1-4 eV for visible light.

The electronic band structure, the electron-light coupling matrix element in the dipole approximation, and the electron-phonon coupling matrix elements were obtained on a coarse 12$\times$12$\times$1 $\mathbf{k}$-point mesh and then interpolated to a finer 480$\times$480$\times$1 $\mathbf{k}$-point grid with the help of maximally localized Wannier functions using modified forms of the Wannier90 Mostofi et al. (2014) and EPWGiustino et al. (2007); Poncé et al. (2016) codes. We only consider the $\pi$ and $\pi^{*}$ bands of graphene in the interpolation as relevant optical transitions to the $\sigma$-bands are forbidden by selection rules.

Obtaining the above mentioned quantities on a very fine $\mathbf{k}$-point grid is necessary to obtain converged results for the Raman matrix element. We checked the convergence of our calculations with respect to the fine $\mathbf{k}$-point mesh by testing mesh sizes ranging from 12$\times$12$\times$1 to 960$\times$960$\times$1 $\mathbf{k}$-points.

II.2 Diagrammatic calculation of Raman intensities

The intensity of the one-phonon Raman peak (the $G$ peak) per solid angle $\mathrm{d}\Omega$ can be calculated via a generalization of Fermi’s golden rule Basko (2009):

[TABLE]

where

[TABLE]

is a short-hand notation for the scattering matrix element. We focus on the case of Stokes scattering and factored out the energy-conserving $\delta$-function from the matrix element. For a phonon with finite lifetime, the $\delta$-function is to be replaced by a Lorentzian with a full width at half maximum of $\gamma_{\mathrm{ph}}$. In the expression above, the state $|\omega_{\mathrm{in}},\mu\rangle$ represents an incoming photon with frequency $\omega_{\mathrm{in}}$ and polarization $\mu$, while the state $|\omega_{\mathrm{out}},\nu;\omega_{\mathrm{ph}},\lambda\rangle\equiv|\omega_{\mathrm{out}},\nu\rangle\otimes|\omega_{\mathrm{ph}},\lambda\rangle$ is comprised of one outgoing photon of frequency $\omega_{\mathrm{out}}$ and polarization $\nu$ and one outgoing (in-plane) optical phonon with frequency $\omega_{\mathrm{ph}}$ and polarization $\lambda$. Since the photon momentum is negligible, we use its frequency rather than its momentum to label the state. All photons are assumed to have momentum along the $z$-direction, i.e., to be polarized in-plane, which is reasonable as the photon detector in an experiment is positioned on the axis of the incoming light perpendicular to the sample.

The operator $i\hat{T}=\hat{S}-1$ denotes the non-trivial part of the scattering matrix ($S$-matrix). Within the framework of diagrammatic, time-dependent perturbation theory, the scattering matrix element is given by Peskin and Schroeder (1995)

[TABLE]

where the subscript text indicates that, in an expansion of the exponential in terms of Feynman diagrams, only the fully connected diagrams are to be retained and the propagators for the three external lines (one each for the in- and outgoing photon and one for the phonon) are to be removed and replaced by the corresponding polarization vectors. The operator $\hat{H}_{\mathrm{int}}$ in the formula above represents the sum of the interaction parts of the Hamiltonian (technically taken to be in the interaction picture), i.e., it is the sum of the electron-photon and electron-phonon interaction Hamiltonian and $\mathcal{T}$ denotes the time-ordering symbol.

The leading order contributions in a diagrammatic expansion of the $T$-matrix element are shown in Fig. 1.

These diagrams represent the leading order terms in the independent particle picture, i.e., correlation effects beyond those included in the one-particle Green’s functions are neglected. In particular, we neglect excitonic effects and correlation contributions to the electron-phonon coupling besides those included already on the DFT level. For optical transitions in the visible regime, excitonic effects are assumed to be negligible in graphene. Electronic correlation effects were shown to give a large contribution to the electron-phonon coupling for phonons with momenta near $\Gamma$ and $K$ Lazzeri et al. (2008). Since here we are only interested in processes involving the phonons at $\Gamma$, these corrections can, in a first approximation, be taken into account by a simple rescaling of the electron-phonon coupling Lazzeri et al. (2008). However, since we will always compare the calculated Raman intensities or matrix elements to a calculated reference value (e.g., to an intensity at a certain $\omega_{\mathrm{in}}$), this correction factor will always cancel in the ratio and is thus not included here.

Explicitly, the contributions from the leading order Feynman diagrams read

[TABLE]

Here, $\gamma_{\mathbf{k}}^{\mu}$ and $g_{\mathbf{k}}^{\mu}$ are 2$\times$2 matrices in electronic band space whose entries are the electron-light and electron-phonon coupling matrix elements between two electronic states at $\mathbf{k}$ for a photon/phonon with polarization $\mu$. The matrix form of the electronic Green’s function is denoted by $G_{\mathbf{k}}(\omega)$, while “$\mathrm{tr}$” denotes the trace in the space of electronic bands. Since all coupling matrices and Green’s functions are diagonal in spin, the spin sum just yields a factor of 2. The sum over $\mathbf{k}$-points is understood to represent an integration over the entire first BZ, which we carry out by sampling the BZ with a regular mesh. Satisfactory convergence was achieved with a 480$\times$480$\times$1-mesh. Due to the fact that the total matrix element is a sum over $\mathbf{k}$-points, contributions from different parts of the BZ can interfere constructively or destructively. We will revisit this important point in section III, when we discuss the interpretation of our results.

The frequency integral can be evaluated via contour integration, where the factor $0^{+}$ prescribes the closure of the contour in the upper half of the complex plane. The application of the residue theorem leads to a sum of three expressions for each Feynman diagram, i.e., six terms in total, which correspond to the six terms obtained from an expansion in terms of time-ordered Goldstone diagrams (compare, e.g., Ref. Yu and Cardona, 2010 for explicit expressions of the six terms). After the frequency integration has been carried out, the contributions to the matrix element at a specific $\mathbf{k}$-point can conveniently be grouped into three classes:

[TABLE]

where the superscripts “aDR”, “SR”, and “NR” label the almost double-resonant, the single-resonant, and the non-resonant contribution to $\mathcal{M}_{\mathbf{k}}$, respectively. Here, “resonant” refers to resonance of an electronic transition with the frequency of the ingoing or outgoing light. For instance, the aDR-contribution reads (in the case of zero doping):

[TABLE]

This amplitude describes a process in which an electron is excited from the $\pi$- to the $\pi^{*}$-band and is subsequently scattered to an intermediate state in the $\pi^{*}$-band accompanied by emission of a phonon before it finally radiatively recombines with the hole it left behind. This contribution is maximal when the energy of an electronic transition $\Delta\varepsilon_{\mathbf{k}}^{\pi^{*},\pi}\equiv\varepsilon_{\mathbf{k}}^{\pi^{*}}-\varepsilon_{\mathbf{k}}^{\pi}$ equals the incoming or outgoing light frequency. In particular, for laser energies $\omega_{\mathrm{in}}\gg\omega_{\mathrm{ph}}$, we have $\omega_{\mathrm{in}}\approx\omega_{\mathrm{out}}$ and the aDR-term becomes quasi-double-resonant. By contrast, the SR-term

[TABLE]

can only become single-resonant as the second factor in the denominator involves the phonon frequency. In the expressions above, the quantity $\bar{\gamma}_{\mathbf{k}}^{s,s^{\prime}}$ is the average of the decay widths of the states $|s,\mathbf{k}\rangle$ and $|s^{\prime},\mathbf{k}\rangle$. For simplicity, we use a constant value of 100 meV for the decay width of one state.

All other terms appearing in $i\mathcal{M}_{\mathbf{k}}$ are summarily included in $i\mathcal{M}_{\mathbf{k}}^{\mathrm{NR}}$. These include terms that show resonant behavior with respect to the phonon frequency. However, the set of electronic transitions form an area in $\mathbf{k}$-space that is much smaller than the corresponding area for the transitions that are in resonance with respect to the incoming or outgoing light and its weight in the sum over all $\mathbf{k}$-points is negligible and hence we do not include these terms in the single-resonant category. We also want to point out that, at least in the case of graphene, the sum of the non-resonant contributions is by no means negligible, as demonstrated in the next section.

III Results and discussion

The calculation of the Raman intensity using Feynman diagrams allows us to investigate the origin of various experimentally established facts about the $G$ line intensity. We will first confirm that our model correctly predicts the strong dependence of $I_{G}$ on the laser frequency $\omega_{\mathrm{L}}$ and show that this strong dependence cannot be explained by the energy dependence of the joint density of states (JDOS) alone.

Secondly, we focus on the combined Fermi energy and $\omega_{\mathrm{L}}$ dependence of $I_{G}$. As measured and conceptually explained in Ref. Chen et al., 2011, the intensity of the $G$ line is greatly enhanced when $2\varepsilon_{\mathrm{F}}$ approaches a value of $\omega_{\mathrm{L}}-\omega_{\mathrm{ph}}/2$. Our first principles approach allows us to investigate this behavior as a function of $\omega_{\mathrm{L}}$ as well, going beyond the low-energy approximation used in the analytic work by Basko Basko (2009). We find that the relative increase of $I_{G}$ at the “critical” values of $\varepsilon_{\mathrm{F}}$ depends non-monotonically on the laser energy.

Finally, we demonstrate that the common textbook approximation Yu and Cardona (2010) of retaining only the almost double-resonant term in the Goldstone diagram expansion of the Raman amplitude fails in graphene. Instead, at low laser energies, the full set of leading order Goldstone diagrams must be included to get quantitative agreement with the full result. At higher energies, however, it is sufficient to approximate the full Raman amplitude with the sum of the almost double-resonant and the single-resonant terms. Moreover, we demonstrate that the simple picture of only taking into account the resonant states (within a band of $\pm$ the electronic decay width), fails. We find that the amount of states which contribute to the total $G$ peak intensity depends strongly on the laser energy. For low $\omega_{\mathrm{L}}$, we can confirm the prediction of the analytical work of Basko Basko (2009) that electronic states from almost the entire BZ contribute. For larger values of $\omega_{\mathrm{L}}$ (\omega_{\mathrm{L}}$$\gtrsim2.5 eV), however, we find that one needs to consider only a broad energy band of width $\approx$2 eV around the resonance energy window $[\omega_{\mathrm{L}}-\omega_{\mathrm{ph}},\omega_{\mathrm{L}}]$ to obtain quantitative agreement with the full Raman amplitude. By contrast, considering only the resonant energy band of width $\pm$100 meV severely overestimates the matrix element as destructive quantum interference effects are neglected.

III.1 Laser energy dependence

The calculated laser energy dependence of $I_{G}$ is shown in Fig. 2a (blue line). We find a strong dependence of $I_{G}$ on $\omega_{\mathrm{L}}$. This behavior cannot be simply understood in terms of the joint density of states (JDOS) (red line). The Raman intensity strongly peaks at a laser energy corresponding to the van Hove singularity of the JDOS at around 4.1 eV, connected to the flatness of the transition energy band structure (see inset of Fig. 2), as generally expected for optical spectra. At low excitation energies, however, the intensity is suppressed even though a sizable number of bright electronic transitions is available. This is in contrast to a simple absorption spectrum, which shows finite intensity whenever dipole-allowed electronic transitions are available. To understand this more complex behavior of $I_{G}$ as a function of laser energy, we divide the excitation energy range of interest into three different regimes, as indicated by the three colored shades in Fig. 2a.

In the first regime, corresponding to laser energies \omega_{\mathrm{L}}$$\lesssim1.5 eV, the Raman intensity is suppressed because of angular momentum conservation. In the low-energy limit, the band structure appears to be circularly symmetric (“Dirac cone” around the $K$ point, see inset of Fig. 2). Due to the continuous in-plane rotation symmetry in this regime, the $z$-component of the angular momentum is conserved and thus the initial and final states must carry the same total angular momentum. Since the electronic system is in its ground state both before and after the scattering event, its contribution to the total angular momentum can be ignored. A (circularly polarized) photon carries angular momentum of $\pm\hbar$ and so does the created (circularly polarized) $E_{2g}$ phonon, as it transforms like a two-dimensional vector under rotations. Therefore the final state, consisting of one photon and one phonon can possess total angular momentum $+2\hbar$, $0\hbar$, or $-2\hbar$. This differs from the angular momentum of the initial state (one photon) of $\pm\hbar$ and hence this process is disallowed by angular momentum conservation in the low-energy regime.

If one increases the laser energy to the second regime (1.5\lesssim$$\omega_{\mathrm{L}}$$\lesssim2.5 eV), the band structure loses its circularly symmetric shape and becomes trigonally warped. The continuous rotation symmetry is broken down to the 120*∘*-rotation symmetry of the lattice and exact angular momentum conservation is consequently reduced to angular momentum conservation up to integer multiples of $3\hbar$ only. Thus an initial state with angular momentum $\pm\hbar$ can be scattered to a final state with total angular momentum $\mp 2\hbar$. As a consequence, incoming light with polarization $\sigma^{\pm}$ has opposite polarization after undergoing one-phonon Raman scattering. The selection rules for linearly polarized light can be obtained from those for circularly polarized light by taking appropriate linear combinations and can be conveniently expressed in terms of Cartesian Raman tensors Loudon (1964):

[TABLE]

where the argument in parenthesis refers to the polarization of the $E_{2g}$ phonon created and the rows (columns) refer to the $x$- or $y$- components of the incoming (outgoing) light polarization. We verified that our calculations respect these selections rules up to a relative factor of the order of $10^{-6}$.

To understand that the Raman intensity is still relatively small in the second regime, we focus on the contribution of the different $\mathbf{k}$-points to the total Raman matrix element $\mathcal{M}_{\mathbf{k}}=\sum_{\mathbf{k}}\mathcal{M}_{\mathbf{k}}$. To this end we look at the absolute value and phase of $\mathcal{M}_{\mathbf{k}}$ on the high-symmetry line $\Gamma$-$K$-$M$-$\Gamma$ for three different values of $\omega_{\mathrm{L}}$, as shown in Fig. 3.

The absolute value of $\mathcal{M}_{\mathbf{k}}$ is given on the vertical axis and its relative phase is color-encoded. The top horizontal axis displays the corresponding electronic transition energy $\Delta\varepsilon_{\mathbf{k}}$, while the shaded background shows the value of the JDOS at this transition energy. The latter is a measure for the weight of the corresponding $\mathbf{k}$ point in the total matrix element when a full 2D integration over the BZ is performed.

One can clearly identify the resonant states from the figure, which correspond to the largest values of $|\mathcal{M}_{\mathbf{k}}|$ (compare zoom-in shown in lower panels). For both $\omega_{\mathrm{L}}$=1.5 and 2.5 eV, the resonant states are sharply centered around one point on each side of $K$. By looking at the phase of the resonant $\mathbf{k}$-points, one notices that for each “resonance peak” the phase continuously undergoes a change of $\pi$ when passing over the peak, i.e., contributions from $\mathbf{k}$-points to both sides of the resonant $\mathbf{k}$-point have opposite phase and cancel each other (see blue to yellow dots in the first two panels of Fig. 3). This change of phase by $\pi$ is typical for a driven system and is well-known from the simple system of a driven and damped harmonic oscillator. It is of particular importance for the case of Raman scattering however, as the total Raman amplitude is a sum over contributions from different $\mathbf{k}$-points and a relative phase of $\pi$ leads to destructive interference, as seen for the states centered around the resonant $\mathbf{k}$-point. While the contributions from the vicinity of the resonant $\mathbf{k}$-points to both sides of $K$ cancel separately, it should also be noted that the “inner” flanks of the two resonance peaks in the $\Gamma$-$K$ and $K$-$M$ directions have opposite phase and also mostly cancel each other as they have similar amplitude and weight. This is related to the residual continuous rotation symmetry. A final point to note is the contribution of the non-resonant states near the van Hove singularity in the vicinity of the $M$ point. While the amplitude of these contributions is two orders of magnitude smaller than those of the resonant $\mathbf{k}$-points, the corresponding $\mathbf{k}$-points are broadly spread along the $K$-$M$ direction of the high-symmetry line, i.e., there are a lot of states along the high-symmetry lines with similar amplitude, which are also all in phase (red dots around $M$ in the first two panels of Fig. 3). Furthermore, when performing a 2D-integration over the entire first BZ, these states enter with a lot of weight as the JDOS peaks at the corresponding transition energies. Thus, the total contribution of the region around the van Hove singularity at $M$, while non-resonant, is still far from negligible. However, it turns out that this contribution is mostly canceled by the sum of the non-resonant contributions from the rest of the BZ (cyan dots in the $\Gamma$-$K$ an $M$-$\Gamma$ directions in the first two panels of Fig. 3), which have opposite phase. This fact will be demonstrated also (and more conclusively) in the last section of this paper.

This picture of dominating destructive quantum interference effects changes significantly for higher laser energy (\omega_{\mathrm{L}}$$\gtrsim2.5 eV, see third shaded region in Fig. 2a). As seen in the third panel of Fig. 3, the contributions of the flanks of the two resonance peaks still destructively interfere with each other (blue to yellow dots). However, the resonance peak in the $K$-$M$ direction becomes very broad along the high-symmetry line and the contributions are mostly in-phase (cyan to green dots). Combined with the high JDOS at the corresponding electronic transition energies, the contributions from these resonant states dominate the total Raman matrix element $\mathcal{M}$ for higher laser energies. This dominance is ultimately a consequence of the flatness of the transition band structure near the van Hove singularity at the $M$ point (see inset of Fig. 2). Compared to the lower laser energy regime, here, the contributions from around the van Hove singularity cannot be canceled by the sum of the non-resonant contributions from the bulk of the BZ, as the relative weight and amplitude of the resonant states is much too large compared to those of the non-resonant states.

We can summarize the discussion of the laser energy dependence of the $G$ peak intensity as follows: In the very low energy regime, the intensity is suppressed due to angular momentum conservation associated with the continuous rotation symmetry in the low-energy regime. As the excitation energy increases, the intensity still remains low due to destructive quantum interference effects, which leads to separate cancellations among both the resonant and non-resonant contributions to the Raman matrix element. Finally, for larger laser energies (\omega_{\mathrm{L}}$$\gtrsim2.5 eV), the resonant, in-phase contributions from the $K$-$M$ direction dominate the Raman matrix element. Due to their increased weight as one approaches the van Hove singularity at $M$, they cannot be suppressed by the sum of the non-resonant contributions.

III.2 Fermi energy dependence

Next, we study the Fermi energy dependence of the $G$ peak intensity. We varied the Fermi level $\varepsilon_{\mathrm{F}}$ from -3 to +3 eV relative to that of pristine graphene and calculated the resulting $G$ peak intensity for different values of the laser energy $\omega_{\mathrm{L}}$. We use the rigid band approximation, i.e., we do not take into account the change of the electronic bands with $\varepsilon_{\mathrm{F}}$. The results of the calculation, normalized to the respective intensity for pristine graphene, are shown in Fig. 4a.

The different values of $\omega_{\mathrm{L}}$ are represented by color, from $\omega_{\mathrm{L}}$=1.5 eV (red) to $\omega_{\mathrm{L}}$=4 eV (violet). For both electron (right half) and hole doping (left half of the plot), the $G$ peak intensity shows a strong increase when $\varepsilon_{\mathrm{F}}$ approaches a critical value. This critical value corresponds well to the conduction or valence band energy of those states in the $K$-$M$ direction that are in resonance with the average of the incoming and outgoing light energy, $(\omega_{\mathrm{in}}$+$\omega_{\mathrm{out}})/2$, as shown in Fig. 4b. The strong behavior of $I_{G}$ as a function of $\varepsilon_{\mathrm{F}}$ has been observed in experiment Chen et al. (2011), but so far only for one fixed laser energy. As can be seen from Fig. 4a, the strength of the relative increase of $I_{G}$ is expected to strongly vary with laser energy. Furthermore, the electron-hole asymmetry (full black line in Fig. 4b) of the $\varepsilon_{\mathrm{F}}$ values at which $I_{G}$ peaks increases with excitation energy $\omega_{\mathrm{L}}$, which can be understood from the electron-hole asymmetry of the band structure at the resonant $\mathbf{k}$-points in the $K$-$M$ direction (dashed black line in Fig. 4b).

A conceptual explanation for the observed sensible dependence of $I_{G}$ on $\varepsilon_{\mathrm{F}}$ has already been suggested in the original paper by Chen et al. Chen et al. (2011). The authors attributed the strong increase of $I_{G}$ at the critical values of $\varepsilon_{\mathrm{F}}$ to blocking of destructive quantum interference due to the Pauli principle. When the Fermi level is increased (lowered), increasingly more transitions from the $\pi$ to the $\pi^{*}$ band are blocked, as the corresponding electronic states in the conduction (valence) band become occupied (unoccupied). As a result, these transitions become blocked by the Pauli principle. In Fig. 4c, this concept is illustrated for the case of electron doping and for three different laser energies. The shaded regions mark those $\mathbf{k}$ points that do not contribute to the total Raman matrix element as the corresponding electronic transitions are Pauli blocked. In particular, if the Fermi level is tuned to certain values, a big part of destructive quantum interference effects is switched off and as a result the $G$ peak intensity increases. The strength of this increase, however, strongly depends on the laser energy. To understand this last point, we consider the two extreme values of $\omega_{\mathrm{L}}$ considered in our calculation, $\omega_{\mathrm{L}}$=1.5 and 4 eV, in more detail.

For the case of $\omega_{\mathrm{L}}$=1.5 eV, we see from Fig. 4a that, on either side of charge neutrality point, the intensity of the $G$ peak first goes through a strong increase, followed by a minimum, and finally a smaller second peak that trails out as $\varepsilon_{\mathrm{F}}$ is increased even more. This sequence of peaks and dips is easily understandable in terms of the picture of quantum interference discussed for the “mid-level” energy regime in the previous section. There we demonstrated the fact that in this $\omega_{\mathrm{L}}$-regime, the total Raman matrix element is suppressed because, on the one hand, almost-resonant contributions of states from $\mathbf{k}$-points around the resonance points mostly cancel each other, while, on the other hand, non-resonant contributions from states around the van Hove singularity are mostly canceled by the sum of non-resonant contributions from the bulk of the BZ (see also the left-most panels of Fig. 3 for illustration). When the Fermi level is now increased from the charge neutrality point at $K$, more and more $\mathbf{k}$-points around $K$ do not contribute anymore to the total matrix element because of the Pauli principle. In terms of the “resonance peak” picture in the $\mathcal{M}_{\mathbf{k}}$ vs. $\mathbf{k}$-plot from Fig. 3, this means that at first, an increase in $\varepsilon_{\mathrm{F}}$ blocks contributions from the “inner” flanks of the two resonance peaks, which destructively interfere with the “outer” flanks. This results in the observed strong increase of $I_{G}$. The increase of $I_{G}$ continues until the Fermi level reaches the point where it blocks the resonant transitions, at which point $I_{G}$ reaches a peak value (at approximately $\varepsilon_{\mathrm{F}}$=0.71 eV for electron doping). An even further increase of $\varepsilon_{\mathrm{F}}$ results in a blocking of contributions from the constructively interfering “outer” flanks and thus to a decrease of $I_{G}$ until a minimum intensity is reached (compare dip in red curve of Fig. 4a at around $\varepsilon_{\mathrm{F}}$=1.55 eV), corresponding to all almost-resonant transitions being Pauli blocked. At this value of $\varepsilon_{\mathrm{F}}$, the only contributions to the total Raman matrix element come from the non-resonant states from around the van Hove singularity and from the bulk of the BZ, which cancel each other. A further increase of $\varepsilon_{\mathrm{F}}$, however, leads to Pauli blocking of the van Hove singularity contribution, which can now no longer destructively interfere with the non-resonant contributions from the bulk of the BZ. This leads to another increase of $I_{G}$ again, which peaks for a second time when $\varepsilon_{\mathrm{F}}$ reaches the van Hove singularity at $M$ (see small bump in red curve of Fig. 4a at around 1.97 eV). An even further increase of $\varepsilon_{\mathrm{F}}$ only leads to more and more parts of the bulk of the BZ being Pauli blocked and consequently $I_{G}$ trails out with increasing $\varepsilon_{\mathrm{F}}$.

By contrast, in the case of $\omega_{\mathrm{L}}$=4 eV the $G$ peak intensity shows a much simpler and less strong behavior as a function of $\varepsilon_{\mathrm{F}}$. When $\varepsilon_{\mathrm{F}}$ is increased away from the charge neutrality point, $I_{G}$ shows almost no response at first before lightly peaking at a value of approx. $\varepsilon_{\mathrm{F}}$=1.65 eV. After that, $I_{G}$ quickly reduces to zero when $\varepsilon_{\mathrm{F}}$ is increased even more. This is again consistent with our earlier observation that for higher laser energies, the $G$ peak is mostly carried by resonant transitions and quantum interference effects only play a minor role. When $\varepsilon_{\mathrm{F}}$ is first increased, $I_{G}$ remains largely insensitive as only non-resonant transitions from around $K$ are Pauli blocked. Only when $\varepsilon_{\mathrm{F}}$ approaches the resonant states near the van Hove singularity at $M$ does $I_{G}$ feature a response, when the destructive quantum interference effects from the “inner” flank (see right-most panel in Fig. 3) are switched off. However, since the “resonance peak” in $\mathbf{k}$-space is broad and mostly in-phase, the overall response of the $G$ peak intensity is not as strong for $\omega_{\mathrm{L}}$=4 eV as for lower laser energy, where the resonance peaks in $\mathbf{k}$ space are very narrow and quantum interference between the opposing flanks plays a much bigger role. After $I_{G}$ for $\omega_{\mathrm{L}}$=4 eV reached a peak value at around $\varepsilon_{\mathrm{F}}$=1.65 eV, any further increase of the Fermi level only leads to the blocking of more and more resonant states, leading to a sharp decrease of $I_{G}$ (compare violet curve in Fig. 4a). After the entire broad “resonance peak” in $\mathbf{k}$-space has been Pauli blocked, the relative value of $I_{G}$ compared to the zero doping case remains insensitive to any further change of $\varepsilon_{\mathrm{F}}$, confirming our previous finding that non-resonant states from the bulk of the BZ only play a minor role at higher laser energies.

To summarize the Fermi energy dependence of the $G$ peak intensity, our ab initio calculations confirm the suggestion of Chen et al. that the blocking of quantum interference effects by shifting the Fermi level is the driving mechanism behind the strong increase of $I_{G}$ as $\varepsilon_{\mathrm{F}}$ reaches a critical value. Going beyond this, we demonstrated that the relative increase of $I_{G}$ at these critical Fermi levels strongly depends on $\omega_{\mathrm{L}}$, as destructive quantum interference effects play an increasingly less dominant role with increasing $\omega_{\mathrm{L}}$. We also predict that for small laser energies, the $G$ peak intensity shows a small resurgence even after all resonant states have already been blocked, when $\varepsilon_{\mathrm{F}}$ approaches the van Hove singularity, thus blocking its destructive influence on the non-resonant contributions from the bulk of the BZ.

III.3 Relevant states for the $G$ peak process

Finally, we address a question first discussed for the low-energy case by Basko Basko (2009), namely which electronic transitions are relevant for the $G$ peak, i.e., which electronic states need to be taken into account in a theoretical description. In the preceding two sections we demonstrated that the influence of quantum interference effects and of non-resonant states strongly depends on the laser energy. It is thus not surprising that the states that need to be considered for a quantitative description of the $G$ peak also vary with $\omega_{\mathrm{L}}$.

To obtain a clear picture of which states contribute to the $G$ peak, we calculated the $G$ peak intensity for different $\omega_{\mathrm{L}}$ with increasingly more states, starting with the ones that are in resonance with the incoming and/or outgoing light. To be more precise, we introduce a transition energy window width $\varepsilon_{\mathrm{cut}}$ and include only those electronic transitions that obey the condition

[TABLE]

i.e., the included states lie within a window of width $2\varepsilon_{\mathrm{cut}}$ around the average resonance frequency $(\omega_{\mathrm{in}}$+$\omega_{\mathrm{out}})/2$, as illustrated in the inset of Fig. 5a.

We then proceed to calculate $I_{G}$ with only the states inside the selected window. The resulting function $I_{G}(\varepsilon_{\mathrm{cut}})$ is plotted for different laser energies $\omega_{\mathrm{L}}$ in Fig. 5a. For each $\omega_{\mathrm{L}}$, we normalized with respect to the result obtained with all states included, so that the all curves tend to one for large $\varepsilon_{\mathrm{cut}}$.

In Fig. 5a, we observe the same behavior as discussed in the previous two sections. When $\varepsilon_{\mathrm{cut}}$ is first increased from zero, the first states to be included are the states that are in resonance with the ingoing and/or outgoing light, which contribute with a large amplitude to the Raman matrix element. As a result, $I_{G}(\varepsilon_{\mathrm{cut}})$ displays a strong increase when the resonant states are first included. The next states to be included stem from the flanks of the “resonance peaks” in $\mathbf{k}$-space (compare once more Fig. 3). These states differ in phase from the resonant states and thus partially destructively interfere, which leads to a decrease of $I_{G}(\varepsilon_{\mathrm{cut}})$ when these states are included.

For low values of the laser energy $\omega_{\mathrm{L}}$ (red to green curves in Fig. 5a), $I_{G}(\varepsilon_{\mathrm{cut}})$ increases once more when $\varepsilon_{\mathrm{cut}}$ is increased. This second peak corresponds to the inclusion of the states from around the van Hove singularity, which are all in-phase and thus constructively interfering. Indeed, the point of steepest ascend in the red to green curves of Fig. 5a corresponds precisely to those values of $\varepsilon_{\mathrm{cut}}$ at which the energy window encompasses the $M$ point. When $\varepsilon_{\mathrm{cut}}$ is increased even more, the destructively interfering states from the bulk of the BZ are taken into account and hence $I_{G}(\varepsilon_{\mathrm{cut}})$ decreases again. When enough of these non-resonant states are taken into account, $I_{G}(\varepsilon_{\mathrm{cut}})$ finally converges to the value of $I_{G}$ calculated with all states.

For higher values of the laser energy though, $I_{G}(\varepsilon_{\mathrm{cut}})$ does not change any more after the resonant states have been taken into account, as seen in the behavior of the blue to violet lines in Fig. 5a, which converge straight to the correct value of $I_{G}$. This once more shows that destructive interference from non-resonant states does not play a significant role in this $\omega_{\mathrm{L}}$-regime.

To further prove this last point, we can consider the resonant and non-resonant terms separately and calculate the $G$ peak intensity only with some of them. Recall from Eq. 5 that the total matrix element can be written as sum of an almost double-resonant (aDR), a single-resonant (SR), and a non-resonant contribution (NR), with each one of them already implicitly summed over $\mathbf{k}$-points. If we calculate $I_{G}$ with only the aDR-, the SR-, or the NR-contribution as a function of $\omega_{\mathrm{L}}$ and normalize to the value of $I_{G}$ including all contributions, we obtain the dashed lines shown in Fig. 5b. Doing the calculation only with the aDR- or the NR-terms severely overestimates the $G$ peak intensity for low excitation energies, while it underestimates the correct result for $I_{G}$ for larger $\omega_{\mathrm{L}}$. A calculation with only the single-resonant terms, on the other hand, always underestimates the correct value for $I_{G}$.

More interesting, however, is the result we obtain when we consider the sum of the aDR- and SR-contributions in the total matrix element, i.e., we include all terms in the total matrix element that can go resonant with either the incoming or outgoing light. The resulting curve is shown as a green line in Fig. 5b. The plot confirms our previous analysis in that, for higher values of the laser energy (\omega_{\mathrm{L}}$$\gtrsim2.5 eV), the resonant contributions are the dominating ones. Indeed, a calculation that only takes into account the resonant terms in the Raman matrix element yields very good results in this laser energy regime.

As a last point, we can finally visualize the required states for different values of $\omega_{\mathrm{L}}$ by reconsidering our result obtained with a transition energy window width $\varepsilon_{\mathrm{cut}}$ shown in Fig. 5a. If we start from the maximum considered value of the cutoff (the full $\pi$-band width), and gradually lower the value of $\varepsilon_{\mathrm{cut}}$ until $I_{G}(\varepsilon_{\mathrm{cut}})$ differs from the correct value of $I_{G}$ by more than 2%, we obtain a minimum transition energy window width around the resonant states that are needed to achieve 2% accuracy. For this minimal value of $\varepsilon_{\mathrm{cut}}$, we can visualize the corresponding $\mathbf{k}$-points within this transition energy window in a plot of the band structure. This is shown for the two laser energies $\omega_{\mathrm{L}}$=1.5 and $\omega_{\mathrm{L}}$=4 eV in the inset of Fig. 5b. In this representation it becomes immediately clear, that for low laser energies the relevant electronic transitions for the $G$ peak come from a large part of the BZ, in agreement with the findings of Basko Basko (2009) in his analytical work in the low-energy limit. For larger values of $\omega_{\mathrm{L}}$, this picture changes, however, in that the relevant states are localized in a broad band of total width $\sim$2 eV around the resonance energy. It is important to note that in neither $\omega_{\mathrm{L}}$-regime is it sensible to approximate the Raman amplitude with an expression that only takes into account resonant electronic transitions (in an energy window of $\pm$ the inverse electronic life time).

IV Conclusions

In conclusion, we developed a fully ab initio approach to the calculation and study of one-phonon Raman intensities and applied it to the case of monolayer graphene. We demonstrated explicitly that quantum interference effects play a crucial role for the understanding and interpretation of the Raman process leading to the $G$ peak of graphene. In particular, quantum interference effects from non-resonant states are very important for a correct low-excitation energy description of the $G$ peak. We confirmed the conceptual picture of Chen et al. Chen et al. (2011) that Pauli blocking of destructive quantum interference by tuning the Fermi level to a critical value can lead to a large enhancement of the $G$ peak intensity. Furthermore, we predict that the relative increase of the $G$ peak intensity at the critical Fermi level strongly depends on the laser energy. Finally, our general approach enabled us to explicitly demonstrate which electronic states need to be taken into account for an accurate description of the $G$ peak for a wide range of laser energies, going beyond the previous analytical analysis by Basko Basko (2009) in the low-energy regime. In particular, we showed that for low laser energy, transitions from almost the entire first Brillouin zone need to be considered, whereas at higher laser energy, only states in a broad band of full width $\sim$2 eV around the resonance energy need to be taken into account.

While in this work we focused on a study of one-phonon Raman intensities in monolayer graphene, our implementation and method of analysis is completely general. As such, it can be applied to any material and enables the study of other materials as well.

V Acknowledgments

The authors would like to thank H.P.C. Miranda for helpful discussions. S.R. and L.W. acknowledge financial support by the National Research Fund (FNR) Luxembourg (projects RAMGRASEA and TMD-nano).