Low-energy theory for strained graphene: an approach up to second-order in the strain tensor

TL;DR

This paper develops a second-order low-energy theory for uniformly strained graphene, revealing anisotropic Fermi velocities, optical properties, and pseudomagnetic effects, advancing understanding of strain-induced electronic phenomena.

Contribution

It introduces a second-order effective Dirac Hamiltonian for strained graphene, capturing trigonal symmetry and strain direction dependence absent in first-order models.

Findings

Anisotropic Fermi velocity tensor derived

Optical conductivity tensor calculated

Dirac point shift related to pseudomagnetic fields

Abstract

An analytical study of low-energy electronic excited states in an uniformly strained graphene is carried out up to second-order in the strain tensor. We report an new effective Dirac Hamiltonian with an anisotropic Fermi velocity tensor, which reveals the graphene trigonal symmetry being absent in low-energy theories to first-order in the strain tensor. In particular, we demonstrate the dependence of the Dirac-cone elliptical deformation on the stretching direction respect to graphene lattice orientation. We further analytically calculate the optical conductivity tensor of strained graphene and its transmittance for a linearly polarized light with normal incidence. Finally, the obtained analytical expression of the Dirac point shift allows a better determination and understanding of pseudomagnetic fields induced by nonuniform strains.

Equations

$$H(\bm{k})=-\sum_{n=1}^{3}t_{n}\left(\begin{array}[]{cc}0&e^{-i\bm{k}\cdot\bm{\delta}_{n}^{\prime}}\\ e^{i\bm{k}\cdot\bm{\delta}_{n}^{\prime}}&0\end{array}\right),$$

$${\mathbf{\delta }}_1=\frac{a_0}{2}(\sqrt{3},1),{\mathbf{\delta }}_2=\frac{a_0}{2}(-\sqrt{3},1),{\mathbf{\delta }}_3=a_0(0,-1),$$

$$t_n=t_0{e}^{-\beta (|{\mathbf{\delta }}_n^{\prime }|/a_0-1)},$$

$$E(\mathbf{k})=\pm |t_1{e}^{i\mathbf{k}\cdot {\mathbf{\delta }}_1^{\prime }}+t_2{e}^{i\mathbf{k}\cdot {\mathbf{\delta }}_2^{\prime }}+t_3{e}^{i\mathbf{k}\cdot {\mathbf{\delta }}_3^{\prime }}|,$$

$$\sum _{n=1}^3{t}_n{e}^{i{\mathbf{K}}_D\cdot {\mathbf{\delta }}_n^{\prime }}=\sum _{n=1}^3{t}_n{e}^{i{\mathbf{K}}_D\cdot (\bar{\mathbf{I}}+\overline{\mathbf{ϵ}})\cdot {\mathbf{\delta }}_n}=\sum _{n=1}^3{t}_n{e}^{i\mathbf{G}\cdot {\mathbf{\delta }}_n}=0,$$

$$\mathbf{G}={\mathbf{K}}_0+{\mathbf{A}}^{(1)}+{\mathbf{A}}^{(2)}+\mathcal{O}({\overline{\mathbf{ϵ}}}^3),$$

$$A_x^{(1)}+i{A}_y^{(1)}=\frac{\beta }{2a_0}({ϵ}_{xx}-{ϵ}_{yy}-2i{ϵ}_{xy}),$$

$$A_x^{(2)}+i{A}_y^{(2)}=\frac{\beta (4\beta +1)}{16a_0}({ϵ}_{xx}-{ϵ}_{yy}+2i{ϵ}_{xy}{)}^2.$$

$$H=\hslash v_0\mathbf{\tau }\cdot (\bar{\mathbf{I}}+\overline{\mathbf{ϵ}}-\beta \overline{\mathbf{ϵ}}-\beta {\overline{\mathbf{ϵ}}}^2-\beta {\bar{\mathbf{\kappa }}}_1+{\beta }^2{\bar{\mathbf{\kappa }}}_2)\cdot \mathbf{q},$$

$$\bar{\bm{\kappa}}_{1}=\frac{1}{8}\left(\begin{array}[]{cc}(\epsilon_{xx}-\epsilon_{yy})^{2}&-2\epsilon_{xy}(\epsilon_{xx}-\epsilon_{yy})\\ -2\epsilon_{xy}(\epsilon_{xx}-\epsilon_{yy})&4\epsilon_{xy}^{2}\end{array}\right),$$

$$\bar{\bm{\kappa}}_{2}=\frac{1}{4}\left(\begin{array}[]{cc}\epsilon_{xx}^{2}-\epsilon_{yy}^{2}+2\epsilon_{xx}\epsilon_{yy}+2\epsilon_{xy}^{2}&4\epsilon_{xx}\epsilon_{xy}\\ 4\epsilon_{xx}\epsilon_{xy}&2(\epsilon_{yy}^{2}-\epsilon_{xy}^{2})\end{array}\right).$$

$$\bar{\mathbf{v}}=v_0(\bar{\mathbf{I}}+\overline{\mathbf{ϵ}}-\beta \overline{\mathbf{ϵ}}-\beta {\overline{\mathbf{ϵ}}}^2-\beta {\bar{\mathbf{\kappa }}}_1+{\beta }^2{\bar{\mathbf{\kappa }}}_2),$$

$$\bar{\bm{\epsilon}}(\theta)=\epsilon\left(\begin{array}[]{cc}\cos^{2}\theta-\nu\sin^{2}\theta&(1+\nu)\cos\theta\sin\theta\\ (1+\nu)\cos\theta\sin\theta&\sin^{2}\theta-\nu\cos^{2}\theta\end{array}\right),$$

$$\tan 2{\theta }_v=\frac{2v_{xy}}{v_{xx}-v_{yy}},$$

$$△\theta$$

$$H=\hslash v_0\mathbf{\tau }\cdot (\bar{\mathbf{I}}+\bar{\mathbf{\Delta }})\cdot \mathbf{q},$$

$$\bar{\mathbf{\sigma }}(\omega )$$

$$\bar{\mathbf{\Delta }}=\overline{\mathbf{ϵ}}-\beta \overline{\mathbf{ϵ}}-\beta {\overline{\mathbf{ϵ}}}^2-\beta {\bar{\mathbf{\kappa }}}_1+{\beta }^2{\bar{\mathbf{\kappa }}}_2,$$

$$\bar{\mathbf{\sigma }}(\omega )$$

$$T={\left(1+\frac{\text{mbox}Re[{\sigma }_{xx}{\cos }^2{\theta }_i+{\sigma }_{yy}{\sin }^2{\theta }_i+{\sigma }_{xy}\sin 2{\theta }_i]}{2{\epsilon }_{0}c}\right)}^{-2}$$

$$T$$

$$△T(\theta )=2\pi \alpha \tilde{\beta }(1+\nu )ϵ+\pi \alpha (1+\gamma \cos 6\theta )(1+\nu {)}^2{ϵ}^2.$$

$$\sum _{n=1}^3{t}_n{e}^{i{\mathbf{K}}_D\cdot {\mathbf{\delta }}_n^{\prime }}=\sum _{n=1}^3{t}_n{e}^{i{\mathbf{K}}_D\cdot (\bar{\mathbf{I}}+\overline{\mathbf{ϵ}})\cdot {\mathbf{\delta }}_n}=\sum _{n=1}^3{t}_n{e}^{i\mathbf{G}\cdot {\mathbf{\delta }}_n}=0,$$

$$\mathbf{G}={\mathbf{K}}_0+{\mathbf{A}}^{(1)}+{\mathbf{A}}^{(2)}+\mathcal{O}({\overline{\mathbf{ϵ}}}^3),$$

$$t_n=t_0[1+{\Delta }_n^{(1)}+{\Delta }_n^{(2)}+\mathcal{O}({\overline{\mathbf{ϵ}}}^3)],$$

$$\sum _{n=1}^3[{\Delta }_n^{(1)}+i{\mathbf{A}}^{(1)}\cdot {\mathbf{\delta }}_n]e^{i{\mathbf{K}}_0\cdot {\mathbf{\delta }}_n}=0.$$

$$\sum _{n=1}^3$$

$$t_n/t_0$$

$${\Delta }_n^{(1)}=-\frac{\beta }{a_0^2}{\mathbf{\delta }}_n\cdot \overline{\mathbf{ϵ}}\cdot {\mathbf{\delta }}_n,$$

$${\Delta }_n^{(2)}=-\frac{\beta }{2a_0^4}(\overline{\mathbf{ϵ}}\cdot {\mathbf{\delta }}_n{)}^2+\frac{\beta (\beta +1)}{2a_0^4}({\mathbf{\delta }}_n\cdot \overline{\mathbf{ϵ}}\cdot {\mathbf{\delta }}_n{)}^2.$$

Full text

I Introduction

Given the striking interval of elastic response of graphene Galiotis et al. (2015); Daniels et al. (2015), can withstand a reversible stretching up to $25\,\%$, strain engineering has been widely used to improve and/or to tune its electronic, thermal, chemical and optical properties Bissett et al. (2014); Zhang and Zhang (2015); Si et al. (2016); Amorim et al. (2016); Naumis et al. (2016). For instance, theoretical predictions have been made of a band-gap opening by large uniaxial strains from both tight-binging approach Pereira et al. (2009) and density functional theory Choi et al. (2010a), whenever the strain produces such a Hamiltonian modification beyond the inequalities obtained by Hasegawa, et al. Hasegawa et al. (2006) The emergence of the pseudomagnetic field caused by a nonuniform strain is possibly the most interesting strain-induced electronic effect, due to the possibility of observing a pseudoquantum Hall effect under zero external magnetic fields Guinea et al. (2010a, b). Nowadays, the transport signatures of the such fictitious fields are actively investigated Gradinar et al. (2013); Bahamon et al. (2015); Burgos et al. (2015); Settnes et al. (2016a); Stegmann and Szpak (2016); Settnes et al. (2016b); Carrillo-Bastos et al. (2016); Georgi et al. (2016). Moreover, from a view point of basic research, strained graphene opens an opportunity to explore mixed Dirac–Schrödinger Hamtiltonian de Gail et al. (2012), fractal spectrum Naumis and Roman-Taboada (2014), superconducting states Kauppila et al. (2016), magnetic phase transitions López-Sancho and Brey (2016), metal-insulator transition Tang et al. (2015), among others exotic behaviours.

The concept of strain engineering has been also extended to the optical context Bae et al. (2013); Ni et al. (2014); Chen et al. (2016); Rakheja and Sengupta (2016). The optical properties of graphene are ultimately provided by its electronic structure, which can be modified by strain. For example, pristine graphene presents a transparency defined by fundamental constants, around $97.7\%$, over a broad band of frequencies Nair et al. (2008). This remarkable feature is essentially a consequence of its unusual low-energy electronic band structure around the Dirac points. Under uniform strain, such conical bands are deformed which produces anisotropy in the electronic dynamics Oliva-Leyva and Naumis (2013). Accordingly, this effect gives rise an anisotropic optical conductivity of strained graphene Pellegrino et al. (2010); Pereira et al. (2010); Oliva-Leyva and Naumis (2014a, b) and, therefore, a modulation of its transmittance as a function of the polarization of the incident light, as experimentally observed Ni et al. (2014). From a theoretical viewpoint, this optoelectronic behaviour of strained graphene has been quantified by continuum approaches up to first-order in the strain tensor Pereira et al. (2010); Oliva-Leyva and Naumis (2014a, b, 2015a). However, nowadays there are novel methods for applying uniaxial strain larger than $10\%$ in a nondestructive and controlled manner Pérez Garza et al. (2014). So, a low-energy continuum theory for the electronic and optical properties of strained graphene, up to second-order in the strain tensor, seems to be needed Masir et al. (2013); Li (2014); Crosse (2014); Ray et al. (2016).

In this paper, we derive the effective Dirac Hamiltonian for graphene under uniform strain up to to second-order in the strain tensor. For this purpose, we start from a nearest-neighbor tight-binding model and carry out an expansion around the real Dirac point. Unlike previous approaches to the first-order in strain, we show how the obtained low-energy Hamiltonian reveals the trigonal symmetry of the graphene. Also, we calculate the optical conductivity of strained graphene and characterize its transmittance for a uniaxial strain up to second-order in the stretching magnitude. These findings describe in a more accurate form the electronic and optical properties of strained graphene and, hence, can be potentially utilized towards novel optical characterizations of the strain state of graphene.

II Tight-binding model as starting point

Strain effects on electronic properties of graphene are usually captured by using a nearest-neighbor tight-binding model Pereira et al. (2009); Qi et al. (2013); Sloan et al. (2013); Settnes et al. (2016b). Within this approach, one can demonstrate that the Hamiltonian in momentum space for graphene under a uniform strain is given by Pereira et al. (2009); Oliva-Leyva and Naumis (2013)

[TABLE]

where the strained nearest-neighbor vectors are obtained by $\bm{\delta}_{n}^{\prime}=(\bar{\bm{I}}+\bar{\bm{\epsilon}})\cdot\bm{\delta}_{n}$, being $\bar{\bm{I}}$ the $(2\times 2)$ identity matrix and $\bar{\bm{\epsilon}}$ the rank-two strain tensor, whose components are independent on the position. Here, we choose the unstrained nearest-neighbor vectors as

[TABLE]

where $a_{0}$ is the intercarbon distance for pristine graphene. Thus, the $x$ $(y)$ axis of the Cartesian coordinate system is along the zigzag (armchair) direction of the honeycomb lattice. Owing to the changes in the intercabon distance, the nearest-neighbor hopping parameters are modified. Here we consider this effect by means of the commonly used model Papaconstantopoulos et al. (1998); Pereira et al. (2009); Ribeiro et al. (2009); Settnes et al. (2016a)

[TABLE]

where $t_{0}=2.7\mbox{eV}$ is the hopping parameter for pristine graphene and $\beta\approx 3$.

From equation (1) follows that the dispersion relation near the Fermi energy of graphene under uniform strain is given by two bands,

[TABLE]

which remains gapless as long as the triangular inequalities, $|t_{1}-t_{2}|\leq|t_{3}|\leq|t_{1}+t_{2}|$, are satisfied Hasegawa et al. (2006). Evaluating equation (4) for uniaxial strains, V. Pereira, et al., found the minimum uniaxial deformation that leads to the gap opening is about $23\%$ Pereira et al. (2009). This result is confirmed by the ab initio calculations, finding that this gap in strained graphene requires deformations larger than $20\%$ Choi et al. (2010b); Ni et al. (2009). Therefore, the use of an effective Dirac Hamiltonian obtained from equation (1) is justified for uniform deformations up to the order of $10\%$.

For this purpose, it is important to take into account a crucial detail: the strain-induced shift of the Dirac points in momentum space. In absence of deformation, the Dirac points $\bm{K}_{D}$ (determined by condition $E(\bm{K}_{D})=0$) coincide with the corners of the first Brillouin zone. Then, to obtain the effective Dirac Hamiltonian in this case, one simply expand the Hamiltonian (1) around such corners, e.g., $\bm{K}_{0}=(\frac{4\pi}{3\sqrt{3}a_{0}},0)$. However, in presence of deformations, the Dirac points do not coincide even with the corners of the strained first Brillouin zone Pereira et al. (2009); Li et al. (2010). Thus, to obtain the effective Dirac Hamiltonian, one should no longer expand the Hamiltonian (1) around $\bm{K}_{0}$. As demonstrated Oliva-Leyva and Naumis (2015b); Volovik and Zubkov (2014), such expansion around $\bm{K}_{0}$ yields an incorrect derivation of the anisotropic Fermi velocity. The appropriate procedure is to find first the new positions of the Dirac points and then carry out the expansion around them Yang (2011); Oliva-Leyva and Naumis (2015b); Volovik and Zubkov (2014); Bahat-Treidel et al. (2010); Volovik and Zubkov (2015).

III Effective Dirac Hamiltonian

As first step, we determine the new positions of Dirac points from the condition $E(\bm{K}_{D})=0$, up to second order in the strain tensor, which is the leading order used throughout the rest of the paper. Essentially, we calculate the strain-induced shift of the Dirac point $\bm{K}_{D}$ from the corner $\bm{K}_{0}$ of the first Brillouin zone by using equation $E(\bm{K}_{D})=0$, which leads to

[TABLE]

where $\bm{G}\equiv(\bar{\bm{I}}+\bar{\bm{\epsilon}})\cdot\bm{K}_{D}$ is the effective Dirac point. As demonstrated in Appendix A, $\bm{G}$ can be expressed as

[TABLE]

where

[TABLE]

and

[TABLE]

Notice that the correction up to first order, $\bm{A}^{(1)}$, coincides with the value previously reported Oliva-Leyva and Naumis (2013), which is interpreted as a gauge field for nonuniform deformations Yang (2011); Oliva-Leyva and Naumis (2015b); Volovik and Zubkov (2014). On the other hand, the expression (8) for the second-order correction $\bm{A}^{(2)}$ is one of the main contributions of this work. To demonstrate its relevance, we numerically calculate the positions of $\bm{G}$ for two deformations and compare them with the analytical results given by (6-8). As illustrated in Fig. 1(a) for a uniaxial strain along zigzag direction, the values of $G_{x}$ estimated up to first order in the strain magnitude $\epsilon$ (blue solid circles) clearly differ from the exact numerical values of $G_{x}$ (gray line) as $\epsilon$ increases, while the values of $G_{x}$ estimated up to second order (red open circles) show a significantly better approximation. The case of a shear strain is an even more illustrative example of the relevance of $\bm{A}^{(2)}$. According to the first-order correction, $G_{x}$ does not change under a shear strain, which is at variance with the exact numerical result displayed in Fig. 1(b). In contrast, the values of $G_{x}$ estimated up to second order present a good agreement with the numerical values over the studied range of $\epsilon$. Beyond the present work, the second-order correction $\bm{A}^{(2)}$ for nonuniform strain could be relevant to a more complete analysis of the strain-induced pseudomagnetic fields. For example, in presence of a deformation field given by $\bm{u}=(u(y),0)$, for which $\bar{\epsilon}_{xx}=\bar{\epsilon}_{yy}=0$ and $\bar{\epsilon}_{xy}=\partial_{y}u(y)/2$, the pseudomagnetic field $B_{ps}$, derived from the standard expression $B_{ps}=\nabla\times\bm{A}^{(1)}$, results equal to zero. However, if $\bm{A}^{(2)}$ is taken into account by means of the possible generalized expression $B_{ps}=\nabla\times\bm{A}^{(1)}+\nabla\times\bm{A}^{(2)}$, one can demonstrate that the resulting pseudomagnetic field $B_{ps}$ is not zero. The implications of this issue will be discussed with details in an upcoming work.

Knowing the position of the Dirac point $\bm{K}_{D}$, through equation (6), one can now proceed to the expansion of Hamiltonian (1) around $\bm{K}_{D}$, by means of $\bm{k}=\bm{K}_{D}+\bm{q}$, to obtain the effective Dirac Hamiltonian. Following this approach up to second order in the strain tensor $\bar{\bm{\epsilon}}$, the effective Dirac Hamiltonian can be written as (see Appendix B)

[TABLE]

where $v_{0}=3t_{0}a_{0}/2\hbar$ is the Fermi velocity for pristine graphene, $\bm{\tau}=(\tau_{x},\tau_{y})$ is a vector of ($2\times 2$) Pauli matrices describing the pseudospin degree of freedom,

[TABLE]

and

[TABLE]

It is important to emphasize that the explicit form of equations (10) and (11) is a consequence of the Cartesian coordinate system $xy$ chosen. For an arbitrary coordinate system $\tilde{x}\tilde{y}$, rotated by an angle $\vartheta$ respect to the system $xy$, the new expressions for $\bar{\bm{\kappa}}_{1}$ and $\bar{\bm{\kappa}}_{2}$ should be found by means of the transformation rules of a second order Cartesian tensor Barber (2003).

From equation (9) one can recognize the Fermi velocity tensor as

[TABLE]

which generalizes the expression, $v_{0}(\bar{\bm{I}}+\bar{\bm{\epsilon}}-\beta\bar{\bm{\epsilon}})$, for the Fermi velocity tensor up to first-order in the strain tensor reported in Refs. [Oliva-Leyva and Naumis, 2013; Volovik and Zubkov, 2014].

As a consistency test, let us consider an isotropic uniform strain of the graphene lattice, which is simply given by $\bar{\bm{\epsilon}}=\epsilon\bar{\bm{I}}$. Under this deformation, the new intercarbon distance $a$ is rescaled as $a=a_{0}(1+\epsilon)$, whereas the new hopping parameter $t$, expanding equation (3) up to second order in strain, results $t=t_{0}(1-\beta\epsilon+\beta^{2}\epsilon^{2}/2)$. Therefore, the new Fermi velocity, $v=3ta/2\hbar$, obtained straight away from the nearest-neighbor tight-binding Hamiltonian, takes the value $v=v_{0}(1-\beta\epsilon+\epsilon-\beta\epsilon^{2}+\beta^{2}\epsilon^{2}/2)$. This result can be alternatively obtained by evaluating our tensor (12) for $\bar{\bm{\epsilon}}=\epsilon\bar{\bm{I}}$.

The tensorial character of $\bar{\bm{v}}$ is due to the elliptic shape of the isoenergetic curves around $\bm{K}_{D}$. Notice that the principal axes of the Fermi velocity tensor up to first-order in the strain tensor, $v_{0}(\bar{\bm{I}}+\bar{\bm{\epsilon}}-\beta\bar{\bm{\epsilon}})$, are collinear with the principal axes of $\bar{\bm{\epsilon}}$. Therefore, within the effective low-energy Hamiltonian up to first-order in the strain tensor, the anisotropic electronic behaviour is only originated from the strain-induced anisotropy. Nevertheless, the terms $\bar{\bm{\kappa}}_{1}$ and $\bar{\bm{\kappa}}_{2}$ in equation (12) suggest that the second-order deformation theory might reveal the anisotropy (trigonal symmetry) of the underlying honeycomb lattice.

To clarify this issue, let us to consider graphene subjected a uniaxial strain such that the stretching direction is rotated by an arbitrary angle $\theta$ respect to the Cartesian coordinate system $xy$ (see Fig. 2). In this case, the strain tensor ($\bar{\bm{\epsilon}}$) in the reference system $xy$ reads

[TABLE]

where $\epsilon$ is the strain magnitude. Note that both $\bar{\bm{\epsilon}}(\theta)$ and $\bar{\bm{\epsilon}}(\theta+180^{\circ})$ represent physically the same uniaxial strain, which can be confirmed in equation (13). It is important to mention that for $\theta=n60^{\circ}$ ($\theta=90^{\circ}+n60^{\circ}$), being $n$ an integer, the stretching is along a zigzag (armchair) direction of graphene lattice.

As discussed above, under the strain (13), the Fermi velocity tensor up to first-order in the strain tensor, $v_{0}(\bar{\bm{I}}+\bar{\bm{\epsilon}}-\beta\bar{\bm{\epsilon}})$, is diagonal in the coordinate system $x^{\prime}y^{\prime}$, rotated by the angle $\theta$ respect to the coordinate system $xy$. However, the Fermi velocity tensor (12), up to second-order in the strain tensor, is diagonal in a coordinate system $x^{\prime\prime}y^{\prime\prime}$, rotated by an angle $\theta_{v}$ such that

[TABLE]

which determines the direction of lower electronic velocity. In the reciprocal space, the angle $\theta_{v}$ characterizes the pulling direction of isoenergetic curves, i.e., the principal axis of the isoenergetic ellipses, as illustrated in Fig. 2(d).

In Fig. (3), we show the difference $\triangle\theta=\theta_{v}-\theta$, numerically calculated from equation (14), as a function of the stretching direction $\theta$ for two different strain magnitudes $\epsilon=5\%$ and $10\%$. The observed six-fold behaviour of $\triangle\theta$ can be analytically evaluated by

[TABLE]

in good agreement with the numerical values, as shown in Fig. (3). From the last expression, it follows that the principal axes of the Fermi velocity tensor (12) are only collinear with the principal axes of $\bar{\bm{\epsilon}}(\theta)$ for $\theta=n30^{\circ}$, i.e., when the stretching is along the zigzag or armchair crystallographic directions. This result demonstrates that our Hamiltonian (9), a second-order deformation theory, reveals the trigonal symmetry of underlying honeycomb lattice.

IV Optical properties

An anisotropic Dirac system described by the effective Hamiltonian

[TABLE]

being $\bar{\bm{\Delta}}$ a symmetric ($2\times 2$) matrix such that $\Delta_{ij}\ll 1$, presents an anisotropic optical response captured by the conductivity tensor (see Appendix C):

[TABLE]

where $\omega$ is the frequency of the external electric field and $\sigma_{0}(\omega)$ is the optical conductivity of the unperturbed Dirac system, i.e., the optical conductivity of unstrained graphene. Equation (17) is a generalization up to second-order in $\bar{\bm{\Delta}}$ of previous expression until first-order in $\bar{\bm{\Delta}}$ for the optical conductivity of an anisotropic Dirac system, as it can be seen in equation (17) of Ref. [Oliva-Leyva and Naumis, 2016].

Now, comparing equations (9) and (16), the optical conductivity tensor $\bar{\bm{\sigma}}(\omega)$ of strained graphene is straightforward obtained by making the replacement:

[TABLE]

into equation (17). Regarding terms up to second-order in the strain tensor, it results

[TABLE]

where $\tilde{\beta}=\beta-1$. This equation generalizes previous works Pereira et al. (2010); Pellegrino et al. (2010); Oliva-Leyva and Naumis (2014a, b), in which the optical conductivity of graphene under uniform strain was reported up to first-order in the strain tensor.

Let us make a proof about the consistency of equation (19). When graphene is at half filling, i.e., the chemical potential equals to zero, the optical conductivity $\sigma_{0}(\omega)$ is frequency-independent and is given by the universal value $e^{2}/(4\hbar)$ Ziegler (2007); Gusynin et al. (2007). It is important to emphasize that this result is independent on the value $v_{0}$ of the Fermi velocity Stauber et al. (2015). Therefore, under an isotropic uniform strain $\bar{\bm{\epsilon}}=\epsilon\bar{\bm{I}}$, which only leads to a new isotropic Fermi velocity $v=v_{0}(1-\beta\epsilon+\epsilon-\beta\epsilon^{2}+\beta^{2}\epsilon^{2}/2)$, the optical conductivity does not change and remains equal to $\sigma_{0}=e^{2}/(4\hbar)$, at least within the Dirac cone approximation Stauber et al. (2015). In other words, any expression reported as optical conductivity tensor for uniformly strained graphene, as a function on the strain tensor, to be evaluated for $\bar{\bm{\epsilon}}=\epsilon\bar{\bm{I}}$ must give rise $\sigma_{0}\bar{\bm{I}}$, as occurred when one evaluates the tensor (19).

The optical conductivity up to first-order in the strain tensor, $\sigma_{0}[\bar{\bm{I}}+\tilde{\beta}\mbox{tr}(\bar{\bm{\epsilon}})\bar{\bm{I}}-2\tilde{\beta}\bar{\bm{\epsilon}}]$, under a uniaxial strain (13) can be characterized by $\sigma_{\parallel}=\sigma_{0}[1-\tilde{\beta}\epsilon(1+\nu)]$ and $\sigma_{\perp}=\sigma_{0}[1+\tilde{\beta}\epsilon(1+\nu)]$, where $\sigma_{\parallel}\ (\sigma_{\perp})$ is the optical conductivity parallel (perpendicular) to the stretching direction (see blue lines in Fig. 4). Within the first-order approximation, the optical conductivity along the stretching direction decreases by the same amount that the transverse conductivity increases, independently of $\theta$. This behaviour is modified when second-order terms are taken into account.

In Fig. 4(a), we plot the components of the optical conductivity tensor (19) versus the stretching magnitude $\epsilon$ for a uniaxial strain along the armchair direction. The perpendicular conductivity to the stretching direction, $\sigma_{xx}(\epsilon)$, does not have appreciable difference respect to the lineal approximation $\sigma_{0}[1+\tilde{\beta}\epsilon(1+\nu)]$ whereas the parallel conductivity, $\sigma_{yy}(\epsilon)$, noticeably differs from $\sigma_{0}[1-\tilde{\beta}\epsilon(1+\nu)]$ with increasing strain. On the other hand, Fig. 4(b) displays a contrary behaviour of the optical conductivity for a uniaxial strain along the zigzag direction. For this case, the parallel conductivity, $\sigma_{xx}(\epsilon)$, looks slight different from $\sigma_{0}[1-\tilde{\beta}\epsilon(1+\nu)]$ whereas the perpendicular conductivity, $\sigma_{yy}(\epsilon)$, is noticeably greater than $\sigma_{0}[1+\tilde{\beta}\epsilon(1+\nu)]$ with increasing strain. This increase of $\sigma_{yy}(\epsilon)$ respect to $\sigma_{0}[1+\tilde{\beta}\epsilon(1+\nu)]$ might help to give a better understanding of the change in the transmission of hybrid graphene integrated microfibers elongated along their axial direction Chen et al. (2016). For example, in Figure 2(b) of Ref. [Chen et al., 2016], it is possible to appreciate that the experimental data of this change gradually differ, with increasing strain, from the theoretical calculation using the first-order linear approximation $\sigma_{0}[1+\tilde{\beta}\epsilon(1+\nu)]$, which can be improved by considering the second-order contribution as shown in Fig. 4(b).

To complete our discussion about the emergence of the trigonal symmetry of graphene in the continuum approach presented here, we now study the transmittance of linearly polarized light on strained graphene. Considering graphene as a two-dimensional sheet with conductivity $\bar{\bm{\sigma}}$ and from the boundary conditions, vacuum-graphene-vacuum, for the electromagnetic field on the interfaces, the transmittance for normal incidence reads as Stauber et al. (2008); Oliva-Leyva and Naumis (2015a)

[TABLE]

where $\varepsilon_{0}$ is the vacuum permittivity, $c$ is the speed of light in vacuum and $\theta_{i}$ is the incident polarization angle. Note that for a pristine graphene with $\bar{\bm{\sigma}}=\sigma_{0}\bar{\bm{I}}$, equation (19) reproduces the experimentally observed constant transmittance $T=(1+\pi\alpha/2)^{-2}\approx 1-\pi\alpha$ over visible and infrared spectrum Nair et al. (2008), being $\alpha\approx 1/137$ the fine-structure constant. From equation (20) it can be seen that an anisotropic absorbance yields a periodic modulation of the transmittance as a function of the polarization direction $\theta_{i}$ Pereira et al. (2010); Oliva-Leyva and Naumis (2015a); Ni et al. (2014).

For the case of a uniaxial strain, and assuming the chemical potential equal to zero, from equations (13), (19) and (20) it follows that the transmittance up to second-order in the strain magnitude $\epsilon$ is given by

[TABLE]

where $\gamma=\beta(2\beta+1)/4$. Expression (21) reveals two new remarkable features in comparison with the first-order theory. As illustrated in Fig. 5, the transmittance mean value, $\langle T\rangle=1-\pi\alpha-\pi\alpha\tilde{\beta}^{2}(1+\nu)^{2}\epsilon^{2}/2$, has a negative shift with respect to the first-order average value $T_{0}=1-\pi\alpha$. Second, the transmittance oscillation amplitude ($\triangle T$) is determined by

[TABLE]

While the first-order expression for the transmittance oscillation amplitude, $2\pi\alpha\tilde{\beta}(1+\nu)\epsilon$, is independent on the stretching direction $\theta$, $\triangle T$ of equation (22) depends on $\theta$. For example, for a uniaxial strain along the zigzag (armchair) direction with $\theta=n60^{\circ}$ ($\theta=90^{\circ}+n60^{\circ}$), $\triangle T$ takes its highest (lowest) value, as displayed in Fig. 5. This strectching direction dependent $\triangle T$ might be used to confirm experimentally the present theory up to second-order in the strain tensor, as done for small strain less than $1\%$ Ni et al. (2014).

V Conclusion

We have analytically deduced a new effective Dirac Hamiltonian of graphene under a uniform deformation up to second-order in the strain tensor, including new Dirac-point positions that are qualitatively different from those predicted by first-order approaches, as occurred for the shear strain. Moreover, based on a detailed analysis about the anisotropic Fermi velocity tensor, we demonstrated how our second-order deformation theory reveals the trigonal symmetry of graphene unlike the previous first-order results.

We further derived, for the first time, analytical expressions for the high-frequency electric conductivity and light transmittance of a strained graphene up to second-order in the strain tensor. The magnitude of this transmittance oscillates according to the incident light polarization and the oscillation amplitude depends on the stretching direction, in contrast to the first-order prediction. In fact, within the first-order theory, the maximal transmittance occurs when the light polarization coincides to the stretching direction. However, the second-order theory predicts such coincidences only for stretching along zigzag and armchair directions. Therefore, the obtained light transmittance results can be experimentally verified by optical absorption measurements and they would be used for characterizing the deformation states of strained graphene. In general, the analytical study presented in this article has the advantage of being concise and establishes a reference point for upcoming numerical and experimental investigations.

It would be important to stress that the observed absence of lattice symmetry in the optical properties of strained graphene is due to the combination of the low-energy effective Dirac model and first-order approximation in the strain tensor. Such absence can be overcome by carrying out the study within the first-neighbour tight-binding model as occurred for high-energy electron excitations Pereira et al. (2009) or by introducing second-order effects in the strain tensor even within the simplest Dirac model, as done in this article. This finding of trigonal symmetry in optical response reveals the capability of low-energy effective Dirac theory to describe properly anisotropic electron behaviour in graphene under strong uniform deformations. However, a tight-binding model beyond nearest-neoghbour interactions would be required to analyze both the gap opening and the electron-hole spectrum symmetry induced by lattice strain Castro Neto et al. (2009). Finally, the present work can be extended to perform an analytical study of the pseudomagnetic fields induced by nonuniform strains.

Acknowledgements.

This work has been partially supported by CONACyT of Mexico through Project 252943, and by PAPIIT of Universidad Nacional Autónoma de México (UNAM) through Projects IN113714 and IN106317. Computations were performed at Miztli of UNAM. M.O.L. acknowledges the postdoctoral fellowship from DGAPA-UNAM.

Appendix A Dirac point position

Here we provide the derivation of expressions (6–8) of main text. Equation, $E(\bm{K}_{D})=0$, can be rewritten as

[TABLE]

where $\bm{G}\equiv(\bar{\bm{I}}+\bar{\bm{\epsilon}})\cdot\bm{K}_{D}$ is the effective Dirac point associated to a pristine honeycomb lattice with strained nearest-neighbor hopping integrals $t_{n}$. To solve equation (23) in a perturbative manner, we cast the position of $\bm{G}$ as

[TABLE]

where $\bm{A}^{(1)}$ ( $\bm{A}^{(2)}$) is the correction from first (second) order in the strain tensor. Similarly, we consider Taylor expansions of $t_{n}$, up to second order in strain tensor, in the form

[TABLE]

where $\Delta_{n}^{(1)}$ ($\Delta_{n}^{(2)}$) are terms of the first (second) order in the strain tensor.

Substituting equations (24) and (25) into equation (23), the coefficient of the first-order strain tensor should be equal to zero, which leads to

[TABLE]

Analogously, the coefficient of the second-order strain tensor should also be zero, yields

[TABLE]

From equation (26), $\bm{A}^{(1)}$ can be determined and it is used as input of equation (27) to obtain $\bm{A}^{(2)}$. To carry out this procedure, it is necessary to explicitly know $\Delta_{n}^{(1)}$ and $\Delta_{n}^{(2)}$ as functions of the strain tensor.

Expanding $t_{n}$, up to second order in strain tensor, gives

[TABLE]

Then, by comparing equations (25) and (28) one obtains

[TABLE]

and

[TABLE]

Finally, substituting $\Delta_{n}^{(1)}$ into equation (26), we get

[TABLE]

and consequently, using this result and the expression of $\Delta_{n}^{(2)}$, equation (27) can be rewritten as

[TABLE]

Note that equations (31) and (32) are the first- and second-order corrections to the Dirac point position given in equation (5) of the main text.

Appendix B Effective Dirac Hamiltonian

In order to derive the effective Dirac Hamiltonian given by equation (9) in the main text, we start from the tight-binding model in momentum space for graphene under a uniform strain,

[TABLE]

and we consider momenta close to the Dirac point $\bm{K}_{D}$, by means of the substitution $\bm{k}=\bm{K}_{D}+\bm{q}$. Then, expression (33) transforms as

[TABLE]

where $h=-\sum_{n=1}^{3}t_{n}e^{i(\bm{K}_{D}+\bm{q})\cdot(\bm{I}+\bar{\bm{\epsilon}})\cdot\bm{\delta}_{n}}$. Now, using equation (24), $h$ can be expanded up to first-order in $\bm{q}$ and second-order in $\bar{\bm{\epsilon}}$ as

[TABLE]

and substituting expression (25) for $t_{n}$ in equation (35), the expansion of $h$ results

[TABLE]

By taking into account equation (27) and $\sum_{n=1}^{3}e^{i\bm{K}_{0}\cdot\bm{\delta}_{n}}=0$, the $\bm{q}$-independent terms in the last expression are cancelled. Thus, $h$ can be rewritten as

[TABLE]

where

[TABLE]

[TABLE]

[TABLE]

[TABLE]

and

[TABLE]

being $\bm{Q}=\bm{q}\cdot\bar{\bm{\epsilon}}$. To simplify each term of $h$ in equation (37), we have used of $\bm{K}_{0}=(\frac{4\pi}{3\sqrt{3}a_{0}},0)$, equations (2) of the main text for $\bm{\delta}_{n}$ and the expressions obtained in the previous section for $\bm{A}^{(1)},\bm{A}^{(2)},\bm{\Delta}^{(1)}$ and $\bm{\Delta}^{(2)}$. In addition, note the same algebraic form between the initial expression of equation (39) and equation(38) if one defines $\bm{Q}=\bm{q}\cdot\bar{\bm{\epsilon}}$. This similarity is also observed between equations (41) and (40).

In consequence, using equations (38–42) we obtain the contribution of each term of $h$ to equation (34) as

[TABLE]

[TABLE]

[TABLE]

[TABLE]

and

[TABLE]

where $v_{0}=3t_{0}a_{0}/2\hbar$ is the Fermi velocity for pristine graphene, $\bm{\tau}=(\tau_{x},\tau_{y})$ is a vector of ($2\times 2$) Pauli matrices, $\bar{\bm{\kappa}}_{1}$ and $\bar{\bm{\kappa}}_{2}$ are respectively the matrices (10) and (11) of the main text. To obtain the expressions (43–47) in terms of Pauli matrices, we used the identity

[TABLE]

where $\chi_{ij}$ are the elements of an arbitrary ($2\times 2$) matrix $\bar{\bm{\chi}}$, being $\bar{\bm{\chi}}=\bm{I}$, $\bar{\bm{\chi}}=\bar{\bm{\epsilon}}$, $\bar{\bm{\chi}}=-\beta\bar{\bm{\epsilon}}$, $\bar{\bm{\chi}}=-\beta\bar{\bm{\epsilon}}^{2}$ and $\bar{\bm{\chi}}=-\beta\bar{\bm{\kappa}}_{1}+\beta^{2}\bar{\bm{\kappa}}_{2}$ for equations (43) to (47), respectively. It is worth mentioning that expression (43) is the effective Dirac Hamiltonian for pristine graphene, while (44) and (45) are corrections to first-order in $\bar{\bm{\epsilon}}$, previously derived in Ref. [Oliva-Leyva and Naumis, 2013]. The second-order corrections in $\bar{\bm{\epsilon}}$, equations (46) and (47), are among the principal contributions of our work.

Finally, combining equations (34), (37) and (43–47), we obtain the effective Dirac Hamiltonian for graphene under a uniaxial strain, up to second-order in the strain tensor $\bar{\bm{\epsilon}}$, given by

[TABLE]

which is the equation (9) reported in the main text.

Appendix C Optical conductivity of an anisotropic Dirac system

In this section, we derive the optical conductivity tensor $\bar{\sigma}_{ij}(\omega)$ of an anisotropic Dirac system described by the effective Hamiltonian

[TABLE]

where the anisotropic behaviour is expressed through the perturbation $\bar{\bm{\Delta}}$, which is a symmetric ($2\times 2$) matrix such that $\bar{\Delta}_{ij}\ll 1$. Essentially, we now extend, up to second-order in $\bar{\bm{\Delta}}$, a previous calculation of $\bar{\sigma}_{ij}(\omega)$ up to first-order in $\bar{\bm{\Delta}}$ reported in Ref. [Oliva-Leyva and Naumis, 2016].

Assuming that the considered system has linear response to an external electric field of frequency $\omega$, its optical conductivity $\bar{\sigma}_{ij}(w)$ can be calculated by combining the Hamiltonian (50) and the Kubo formula. Following the approach used in Refs. [Ziegler, 2006, 2007], $\bar{\sigma}_{ij}(\omega)$ can be expressed as a double integral with respect to two energies $E$, $E^{\prime}$:

[TABLE]

where $f(E)=(1+\exp[E/(k_{B}T)])^{-1}$ is the Fermi function at temperature $T$, Tr is the trace operator including the summation over the $\bm{q}$-space (as defined in equation (7) of Ref. [Ziegler, 2007]) and $v_{l}=i[H,r_{l}]$ is the velocity operator in the $l$-direction, with $l=x,y$.

To calculate the integral (51) it is convenient to make the change of variables

[TABLE]

which yields that the Hamiltonian (50) becomes $H=\hbar v_{0}\bm{\tau}\cdot\bm{q}^{*}$, corresponding to the case of a unperturbed and isotropic Dirac system, as unstrained graphene. At the same time, the velocity operator components transform as

[TABLE]

and analogously

[TABLE]

where $v_{x}^{*}=(\partial H/\partial q_{x}^{*})$ and $v_{y}^{*}=(\partial H/\partial q_{y}^{*})$ are the velocity operator components for the unperturbed Dirac system.

Then, substituting equations (53) and (54) into equation (51) we find

[TABLE]

and

[TABLE]

where $J$ is the Jacobian determinant of the transformation (52) originated by expressing the trace operator Tr of equation (51) in the new variables $\bm{q}^{*}$ and $\sigma_{0}(\omega)$ is the optical conductivity of the unperturbed Dirac system, i.e., the reported optical conductivity of unstrained graphene Ziegler (2007); Gusynin et al. (2007); Stauber et al. (2008). Note that equations (55-57) can be written in a compact manner as

[TABLE]

Now, if $J$ is expressed up to second-order in $\bar{\bm{\Delta}}$, results

[TABLE]

where $\mbox{tr}(\bar{\bm{\Delta}})=\bar{\Delta}_{xx}+\bar{\Delta}_{yy}$.

Finally, substituting equation (59) into (58) we obtain that the optical conductivity tensor of the anisotropic Dirac system, described by the Hamiltonian (50), is given by

[TABLE]

where the second term is the contribution to second-order in the perturbation $\bar{\bm{\Delta}}$ to the conductivity.