Impurity induced broadening of Drude peak in strained graphene

The recent experimental study of the far-infrared transmission spectroscopy of monolayer graphene has shown that the Drude peak width increases by more than $10\%$ per $1\%$ of applied mechanical strain, while the Drude weight remains unchanged. We study the influence of the strain on the resonant impurity scattering. We propose a mechanism of augmentation of the scattering rate due to the shift in the position of resonance. Using the Lifshitz model of substitutional impurities, we investigate changes in the Drude peak weight and width as functions of the Fermi energy, impurity concentration and magnitude of the impurity potential.


Introduction
This work is dedicated to 80th birthday of our colleague and a remarkable physicist Professor Ihor Stasyuk. It is a pleasure for us to present a recent result on strained graphene. We study intricate physics related to the interplay between the band structure and strain using the methods of the Green's function approach, to which I. Stasyuk has made a significant contribution [1].
Outstanding mechanical strength and stiffness of graphene are marked by its capability to sustain reversible elastic deformations in excess of 20% [2]. These features are related to in-plane σ bonds formed by hybridized 2s, 2p x and 2p y electron orbitals of carbon atoms in the lattice. The remaining 2p z orbitals form valence and conductance bands responsible for the observed transport and optical properties of graphene (for a review, see, e.g. [3]).
In the simplest case, the strain is uniaxial. This strain configuration has been investigated by quite a few authors. In particular, a number of theoretical studies [4][5][6][7][8] predicted that the optical conductivity changes anisotropically with respect to the direction of the strain. These predictions seem to be in agreement with the transparency measurements in the visible range made on large-area chemical vapor deposited (CVD) monolayer graphene pre-strained on a polyethylene terephthalate (PET) [9].
The recent study in optical conductivity response to the strain [10] revealed some unexpected features in the far-infrared transmission spectra of a CVD monolayer graphene on PET substrate. It was found that the Drude peak width increases by more than 10% per 1% of the applied uniaxial strain, while the Drude weight remains unchanged.
The purpose of the present work is to study the impact of strain on the optical conductivity in the presence of point defects. These defects are assumed to be either chemically substituted carbon atoms including their absence, i.e., vacancies, or adsorbed atoms or molecules on a graphene sheet. They may originate as a byproduct of a fabrication method, or can be purposefully superimposed on the graphene sheet by exposure to an active environment. In either case, one should expect a finite amount of defects on graphene [11].
In [10], contribution of the point defects to the Drude peak width was considered within the Born approximation. It was found that the resulting variation of the Drude width under the strain is insignificant. It is known, however, that certain point defects may lead to the appearance of the resonance impurity states. The resonance energy can be located near the Dirac point, in which case the resonance is well defined. This makes properties of the system, such as optical conductivity, very sensitive to the position of the Fermi energy E F relative to the resonance energy. Even a small variation of this parameter can lead to a dramatic increase in the Drude width.
The paper is organized as follows. We begin by presenting in section 2 the Hamiltonian of the electronic subsystem in strained graphene. To describe impurities, we employ the Lifshitz model. In section 3, we use the Green's function formalism to calculate the electron self-energy function within the average T-matrix approximation. Then, we use the self-energy function to calculate the Drude weight and width. In section 4, we present the results for various impurity concentrations and impurity perturbations, and analyze conditions for which the used approximations are valid.
Throughout the paper, units = e = 1 are chosen. Particularly, we will omit e 2 / 2 factor in conductivity and the Drude weight (the last one is thus measured in units of energy).

Model
It is assumed that electrons in impure graphene can be described by the Hamiltonian of separable formĤ with the first term characterizing the electronic subsystem of a clean graphene sheet, and the second term adding impurity perturbation to the model. The crystal lattice for an undeformed sample with impurities is depicted in figure 1. We assume that the substitutional impurities do not change the lattice structure. For a clean graphene crystal, the hexagonal structure has two non-equivalent sublattice positions A and B for each Bravais lattice cell. We use the second quantization formalism, with the creation (annihilation) operatorsĉ + αi (ĉ αi ) which add (remove) the electronic state at the i-th cell on the α = A, B sublattice. Usually, one works in the nearest-neighbour approximation, where the only non-zero hopping integrals are those between the nearest neighbours A

33703-2
and B. Therefore, the resulting form of the Hamiltonian for the hexagonal lattice is: where the sum i j goes over the closest Bravais cells i and j. Here, t i j are the hopping integrals for the nearest-neighbouring atoms with indices (i, A) and ( j, B). In the case of unstrained graphene, the hopping integrals are position and direction independent, t i j = t. The value of the hopping integral is usually estimated as t ≈ 2.7 eV [12]. As we turn on the strain, the distances between the atoms vary, and in general case, t i j no longer takes such a simple form.
Since we consider a simple uniaxial strain, which is uniform along the direction of the applied stress, the components of two-dimensional strain tensorε ε ε are independent of equilibrium positions of atoms r. Accordingly, the displacement vector u(r) reads u(r) =ε ε ε · r. Thus, the actual position of an atom r = r + u(r) can be written as r = (Ī +ε ε ε) · r, whereĪ is the unit 2 × 2 matrix.
The strain tensor acquires the diagonal formε ε ε = diag(ε, −νε) in the coordinate system aligned with the axis of the applied stress. Here, ε is the relative longitudinal expansion, and ν is the Poisson's ratio (ν = 0.16 in graphite and may possibly be even smaller in graphene on substrate [10]). Note that the planar deformation of the hexagonal crystal in the basal plane is determined by the two independent stiffness (compliance) tensor components, viz. it behaves as an isotropic planar solid [13].
The three vectors that connect atom A with the nearest neighbouring atoms B, δ n (n = 1, 2, 3), change from their equilibrium value δ (0) n (see figure 1) independently of the cell index: δ n = (Ī +ε ε ε) · δ (0) n . Accordingly, the hopping integral does not depend on the cell index i of the A atom, but only on the direction to the neighbours. Thus, we end up with the three distinct hopping integrals, which we will denote as t {n} . It is usually assumed that the change in the overlap integrals of p z -orbitals [14] can be adequately described by the first order expansion in strain (see also in [15]): where a ≈ 1.42 Å is the equilibrium distance between the carbon atoms, and β = −d ln t/d ln |δ| ∼ 3−4 is the dimensionless Gruneisen parameter. It can be shown that under uniform uniaxial strain, the density of states (DOS) per unit cell and spin (including the valley degeneracy) acquires the form [8,16,17] where W ε is the strain dependent bandwidth, Here, W 0 = ( √ 3π) 1/2 t 2.33t is the effective bandwidth or the energy cutoff for unstrained graphene that preserves the number of states in the Brillouin zone. It is worth mentioning that energy is counted from the Dirac point. In what follows, we will assume that ν = 0, because if necessary it can be easily restored by replacing ε → ε(1 − ν).
In the simplest model of substitutional impurities, widely referred to as the Lifshitz model [18], we add the constant potential V L to the sites occupied by the impurity atoms. This leads to the diagonal term in the HamiltonianV where η αl equals one on the impurity sites, and zero otherwise. For a lattice with N atoms with the impurity concentration c, the total number of impurity sites is cN, but the distribution of impurities is not specific, and differs from one sample to another. This model can also be used to describe vacancies in a graphene sheet deposited on a substrate. Additionally, it can describe adsorbed atoms, provided that the hopping integral between the additional orbital of the adsorbed atom or molecule, and the p z orbital of the carbon atom of the host is larger than the bandwidth W 0 .

Formalism
The problem of inclusion of impurities in the electronic structure can be treated by using the Green's function method (see [19][20][21] for a general review, or [22] for a graphene-specific application). The Green's function is related to the Hamiltonian (1) as (Î is the unit matrix): We prefer to work in the site representation, where both the Hamiltonian and the Green's function are N × N matrices in the basis of one-particle statesĉ + αl |0 (here |0 is the vacuum state). In a general case, this Green's function matrix has a complicated structure due to irregularity in positions of the impurities. To treat this problem, we are averaging the Green's function over all possible configurations of impurities at a fixed concentration c:Ĝ = Ĝ . It turns out that with an increasing number of atoms in the lattice, the Green's functionĜ approaches the average valueĜ [23]. The averaging reestablishes the translational invariance, and the averaged Green's functionĜ can be related to the host Green's functionĝ by means of the Dyson equation:Ĝ =ĝ +ĝΣĜ.
If we take into consideration only the single-site scattering, the self-energy operatorΣ becomes diagonal.
Omission of the cluster scattering leaves us with the scalar function Σ = Σ(E), expressed viaΣ = ΣÎ. Therefore, the Green's functionĜ can be expressed viaĝ aŝ In the present work we employ the average T-matrix approximation (ATA), in which Σ(E) is expressed via the diagonal element of the host Green's function in the site representation g 0 (E) = g αi,αi (E) as follows: To obtain an analytical expression for Σ to work with, we derive g 0 via the relation between it and the density of states per unit cell ρ 0 (E): This determines the imaginary part; the real part is restored using the Kramers-Kronig relation. The resulting expression for g 0 is: We now discuss how the the calculated self-energy Σ (ATA) is related to the spectroscopy measurements [10]. One of the advantages of the infrared spectroscopy as compared to the DC transport measurements is that it allows one to find independently both the Drude spectral weight and the optical scattering rate. We assume that the real part of the Drude conductivity is of a Lorentzian shape where D is the Drude weight and Γ is the Drude peak width or optical scattering rate. Both these quantities depend on the Fermi energy. The Drude weight can be related to the dc limit, σ dc = Re σ(ω = 0) = D/(πΓ), of the dynamical conductivity (14).
According to the Matthiessen's rule, the total optical scattering rate is Γ = 2 i Γ i , where Γ i are the contributions from the different channels of single particle scattering, e.g., short-range point defects, long-range charged impurities, acoustic phonons, surface phonons in the substrate and grain boundaries.

33703-4
Each of these contributions may be strain dependent. As was mentioned above, in the present work we restrict ourselves to the scattering by point impurities: Although the uniaxial strain makes the Drude weight anisotropic [10], we will only use the Drude weight averaged over strain directions. Thus, using the expression for the dc conductivity for graphene obtained in the bare bubble approximation at zero temperature [24], we obtain For small concentrations of impurities, when inequality |Σ(E)| E holds, equation (16) takes a rather simple form

Strain dependence of the impurity resonance energy
Since a small impurity concentration, c 1, is considered, it is safe to neglect c in the denominator of the ATA expression (11). Substituting the diagonal element of the host Green's function (12) in equation (11), we obtain, for example, that the imaginary part of the self-energy reads The denominator of (18) has two terms, both of which are non-negative. While the last term is a smooth function of energy, the first term can be zero. Therefore, from the last equation it is clear that the impurity scattering rate (15) should have a maximum at some energy. For |V L | W ε the location of this maximum E r can be approximated by the solution of the well-known Lifshitz equation, 1 = V L Re g 0 (E r ), or For V L < 0, the position of the extremum is above the Dirac point, E r > 0, and vice versa. This unusual positioning property holds true for a general spectrum consisting of two symmetric bands touching each other at a single point [25]. Without compromising generality of our consideration, we will consider the case V L < 0, and the electron-doped graphene E F > 0.
To find the solution of equation (19), one can use a recursive procedure: , with E (n) r , n → ∞ converging to the exact solution E r . Alternatively, one can express E r in terms of the inverse function W(z) defined as a solution to the equation z = W exp(W), also known as the Lambert W-function. This function has two branches W 0 (z) and W −1 (z), which represent separate solutions of the equation. Specifically, we will use the W −1 branch, which gives |W −1 (x)| > 1 for x < 0. In terms of this function,

33703-5
ReΣ(E): This extremum at E r can be regarded as a well-defined impurity resonance when its width is much smaller than the difference between E r and the Dirac point energy. The corresponding condition | Im Σ(E r )| |E r | gives us the following According to (19), (20), (21), small |E r | corresponds to large |V L |. Thus, a well-defined resonance in the Lifshitz model is present only for |V L | W ε . It should be understood that the Lifshitz model is, in a certain sense, phenomenological. More complicated impurity models that are more appropriate to the real physical situation solve the apparent problem of the excessively large impurity potential. At the same time, these models yield qualitatively similar results. We would like to stress that the Lifshitz impurity model is capable of providing a semi-quantitative description of experimental ARPES data on the impure graphene [26].
We show in figure 2 both Im Σ(E) and Re Σ(E) computed using equation (11) in the absence of strain and for ε = 5%. A moderate value V L = −6t is chosen and the concentration of impurities is c = 10 −3 . Both the resonance dip in Im Σ(E) and the corresponding S-shaped feature in Re Σ(E) are evident in figure 2. These features are expressed in a narrow energy window. The span of the window corresponds to the above obtained width of the impurity resonance ζ E r . One can also see in figure 2 that a small strain ε slightly shifts the position of the resonance towards the Dirac point. This shift, as we will show later, results in a considerable increase in the scattering rate Γ.
Let us assume that the position of the resonance E r changes under the strain linearly: where E r (0) is the resonance energy in the absence of strain. Substituting equation (23) to the Lifshitz equation (19), one can analytically estimate the value of the coefficient α: In figure 3 (a), we show both the value of E r (dash-dotted black line) and the energy of the minimum E m of the function Im Σ(E) (solid green line) as functions of strain ε. The value V L = −12t is used. One finds that the dependence E m (ε) is significantly off the line E r (ε) discussed above. However, the function E m (ε) has a linear dependence similar to (23), but with a slightly different slope α . To obtain α for the corresponding slope, one should substitute λ in equation (24) is the value of the energy that gives the minimum at ε = 0. This difference between E r and E m is connected to the fact that the condition for ζ (22) is not decisively satisfied. Figure 3 (b) shows the resonance in Γ as a function of the Fermi energy E F . The peak position is independent of the impurity concentration c. This allows one to tune the strength of the potential V L to match the actual position of the resonance energy E r in the experimental sample.
The imaginary part Im Σ(E) has an additional maximum at E = 0, which corresponds to the Dirac point in the clean sample. In this point, equation (18) yields Im Σ(0) = 0 independent of strain. As it turns out, in the interval of the order |cV L | around the Dirac point, our approximation does not hold. An accurate description of the conduction band boundaries should include the renormalization of the propagator.

Small impurity concentrations
First, we analyze the case of those concentrations c for which |Σ(E F )| E F . This condition is satisfied for c much smaller than some characteristic value c χ . The actual value of c χ depends on V L ; the larger |V L | is, the smaller concentrations are required for deviations from the following analysis to show up. This concentration can be roughly estimated as c χ ∼ E 2 r /W 2 0 (see the sub-section 4.3 for a detailed discussion).
With the strain ε, the Drude peak width Γ(E F ) (15) increases for E F residing below the resonance energy. The relative increase of the Drude width can be quite substantial, especially for values of the Fermi energy that lies at the half-width of the peak. Figure 3 (b) evidently demonstrates this feature. The gain turns out to be linear in strain, and can be described by the increment , which gives a relative increase of the scattering rate Γ per one percent of the strain ε: We can easily estimate the maximum gain m for a given impurity perturbation V L . First, let us assume that E F remains constant throughout the process of deformation. The results for this scenario are shown in figure 4 (a). As can be seen in figure 3 (b), the maximum values of (E F ) are reached when the Fermi   energy lies at the half-width of the resonance. This happens when the first and the second term in the denominator of (18) are equal. Assuming that V L can be chosen arbitrarily, the maximum gain for a given E F can be estimated analytically as This can be seen in figure 4 (a), where the curves for various V L < 0 are enveloped by the dashed curve described by (26). It turns out that the relative increase in Γ is independent of impurity concentration. Apart from observing the substantial increase of the Drude width for E F E r , we see that for E F E r one should expect a decrease in Γ as we increase the strain ε. Assuming that the resonance peak is nearlysymmetrical, we conclude that the same estimate as in (26), but taken with a negative sign, applies. The result is shown in figure 4 (a) by the lower dashed enveloping curve.
However, this result should be treated with care. Apart from simplicity of the impurity model, it was assumed that the value of the Fermi energy remains constant in the strained sample. As follows from equations (4), (5), the density of states ρ ε (E) changes as we stretch the graphene sheet. In an isolated sample, the number of charge carriers remains constant. Accordingly, the Fermi energy E ε F becomes strain dependent. Neglecting the effects of impurities, the concentration of charge carriers per unit cell reads N c = (E ε F /W ε ) 2 . This gives us For a positive strain ε, we expect the Fermi energy to reduce to approximately 98% of its initial value per 1% of strain. The Drude weight D ≈ E ε F changes accordingly. The results for constant carrier density are shown in figure 4 (b). We present a plot in terms of W 0 √ N c , which corresponds to the zero-strain value of the Fermi energy. The inclusion of the E ε F shift damps the gain significantly.
Analogously to (26), we can estimate the maximum gain by including both the effect of the narrowing of the bandwidth (5) [already included in equation (26)] and the reduction of the Fermi energy (27). We end up with the following analytical expression: We can see in figure 4 (b) that this estimate is approached closely in precise calculation, both for upper and lower enveloping curves. 33703-8

Moderate concentrations
As the impurity concentration increases, we should pay attention to the deviation of the actual Drude weight from its zeroth-order expansion in the impurity concentration D ≈ E F . Figure 5 (a) shows D(E F ) calculated from equation (16) as a function of the concentration c for various values of E F . This deviation increases with c, and it is most pronounced for E F near the resonance energy. The sign of the deviation changes when the Fermi energy crosses the resonance energy. Figure 5 (b) shows the dependence of D on the Fermi energy for a fixed concentration c = 10 −3 .
This difference can be understood using the first order expansion (17). As we increase the impurity concentration, the absolute value of Re Σ(E F ) becomes comparable to E F near the resonance. The real part Re Σ(E F ) changes its sign at the resonance energy (see figure 2), so D < E F for E F smaller than the resonance energy, and vice versa.
We note that D = D(E F ) has a plateau in the vicinity of the resonance energy, as can be seen in figure 5 (b). This plateau becomes flatter as we increase the impurity concentration c, so one could expect a slow change in D with a large change of E F . This means that one cannot properly determine the Fermi energy E F in the resonance region by measuring the Drude weight D.
When exceeding certain impurity concentration c χ , we expect the slope of the D = D(E F ) curve to become negative. We end up with a non-monotonous function, so one cannot determine E F from D uniquely near the resonance. The curves for different E F in figure 5 (a) cross each other as we approach this concentration. However, for concentrations c ∼ c χ , the applicability criterion for the ATA does not hold.
We can determine c χ by using the first-order expansion in concentration (17). The corresponding estimation can be obtained from the condition of flatness of the function D = D(E) near the solution of the Lifshitz equation E r : This gives us the following estimate: For larger values of the perturbation V L , the value of critical concentration is smaller. One can compare it to the spectrum transformation concentration c ST [25]:

33703-9
As we apply the strain, the resonance shifts towards the Dirac point. The dependence D = D(E F ) reflects this shift, as one can see in figure 5 (b). Thus, given that E F remains constant, one should expect the Drude weight D to change. This difference in D for ε = 0 and ε = 5% can be seen in figures 5 (a), 5 (b). As we can see in figure 5 (b), the most significant difference is observed when the Fermi energy E F approaches the resonance energy. However, the Drude weight shift fades as we move away from the resonance. In addition to that, the apparent flatness of the curve near the resonance increases with increasing the strain.
The results of the previous subsection still hold, but one should not assume that D = E F , as it is commonly accepted for the clean graphene. For E F lower than the resonance energy, where we expect the gain in the Drude width Γ (i.e., > 0), the Drude weight considerably exceeds the corresponding value of the Fermi energy [see figure 5 (b)]. Consequently, the values of D for which the maximum of is reached will exceed the corresponding values of E F .

Conclusion
In this paper, the influence of impurities on the Drude peak parameters was studied. Using the Lifshitz model and the average T-matrix approximation, we studied the self-energy function to establish the resonant behaviour of the partial scattering rate for point defects, treated as a function of the Fermi energy.
Modelling the uniform uniaxial strain by modification of the lattice hopping integrals, we have found that the effect of the strain on the electronic spectrum can be described by introduction of the effective bandwidth. The corresponding bandwidth variation was shown to lead to the shift in the position of the impurity resonance towards the Dirac point with increasing magnitude of strain. Although this shift is rather moderate in absolute value compared to the resonance energy, it is shown to significantly increase the Drude width. Thus, the Drude peak has a resonance behaviour as a function of Fermi energy. This effect is expressed more pronouncedly for the resonances which are located closer to the Dirac point, i.e., for larger impurity potentials.
In addition, it was shown that the Drude weight is not always a linear function of the Fermi energy. In the extreme case, the Drude weight can have a plateau feature as a function of the Fermi energy, which hinders the determination of the Fermi energy from the Drude weight data. Furthermore, the shift in the resonance energy due to a strain leads to a slight shift in the Drude weight at the constant Fermi energy. The magnitude and the direction of a shift was connected to the impurity resonance energy behaviour.
Although the article [10] gave us an impetus to write the present paper, we did not pose the goal to provide an exhaustive explanation of the experimental results presented in it. This was a reason behind the choice for impurity model. Instead, on a simple example, we made an attempt to show that for an impure graphene, the resonant behaviour of the electronic self-energy function leads to the resonant dependency of the Drude weight on the Fermi energy, and to the subsequent intense change of the Drude weight due to the applied strain.