Photogalvanic effect in Weyl semimetals

TL;DR

This paper provides a theoretical analysis of how impurities affect the photogalvanic effect in Weyl semimetals, revealing complex behaviors influenced by impurity scattering and tilt vectors, with potential experimental relevance.

Contribution

It introduces a comprehensive kinetic model for the photogalvanic effect in Weyl semimetals considering impurity types and tilt effects, extending understanding of photocurrent generation mechanisms.

Findings

Nonmonotonic time evolution of circular PGE with intermediate internode scattering.

Tilt vector adds three independent photocurrent components, including a dominant one aligned with the tilt.

Skew scattering from tilt causes a perpendicular photocurrent component.

Abstract

We theoretically study the impact of impurities on the photogalvanic effect (PGE) in Weyl semimetals with weakly tilted Weyl cones. Our calculations are based on a two-nodes model with an inversion symmetry breaking offset and we employ a kinetic equation approach in which both optical transitions as well as particle-hole excitations near the Fermi energy can be taken into account. We focus on the parameter regime with a single photoactive node and control the calculation in small impurity concentration. Internode scattering is treated generically and therefore our results allow to continuously interpolate between the cases of short range and long range impurities. We find that the time evolution of the circular PGE may be nonmonotonic for intermediate internode scattering. Furthermore, we show that the tilt vector introduces three additional linearly independent components to the steady state photocurrent. Amongst them, the photocurrent in direction of the tilt takes a particular role inasmuch it requires elastic internode scattering or inelastic intranode scattering to be relaxed. It may therefore be dominant. The tilt also generates skew scattering which leads to a current component perpendicular to both the incident light and the tilt. We extensively discuss our findings and comment on the possible experimental implications.

Figures

Equations

$$H_{\mathrm{kin}}(\text{p})=\sum _{\zeta =\pm }\zeta [v(\text{p}-\zeta \text{b})\cdot \mathbf{\sigma }+u\hat{t}\cdot (\text{p}-\zeta \text{b})+\gamma ]\frac{1+\zeta {\kappa }_3}{2}.$$

$$H_{\mathrm{dis}}(\text{x})=V(\text{x})[1_{\kappa }+{\kappa }_1],$$

$${\partial }_{t}f(\text{p},t)+\dot{\text{p}}{\nabla }_{\text{p}}f(\text{p},t)=S{t}_{\mathrm{inj}}[f]+S{t}_{\mathrm{dis}}[f].$$

$$S{t}_{\mathrm{inj}}[f]$$

$$S{t}_{\mathrm{dis}}[f]=-{\int }_{\text{p}}{w}_{{\text{p}}^{\prime }\text{p}}f(\text{p},t)-w_{\text{p}{\text{p}}^{\prime }}f({\text{p}}^{\prime },t).$$

$$w_{{\text{p}}^{\prime },\text{p}}^{(s)}$$

$$w_{{\text{p}}^{\prime },\text{p}}^{(a)}$$

$$\times \left(-\frac{\pi u}{2v}\right){\nu }_{\zeta }({ϵ}_{\xi ,\zeta }(\text{p}))n_{\mathrm{imp}}{\mathcal{V}}_0^3\hat{t}\cdot (\hat{p}\times {\hat{p}}^{\prime }).$$

$${\nu }_{\zeta }(ϵ)=\frac{(\zeta ϵ-\gamma {)}^2}{2{\pi }^2{v}^3}.$$

$$f_{\xi }(\text{p})$$

$$I=\frac{2\pi v^2}{{\omega }^2}\xi e^2{\mathcal{E}}_{+}\cdot {\mathcal{E}}_{-},$$

$$I_i=\frac{2\pi v^2}{{\omega }^2}\xi \left\{-i{e}^2({\mathcal{E}}_{+}\times {\mathcal{E}}_{-}{)}_i\right\},$$

$$I_{ij}=\frac{2\pi v^2}{{\omega }^2}\xi \left\{-\text{Re}[e^2{\mathcal{E}}_{+,i}{\mathcal{E}}_{-,j}]\right\}.$$

$${\Gamma }_{\zeta }(ϵ)=\frac{2\pi }{3}{n}_{\mathrm{imp}}{\mathcal{V}}_0^2{\nu }_{\zeta }(ϵ),$$

$${\Gamma }_{\mathrm{sk},\zeta }(ϵ)=\frac{{\pi }^2}{3}{n}_{\mathrm{imp}}{\mathcal{V}}_0^3{\nu }_{\zeta }^2(ϵ),$$

$$\text{j}=-e\sum _{\xi =\pm 1}{\int }_{\text{p}}v\left[\xi \hat{p}+\frac{u}{v}\hat{t}\right]f_{\xi }(\text{p}).$$

$$=$$

$${\tau }_{\mathrm{CPGE}}\simeq \frac{1}{2}\sum _{ϵ={ϵ}_e,{ϵ}_h}1/{\Gamma }_{+},$$

$${\tau }_{\mathrm{skew}}\simeq \frac{1}{2}\sum _{ϵ={ϵ}_e,{ϵ}_h}{\Gamma }_{\text{sk},+}/{\Gamma }_{+}^2,$$

$${\tau }_{\hat{t}}\simeq -\frac{1}{2}\sum _{ϵ={ϵ}_e,{ϵ}_h}\frac{2}{3{\Gamma }_{-}\kappa }\text{sgn}\left(ϵ-\gamma \right),$$

$${\tau }_{\perp }\simeq -\frac{\kappa }{2}\sum _{ϵ={ϵ}_e,{ϵ}_h}\left[\frac{2{\Gamma }_{-}}{5{\Gamma }_{+}^2}+\frac{2}{5{\Gamma }_{-}}\right]\text{sgn}\left(ϵ-\gamma \right).$$

$${\tau }_{\mathrm{CPGE}}\simeq \sum _{ϵ={ϵ}_e,{ϵ}_h}\frac{1}{3[{\Gamma }_{-}+{\Gamma }_{+}]},$$

$${\tau }_{\mathrm{skew}}\simeq \sum _{ϵ={ϵ}_e,{ϵ}_h}(1-\kappa )\frac{2({\Gamma }_{\text{sk},+}{\Gamma }_{-}+{\Gamma }_{\text{sk},-}{\Gamma }_{+})}{9({\Gamma }_{+}+{\Gamma }_{-}{)}^2},$$

$${\tau }_{\hat{t}}\simeq -\sum _{ϵ={ϵ}_e,{ϵ}_h}\frac{9{\Gamma }_{-}+5{\Gamma }_{+}}{15{\Gamma }_{-}({\Gamma }_{+}+{\Gamma }_{-})}\text{sgn}\left(ϵ-\gamma \right),$$

$${\tau }_{\perp }\simeq -\sum _{ϵ={ϵ}_e,{ϵ}_h}\frac{2}{15({\Gamma }_{+}+{\Gamma }_{-})}\text{sgn}\left(ϵ-\gamma \right).$$

$$T(t)\simeq \sum _{ϵ={ϵ}_e,{ϵ}_h}\{\begin{array}{ll}\frac{1}{2{\Gamma }_{+}}\left[1-e^{-{\Gamma }_{+}t}\right],& \kappa \to 0,\\ \frac{1}{3({\Gamma }_{+}+{\Gamma }_{-})}\left[1-e^{-3t({\Gamma }_{+}+{\Gamma }_{-})/2}\right],& \kappa \to 1,\end{array}$$

$$H_{\mathrm{kin}}$$

$$H_{\mathrm{dis}}$$

$${\psi }_l(\text{x})=⟨\text{x}⟩l=e^{i\text{p}\text{x}}|u_{\xi ,\text{p}}⟩{\hat{e}}_{\zeta },$$

$$\ket{u_{\xi,\textbf{{p}}}}=\frac{1}{\sqrt{2(1+\xi\hat{p}_{3})}}\left(\begin{array}[]{c}\xi+\hat{p}_{3}\\ \hat{p}_{1}+i\hat{p}_{2}\end{array}\right)_{\sigma}.$$

$${ϵ}_l=\zeta (u\hat{t}\cdot \text{p}+\xi vp+\gamma ).$$

$$⟨u_{\xi ,\text{p}}⟩u_{{\xi }^{\prime },{\text{p}}^{\prime }}$$

$$|⟨u_{\xi ,\text{p}}⟩u_{{\xi }^{\prime },{\text{p}}^{\prime }}{|}^2$$

$$|u_{\xi ,\text{p}}⟩⟨u_{\xi ,\text{p}}|$$

$${\mathbf{\Omega }}_{\xi }(\text{p})$$

$$\text{p}=\frac{\xi (\zeta ϵ-\gamma )}{v}\left\{\text{k}-\xi \hat{t}\frac{u}{v}\right\}.$$

$$(dp)=\xi (\zeta ϵ-\gamma ){\nu }_{\xi \zeta }(ϵ)\frac{d^{3}k}{4\pi },$$

$$\delta ({ϵ}_l-ϵ)=\frac{1-\xi k_{t}u/v}{\xi (\zeta ϵ-\gamma )}\delta (|\text{k}|-1)$$

$$\int (dp)\delta ({ϵ}_l-ϵ)\cdots ={\nu }_{\xi \zeta }(ϵ){⟨(1-\xi \frac{u}{v}\hat{k}\cdot \hat{t})\dots ⟩}_{\hat{k}}.$$

$$\hat{p}=\hat{k}+\xi \frac{u}{v}[\hat{k}(\hat{k}\cdot \hat{t})-\hat{t}].$$

Full text

I Introduction

In the present days, Weyl semimetals [Burkov, 2016] enjoy significant scientific interest which is generated by a multitude of reasons. From the theoretical viewpoint, this class of materials is special since it represents a solid state realization of Weyl’s theory of chiral relativistic fermions [Weyl, 1929] – an elegant theory which regrettably seems to have lost its connection to fundamental particle physics after the discovery of neutrino oscillations [NeutrinoOscillations, ]. In the context of transport theories, the unique features related to the chiral anomaly [Bertlmann, 2000] are at the basis [Son and Spivak, 2013] of the experimentally observed giant magnetoresistance [Xiong et al., 2015]. From the viewpoint of applications, Weyl semimetals are attractive as their topologically protected band touching associated with a Berry curvature singularity may allow to explore and exploit novel regimes and phenomena of semiconductor physics. As we understand now, the appearance of Weyl nodes in systems lacking inversion or time reversal symmetry is far from being exceptional and the observation of Weyl physics has been already reported in various materials [WeylRealizations, ]. One particularly exciting phenomenon occurring in Weyl semimetals is the photogalvanic effect (PGE) [Cha, ; Wu et al., 2016; Morimoto et al., 2016; Jarillo-Herrero-Gedik2017, ], i.e. the generation of current due to the exposure to light. In contrast to most ordinary semiconductors, Weyl semimetals are susceptible to lowest frequencies and allow for novel technologies in the infrared. This consequence of the protected gapless spectrum comes along with the theoretical prediction of a quantized circular PGE [de Juan et al., 2016] which is awaiting its experimental verification.

The theoretical foundations of PGE go back to the seventies (for review see Refs. [Belinicher and Sturman, 1980; Sturman and Fridkin, 1992; Ivchenko, 2005]) while modern derivations rely on the nonlinear susceptibility framework [Sipe and Shkrebtii, 2000], Floquet theory [Morimoto and Nagaosa, 2016] and the Keldysh quantum kinetic equation approach [Koe, ]. In this paper, we concentrate on the dc response. It is common to distinguish two contributions: the injection current and the shift current. The former of the two represents the current generated by photoelectrons which are excited by optical (vertical) transitions with an anisotropic transition rate. In the steady state, this contribution diverges linearly when the current relaxation time becomes infinite. Differently stated, in the absence of a relaxation mechanism, the injection current grows linearly in time. In contrast, the shift current [von Baltz and Kraut, 1981; Belinicher et al., 1982] is due to the real space displacement that electrons undergo upon an optical transition. This contribution is always finite and thus subleading for long relaxation times.

Therefore, in the present paper we keep only the dc injection current. We focus on the parameter regime with a single photoactive node which materializes the predicted topological quantization [de Juan et al., 2016]. We investigate its relaxation due to impurities taking into account both intra and inter valley scattering in a model of two Weyl nodes. We scrutinize the time evolution of the current in response to a sudden illumination as well as its steady state characteristics. Furthermore we study the effect of a weak tilt in the Weyl spectrum. As we show, the interplay of impurity scattering and tilt leads to skew scattering and to the appearance of current components transversal to the direction of illumination.

It is also worthwhile to emphasize, that the current injection mechanism for the nonlinear PGE response at large frequencies is strongly different from the nonlinear intraband response at small frequencies considered in Refs. [Sodemann and Fu, 2015; Zyuzin and Zyuzin, 2017; Ishizuka et al., 2017]. We explicitly show that the latter exactly vanishes in our model for the case of absent tilt and intervalley scattering and that it is generically subdominant for frequencies exceeding the elastic scattering rate.

This article is structured as follows. In Sec. II we introduce the model under consideration and the assumptions for the calculation. In Sec. III we present a sketch of the technical derivation as well as the main results, which are subsequently discussed in Sec. IV. We conclude with a summary and an outlook about the experimental relevance of our findings, Sec. V. The Appendix contains details on the model (Appendix A), the elastic scattering rates (Appendix B), the relaxation of the photocarriers (Appendix C), the intraband response (Appendix D), and a microscopic tight binding model (Appendix E).

II Model and Assumption

The terms entering the free Hamiltonian $H_{0}(\textbf{{p}},\textbf{{x}})=H_{\rm kin}(\textbf{{p}})+H_{\rm dis}(\textbf{{x}})$ are as follows. The kinetic term is a generalized low energy version of lattice models presented in Refs. [Shapourian and Hughes, 2016; de Juan et al., 2016] and Appendix E

[TABLE]

By construction, the model defined by Eqs. (1b) is not invariant under time reversal, space inversion or rotation symmetries, see App. A for more details. We note that the tilt direction could also be assumed to be equal in both nodes [DetassisGrubinskas2017, ]. In contrast, we chose a model with opposite tilt in opposite cones, so that the tilt preserves inversion symmetry.

Our calculation of the photocurrent in Weyl semimetals relies on the following simplifying assumptions. First, we neglected any spatial dependence of the electric field which is justified for $v/c\rightarrow 0$, where $c$ is the speed of light. In doing so, we omit the photon drag effect. The assumption of a homogeneous electric field treatment within a 3D bulk Hamiltonian also relies on the second assumption of a long Thomas-Fermi screening length. In this limit, the contribution of surface states can be expected to be subdominant. Third, in neglecting all other sources of scattering we assume impurities to be the dominant source of photocurrent relaxation. While this statement is generally expected to be true at sufficiently low temperatures where electron-phonon and electron-electron scattering rates are small, we will explicitly show the limitations of this picture. Fourth, the current relaxation is calculated using a semiclassical kinetic equation approach which is justified if the mean free path of photoelectrons and photoholes exceeds their wavelength. This implies that the frequency of light is much larger than the elastic scattering rate. Formally, we implement this assumption by keeping only the leading $\mathcal{O}(n_{\rm imp}^{-1})$ terms in small impurity density. For the evaluation of quantum transition probabilities we truncate the T-matrix at next to leading (i.e. second) order in powers of weak impurity potentials. We also assume that the impurity potential is short ranged as compared to the wavelength of photocarriers (but it may be long-ranged as compared to $1/|\textbf{{b}}|$). Fifth, our calculation at the lowest temperatures is exponentially accurate if the energy difference between Fermi energy and the energy of photocarriers exceeds the temperature. Finally, all presented formulae are valid to linear order in $u/v\ll 1$, unless stated otherwise.

III Sketch of the calculation and results

In this section we present a sketch of the derivation of the photocurrent for a tilted, disordered Weyl semimetal. For the sake of a clearer presentation, we here restrict ourselves to the case of photocarrier generation at the $\zeta=+1$ cone only (see Fig. 1) and of absent intervalley scattering and relegate details and the more general case of finite intervalley scattering to Appendices B and C.

Our calculation is based on the Boltzmann kinetic equation [LyandaGellerAndreev2015, ; DeyoSpivak2009, ] describing the time evolution of the distribution function $f(\textbf{{p}},t)$

[TABLE]

Here, $\dot{\textbf{{p}}}=e\textbf{{E}}(t)$ and the collision integral contains two contributions. First, there is a term describing the excitation rates of photocarriers

[TABLE]

The second part of the collision integral describes the relaxation of photocarrieres by impurities

[TABLE]

We introduced the shorthand notation $\int_{\textbf{{p}}}=\int d^{3}p/(2\pi)^{3}$. Since impurity scattering is fully elastic, it may only relax anisotropic terms in the photoexcitation rates. The intervalley scattering discussed below and in Appendix C may also relax the valley imbalance of the isotropic contribution to the photoexcitation rates.

In the limit $u/v=0$, time reversal and rotational symmetries dictate that the scattering probability $w_{\textbf{{p}}^{\prime}\textbf{{p}}}$ contains only terms proportional to $1$ and $\textbf{{p}}\cdot\textbf{{p}}^{\prime}$ and thus $w_{\textbf{{p}}^{\prime}\textbf{{p}}}=w_{\textbf{{p}}\textbf{{p}}^{\prime}}$. This statement is fulfilled also when intervalley scattering is present, see Appendix B. However, the symmetry of the scattering probability is lost in the presence of a finite tilt $u/v\neq 0$. Here we present the dominant contributions of symmetric and antisymmetric scattering probability:

[TABLE]

Here, $\mathcal{V}_{\textbf{{p}}}$ is the Fourier transform of the impurity potential and the density of states (DOS) is

[TABLE]

The skew scattering probability obtained in Eq. (4b) relies on the third moment of the distribution function of the disorder potential. We note that, by means of both semiclassical [SinitsynSinova2007, ] and diffractive [AdoTitov2015, ; KoenigLevchenko2016, ] mechanisms, skew scattering may also occur for Gaussian disorder models. However, skew scattering probabilities due to Gaussian disorder contain an additional factor of $n_{\rm imp}$ and are thus subleading in the perturbation scheme employed here.

To obtain the photocurrent in the steady state, the solution for the distribution function can be obtained by equating the full collision integral to zero. We obtain the following nonequilibrium corrections to the Fermi-Dirac distribution function in band $\xi$

[TABLE]

We also introduced the intravalley scattering rates for momentum relaxation and skew scattering

[TABLE]

In Eq. (6a) the scalar part of the distribution function was omitted (it will be recovered shortly). Furthermore, Eq. (6a) was derived using the premise that only $\mathcal{O}(u/v)$ corrections proportional to odd powers of $\hat{p}$ enter the final expression for the photocurrent.

We insert this expression into the definition of the current

[TABLE]

Here the sum over $\xi$ reflects the two types of photocarriers that the injection term generates: electrons in band $\xi=+1$ at energy $\epsilon_{e}=\gamma+\omega/2$ and holes in band $\xi=-1$ at energy $\epsilon_{h}=\gamma-\omega/2$. As a side remark, we note that anomalous velocity terms $\delta\textbf{{v}}\sim\bm{\Omega}_{\xi}\times\textbf{{E}}$ are unimportant for the photocarriers, since their nonequilibrium distribution function Eq. (6a) is quadratic in electric field.

The final steady state current can be written as

[TABLE]

General formulae which interpolate between these two limits are presented in Eq. (70) of Appendix C and are plotted in Fig. 2. Note that all relaxation rates $\Gamma_{\zeta}(\epsilon)$ and $\Gamma_{\rm sk,\zeta}(\epsilon)$ introduced in these formulae are implicitly energy dependent. However, in the photoactive cone, the particle-hole symmetry about the Weyl node implies for both photoelectrons and photoholes equal scattering rates $\Gamma_{+}(\epsilon_{e})=\Gamma_{+}(\epsilon_{h})$, $\Gamma_{\rm sk,+}(\epsilon_{e})=\Gamma_{\rm sk,+}(\epsilon_{h})$. The rate $\tau_{\hat{t}}$ is not defined in a noninteracting model without intervalley scattering as it is determined by the relaxation time of the isotropic part in Eq. (3a). This can only be achieved by inelastic scattering or by means of intervalley scattering. We discuss this issue and other physical implication of the result presented in Eqs. (9), Eq. (70) and Fig. 2 in Sec. IV.

In Appendix C we also present a derivation of the time evolution after a sudden illumination in the limiting case $u/v=0$. Technically, this amounts to adding the solution of the homogenoeous Boltzmann equation (Liouville equation), Eq. (2), to the particular steady state solution, Eq. (6a). In terms of the final result Eq. (9a) this amounts to the following replacement $\tau_{\rm CPGE}\rightarrow T(t)$ with

[TABLE]

where we assumed the light to be switched on at time $t=0$. For a more general formula of $T(t)$ which interpolates between the two limits $\kappa=0$ and $\kappa=1$, we refer the reader to Eq. (73) of Appendix C as well as to Fig. 3. Again, the discussion of this result is relegated to Sec. IV.

Finally, the Appendix D contains a calculation of the intraband rectified response due to particle hole excitations near the Fermi surface. There we show, that this contribution vanishes in the absence of tilt and intervalley scattering and that it is generically suppressed by the parameter $\Gamma_{+}/\omega\ll 1$.

IV Discussion

In this section we analyze and discuss the physical content of the major results of this work, Eqs. (9) and (10) as well as Figs. 2 and 3.

IV.1 Limit $u/v=0$

We begin the discussion considering the limit of absent tilt $u/v=0$. In this case, the dc photocurrent is proportional to $-i(\bm{\mathcal{E}}_{+}\times\bm{\mathcal{E}}_{-})$ and is thus crucially dependent on circularly polarized light. This is the injection mechanism: vertical transitions inject photocarriers moving predominantly in direction $-i(\bm{\mathcal{E}}_{+}\times\bm{\mathcal{E}}_{-})$, see Eq. (3a). This leads to the linearly increasing time-dependent response, Fig. 3 at small $t\ll 1/\Gamma_{+}$. The slope is topologically quantized [de Juan et al., 2016] in the present case of a single photoactive Weyl node. Disorder exponentially relaxes the linear increase and a steady state is formed, see Fig. 3 at large $t\gg 1/\Gamma_{+}$. It is worthwhile to notice that $1/\Gamma_{+}$ is the “transport” mean free time, which is a factor of 3 larger than the “quantum” mean free time responsible for the level broadening. While these features are generic to any disorder model, intervalley scattering has the following three additional effects. First, it increases the scattering rate by increasing the density of available final states. This implies faster relaxation and a smaller steady state value. Second, since the relaxation rates at energies $\epsilon_{e}$ (photoelectrons) and $\epsilon_{h}$ (photoholes) are different in the $\zeta=-1$ node, intervalley scattering introduces two different relaxation rates for the circular PGE. Third, for intermediate $\kappa$ and provided $\Gamma_{-}<\Gamma_{+}$ intervalley scattering can lead to a peculiar nonmonotonic time dependence, see Eq. (73) and the $\kappa=1/3$ and $\kappa=2/3$ curves in Fig. 3. For the parameters chosen in our plot the nonmonotonicity stems from the photoholes living at energy $\epsilon_{h}$. Those scatter from the photoactive $\zeta=+1$ node to the $\zeta=-1$ node where the current is opposite and the decay much slower in view of a smaller DOS.

IV.2 Finite $u/v$

We now turn to the steady state solution at finite $u/v$. The finite tilt allows for the presence of additional contributions to the current, see Figs. 2 and 4.

First, skew scattering leads to a term $-i(\bm{\mathcal{E}}_{+}\times\bm{\mathcal{E}}_{-})\times\hat{t}$. Just as the circular PGE, it is only present for circular polarization of light. As we already mentioned, the momentum of excited photo electrons predominantly points in direction $-i(\bm{\mathcal{E}}_{+}\times\bm{\mathcal{E}}_{-})$. The time-reversal symmetry breaking tilt introduces a finite anisotropy leading to a preferred direction of momentum relaxation. The resultant imbalance between electrons moving towards $-i(\bm{\mathcal{E}}_{+}\times\bm{\mathcal{E}}_{-})\times\hat{t}$ as compared the opposing direction generates the skew scattering term proportional to $\tau_{\rm skew}$. It is in accordance with intuitive expectation that opposite tilts in opposite nodes imply antagonistic skew scattering contributions from opposite Weyl nodes. It is however surprising to find, that the $\mathcal{O}(u/v)$ skew scattering contribution exactly vanishes at strongest intervalley scattering $\kappa=1$, see the blue dot-dashed curve in Fig. 2.

Apart from skew scattering, the finite tilt in the spectrum also introduces terms which do not rely on circular polarization of light and which may be present also for linearly polarized light. We first focus on the term proportional to $\hat{t}(\bm{\mathcal{E}}_{+}\cdot\bm{\mathcal{E}}_{-})$. This term stems from the $(u/v)\hat{t}$ term in the definition of the current, Eq. (8). Therefore, it stems from the scalar (isotropic) part of the distribution function of photocarriers, or differently said, from the scalar term proportional to $\mathbf{1}$ in Eq. (3a). As we already mentioned, without intervalley scattering impurities cannot relax such a term. This is why the photocurrent in direction $\hat{t}$, represented by a red arrow in Fig. 4, is much stronger than all other contributions. The presence of intervalley scattering $\kappa>0$ introduces a finite $\tau_{\hat{t}}$, which however diverges as $1/\kappa$ for small $\kappa$, see Eq. (9d). This is represented in the red dashed curve in Fig. 2 which, for the sake of a better presentation, had to be downscaled by an extra factor of $20$. In the chosen parameter regime, $\tau_{\hat{t}}$ is particularly large, because $\Gamma_{-}$ is relatively small at $\epsilon_{h}$. Even though inelastic scattering is beyond the scope of the present paper, we mention that in practice, the steady state value of the current in $\hat{t}$ direction is determined by the smaller of $\tau_{\rm inel}$ (inelastic scattering time) and $\tau_{\hat{t}}$.

Finally, the steady state result, Eq. (9a) contains a term proportional to $\text{Re}[\bm{\mathcal{E}}_{+}\otimes\bm{\mathcal{E}}_{-}]\hat{t}$. This term stems from the tensor contribution to the photocarrier excitation rate, i.e. the $\hat{p}\otimes\hat{p}$ term in Eq. (3a) and is absent in the absence of intervalley scattering.

V Summary and Outlook

In summary, we have derived and analyzed the disorder induced relaxation of the photogalvanic effect (PGE) in a simple model for a Weyl semimetal. We took into account intra and intervalley scattering as well as leading order corrections due to finite tilt of the Weyl spectrum.

The major findings, which are pictorially summarized in Figs. 2-4, are that (i) intervalley scattering can lead to nonmonotonic time dependence of the circular PGE proportional to $-i(\bm{\mathcal{E}}_{+}\times\bm{\mathcal{E}}_{-})$; (ii) the finite tilt introduces additional current components proportional to (a) $-i(\bm{\mathcal{E}}_{+}\times\bm{\mathcal{E}}_{-})\times\hat{t}$ due to skew scattering, (b) $\hat{t}(\bm{\mathcal{E}}_{+}\cdot\bm{\mathcal{E}}_{-})$, and (c) $\text{Re}[\bm{\mathcal{E}}_{+}\otimes\bm{\mathcal{E}}_{-}]\hat{t}$. As we discussed, the contribution (ii b) in the direction of the tilt vector $\hat{t}$ is particularly inert to elastic scattering as it stems from isotropic generation of photo carriers and may thus be dominant. Contributions (i) and (ii a) involving $-i(\bm{\mathcal{E}}_{+}\times\bm{\mathcal{E}}_{-})$ change sign when the chirality of the photoactive Weyl node is reversed.

The observation of decisive qualitative and quantitative consequences of a weak spectral tilt corroborates the findings of earlier studies on different observables such as the linear conductivity tensor [TrescherBergholtz2015, ; Carbotte2016, ; SteinerPesin2017, ], the polarization function [DetassisGrubinskas2017, ] as well as transport through Weyl tunnel junctions [YesilyurtJalil2017, ]. Further implications on characteristic features of Weyl fermions, e.g. the natural optical activity [MaPesin2015, ], will be the subject of a separate studies [RouPesin2017, ].

We conclude our paper with an outlook on the experimental relevance of our findings. First of all, we comment on the minimal two Weyl node model that we are considering. Such materials are not discovered so far, but they were suggested theoretically in heterostructures [Burkov and Balents, 2011] and alloys [BulmashQi2014, ] of Chern insulator materials as well as, most recently, in certain magnetic Heusler compounds [Wang et al., 2016]. At the same time our simplified model may be applied to present day materials when pairs of Weyl nodes are well separated in momentum space.

Concerning optical experiments, we are aware of only two published results on the PGE in Weyl semimetals. The work of Ref. [Wu et al., 2016] is devoted to the second harmonic generation, i.e. a different observable as compared to the focus of this paper. In the second experiment [Jarillo-Herrero-Gedik2017, ] the chirality of the Weyl fermions, which is given by the sign of the Berry monopole charge, was inferred from the photocurrent response to mid-infrared light. Two key points are important: (i) Experiment reveals a current component perpendicular to the direction of light which strongly depends on the polarization (linear, circular positive, circular negative) and which is claimed to be proportional to the chirality of the photoactive nodes. As we show, skew scattering also induces a term of the very same tensor structure. (ii) Polarization independent photocurrent is observed in a direction perpendicular to both light and to the polarization dependent current. As we show, such terms may also be induced by the tilt. In addition, unpublished THz spectroscopy data on the dc PGE in TaAs presented at the APS march meeting 2017 [Patankar, ] indicates that the signals of radiated electric field in direction perpendicular or parallel to the polar axis show very different behavior and that for the latter case the signal is nearly independent of the polarization of incident light and is relatively long lived. The preliminary interpretation focused on anisotropic scattering and is thus related to the skew scattering mechanism investigated here. Another route for interpretation could be contributions similar to the terms $\hat{t}(\bm{\mathcal{E}}_{+}\cdot\bm{\mathcal{E}}_{-})$ considered in this paper. We repeat that these are indeed long lived and independent on the polarization of incident light. At the same time, we note that the tilt of the Weyl fermions in TaAs is substantial ($v/u>2$) and, therefore, our perturbative calculations should not be expected to quantitatively describe those experiments since our time-reversal breaking two node toy model has limited resemblance with the time reversal conserving 24 Weyl node material TaAs. However, we are convinced that the present results provide a proof of principle for various types of photocurrent contributions which will trigger future investigations.

VI Acknowledgments

We acknowledge useful discussions with A.V. Andreev, P.-Y. Chang, P. Coleman, L. Golub, A. Grushin, Y. Komijani, J. Moore, A. A. Zyuzin, and, in particular, with S. Patankar and L. Wu, who shared insights from ongoing experimental efforts.

The work of E.J.K., H.-Y.X. and A.L. at the University of Wisconsin-Madison was financially supported in part by NSF Grant No. DMR-1606517 and by the Office of the Vice Chancellor for Research and Graduate Education with funding from the Wisconsin Alumni Research Foundation. The work of D.A.P. at the University of Utah was supported by the NSF Grant No. DMR-1409089. At the latest stage of the project (when writing the manuscript), work by E.J.K. was carried out at Rutgers University, where financial support was provided by DOE, Basic Energy Sciences grant DE-FG02-99ER45790.

Appendix A Model

It is useful to perform a gauge transformation of Eq. (1b), such that

[TABLE]

The potential term $\Phi(\textbf{{x}},t)$ remains unaffected from the gauge transformation. The eigenstates of the clean Hamiltonian are characterized by the quantum number (multiindex) $l=(\zeta,\xi,\textbf{{p}})$, with $\zeta=\pm 1$, $\xi=\pm 1$ and $\textbf{{p}}\in\mathbb{R}^{3}$:

[TABLE]

where $\hat{e}_{+}=(1,0)_{\kappa}$ and $\hat{e}_{-}=(0,1)_{\kappa}$ and

[TABLE]

The unperturbed eigenenergies associated to the eigenstates $\ket{l}$ are

[TABLE]

The following relations will be useful ($\hat{p}_{1}+i\hat{p}_{2}=\sqrt{1-\hat{p}_{3}^{2}}e^{i\phi}$):

[TABLE]

A.1 Deformed coordinate system

In view of the tilt, the Fermi surfaces are ellipsoidal. Here we introduce a coordinate transformation which deforms the ellipse into a circle at the expense of certain Jacobians.

We first use that at a given energy $\epsilon$ the equienergy surface of a band characterized by $(\xi,\zeta)$ is an ellipse in the $p_{t}-p_{\perp}$ plane ($p_{t}=\hat{t}\cdot\textbf{{p}}$, $p_{\perp}=\sqrt{p^{2}-p_{t}^{2}}$) with center at $(p_{t},p_{\perp})=\xi(\zeta\epsilon-\gamma)(\frac{-\xi u}{v^{2}-u^{2}},0)$. The radius in $p_{t}$ ($p_{\perp}$) direction is $\xi(\zeta\epsilon-\gamma)\frac{v}{v^{2}-u^{2}}$ ($\xi(\zeta\epsilon-\gamma)\frac{1}{\sqrt{v^{2}-u^{2}}}$).

It is useful to choose the following energy dependent parametrization of the three vector p .

[TABLE]

The integration measure transforms as

[TABLE]

and, since a delta function in energy space can be written as

[TABLE]

we can write

[TABLE]

Note that, when k is a unit vector, we can write

[TABLE]

We will use the following parametrization of $f_{\xi,\zeta}(\textbf{{p}})$

[TABLE]

where the transformation $\textbf{{p}}\rightarrow(\epsilon,\hat{k})$ depends on the cone $\xi,\zeta$ under consideration. We repeat that the band index $\xi$ is uniquely determined by $\epsilon$ and $\zeta$.

Appendix B Elastic scattering rates

The scattering probability for scattering from the state $l$ to the state $l^{\prime}$ is

[TABLE]

where

[TABLE]

are the matrix elements of the retarded transition matrix. We next present restrictions on the scattering probability imposed by the symmetries discussed in the main text.

B.1 Symmetry considerations

As we mentioned in the main text, our model is not invariant under any of time reversal, rotation or inversion symmetries. We find

[TABLE]

We used the representations $\hat{T}=\sigma_{y}K$ (K is the complex conjugation); $\hat{I}=\kappa_{1}\mathcal{I}$ ($\mathcal{I}:\textbf{{x}}\rightarrow-\textbf{{x}}$) and $\hat{R}_{U}=U\mathcal{R}_{\underline{R}_{U}^{T}}$ ($\mathcal{R}_{\underline{R}_{U}}:\textbf{{x}}\rightarrow\underline{R}_{U}\textbf{{x}}$ and $U\bm{\sigma}U^{\dagger}=\underline{R}_{U}\bm{\sigma}$ is the $\text{SU}(2)$-$\text{SO}(3)$ homomorphism).

The symmetry operations on the wave function are as follows

[TABLE]

Here, again we used $\hat{p}_{1}+i\hat{p}_{2}=(1-\hat{p}_{3}^{2})e^{i\phi}$. For the rotational symmetry it is more useful to investigate its action on a projector onto an eigenstate

[TABLE]

We exploit those transformations in the analysis of the scattering probabilities in the limiting case $u/v=0$. For the present model, rotational symmetry plays a major role and we find that, for $u=0$ and after average over disorder configuration,

[TABLE]

which implies that $\overline{w_{l^{\prime}l}}$ is a scalar or pseudo scalar under transformations and can be expanded by the basis functions $1,\hat{p}\cdot\hat{p}^{\prime},\hat{p}\cdot\textbf{{b}},\hat{p}^{\prime}\cdot\textbf{{b}},(\hat{p}\times\hat{p}^{\prime})\cdot\textbf{{b}}$. We can use the fact that after disorder average, the average translational inveriance implies that b enters only in the form of $|\mathcal{V}_{\bm{b}}|^{2}$ under the assumption of sufficiently short ranged impurities. Therefore, rotational invariance implies that

[TABLE]

In particular, $w_{l^{\prime}l}$ is symmetric under the exchange of momenta.

We can further use the implication of TR symmetry and find

[TABLE]

so that

[TABLE]

and thus the absence of skew scattering: ${w_{l^{\prime}l}}={w_{ll^{\prime}}}$ in the limit $u=0$.

B.2 Scattering probabilities

B.2.1 Born approximation

The scattering probability in Born approximation can be readily evaluated from Eq. (22)

[TABLE]

B.2.2 Skew scattering

For the skew scattering we expand $T_{l^{\prime}l}\simeq\braket{l^{\prime}}{\hat{H}_{\rm dis}+\hat{H}_{\rm dis}G_{0}^{R}\hat{H}_{\rm dis}}{l}$ and obtain [Nagaosa et al., 2010]

[TABLE]

For our specific model this leads to

[TABLE]

We use that

[TABLE]

so that we find

[TABLE]

This concludes the derivation of Eq. (4b) of the main text.

B.3 Scattering rates

It will be useful to expand the distribution function for the band $(\xi,\zeta)$ by means of representations of the rotation group. We write

[TABLE]

Furthermore, the kinetic equation and the collision integral shall be written as a vector $\underline{f}=(f_{\zeta=+1},f_{\zeta=-1})^{T}$ in valley space and denote it by an underbar. Note that, at a given energy $\epsilon$ there are two bands involved, one at each valley.

B.3.1 Born approximation

We consider the Born scattering from impurities introduce scattering rates $\Gamma_{\zeta}=2\pi n_{\rm imp}|\mathcal{V}_{0}|^{2}\nu_{\xi,\zeta}(\epsilon)$ ($\xi$ is determined by $\zeta$ and $\epsilon$ uniquely).

We can write the collision integral entering the equation for electrons in pocket $\xi,\zeta$

[TABLE]

The notation $St[\{f\}]$ reflects that in general the collision integral depends on the whole set of distribution functions, i.e. in the present case on the distribution functions of both nodes. All energies are taken at given $\epsilon_{\xi,\zeta}$ which is either $\epsilon_{e}$ or $\epsilon_{h}$. We employ the valley space vector notation for the scalar part of the distribution function to $\mathcal{O}(u/v)$

[TABLE]

Here, the matrix $\underline{\xi}=\text{diag}(\xi_{+},\xi_{-})$. As we shall see later, for the vector part we don’t need $u/v$ corrections which contain an even power of vectors $\hat{k}$. Under this premise (reflected by the sign $\doteq$), only the nontilted contribution survives

[TABLE]

Next, we turn our attention to the matrix distribution function. Here we obtain

[TABLE]

The sum over repeated indices is to be understood. We will find later that the solution for the tensorial part of the distribution function $f^{ijk}$ is first order in $u/v$. Therefore, we only need the zeroth order in $u/v$ from the collision integral, i.e.

[TABLE]

B.3.2 Skew scattering

For the skew scattering the leading order is $(u/v)^{1}$ and we only need to consider the contributions of the scalar part (vanishes), of the vector part (nonvanishing) and of the matrix part (vanishes) of the distribution function.

After angular integration the skew-scattering collision integral for electrons in band $(\xi,\zeta)$ becomes

[TABLE]

We define the band dependent skew scattering rate $\Gamma_{\zeta}^{\rm(sk)}=\frac{\pi^{2}}{3}n_{\rm imp}\mathcal{V}_{0}^{3}\nu_{\xi\zeta}^{2}$. We further use that $\Gamma_{+}/\Gamma_{-}=\nu_{+}/\nu_{-}$ to reexpress the ratios of DOSs. Then

[TABLE]

Note that $\Gamma_{+}^{\rm(sk)}{\Gamma_{-}}/{\Gamma_{+}}=\Gamma_{-}^{\rm(sk)}{\Gamma_{+}}/{\Gamma_{-}}$.

Appendix C Relaxation of photocarriers

We use the notation $\zeta_{0}=\text{sgn}(\epsilon_{F})$ for photoactive node and we further write $\xi_{0}$ for the band index in that node at a given energy $\epsilon_{e}$ or $\epsilon_{h}$.

C.1 Injection term

The injection term, i.e. the collision term in the Boltzmann equation generating photoelectrons/photoholes in band $(\xi_{0},\zeta_{0})$ in vector notation is given by the Fermi’s golden-rule expression

[TABLE]

For the product of transition matrix elements

[TABLE]

we calculate symmetric (round brackets) and antisymmetric (square brackets) components separately

[TABLE]

This leads to Eq. (3a) of the main text. After the change of coordinates from $\hat{p}$ to $\hat{k}$ the injection term in vector notation can be written as

[TABLE]

where we introduced (recall $\mathbf{E}(t)=\sum_{\pm}\bm{\mathcal{E}}_{\pm}e^{\pm i\omega t}$)

[TABLE]

As we mentioned previously, we are only interested in the dc part of these expressions since $\omega/\Gamma_{\zeta_{0}}\gg 1$ implies that second harmonic generation is only weakly affected by disorder.

C.2 Static solution to $\mathcal{O}[(u/v)^{0}]$

We expand $\underline{f}=\underline{f}^{(0)}+(u/v)\underline{f}^{(1)}$ and keep only the zeroth order in $(u/v)$ in this section. Upon equating the full collision integral $\underline{St}_{\rm inj}+\underline{St}^{(2)}+\underline{St}^{(3a)}$ to zero we obtain

[TABLE]

The solutions for matrix and vector part immediately follow

[TABLE]

The symbol ”$\doteq$” means that this equation is valid only if $\underline{f}^{(0)}$ is multiplied by $(1,-1)$ from the left. This is, because the scalar symmetric part $(f_{\zeta=+}+f_{\zeta=-})$ corresponds to a zero mode of the collision integral:

[TABLE]

The inverse of $\underline{\Gamma}_{s}^{-1}$ means (symbolically in the subspace of the antisymmetric part $(f_{\zeta=+}-f_{\zeta=-})$)

[TABLE]

C.3 Static solution to $\mathcal{O}[(u/v)^{1}]$

As we shall show below, for the calculation of the current we only need terms containing an odd number of $\hat{k}$s in the linear corrections to the distribution function. We denote this assumption by a “$\doteq$” in the following equation for the first order terms

[TABLE]

After some algebra, this leads to

[TABLE]

C.4 Static current response

The current is defined as

[TABLE]

Since we found that distribution function is proportional to $\delta(\omega-2vp)$ it is useful to switch to spherical coordinates in p space (i.e. we transform back from the $\hat{k}$ representation). From this equation we readlily see that in the terms of first order in $(u/v)$ of $f(\textbf{{p}})$ we need only terms which contain odd powers of $\hat{p}$.

We will write

[TABLE]

Here we introduced

[TABLE]

as well as

[TABLE]

We used that both $\underline{\Gamma}_{t}$ and $\underline{\xi}$ are diagonal matrices.

In order to present the final result for the current we define $\tilde{\underline{\Gamma}}_{v}=\underline{\xi}{\underline{\Gamma}}_{v}\underline{\xi}$

[TABLE]

Then we obtain using $\tilde{\underline{\Gamma}}_{\rm sk}=\underline{\xi}{\underline{\Gamma}}_{\rm sk}\underline{\xi}$

[TABLE]

Here $\underline{\kappa}_{3}=\text{diag}(1,-1)_{\zeta}$. Bear in mind that the second harmonic term should not be inferred from this term. The prefactor ${{e^{3}}/{12\pi}}$ corresponds to the quantization unit presented in Ref. de Juan et al., 2016. We can readily read off the relaxation times introduced in Eq. (9a)

[TABLE]

The evaulation of these equations in the case $\zeta_{0}=+1$ leads to

[TABLE]

We note that the two relaxation times $\tau_{\hat{t}}$ and $\tau_{\perp}$ contain additional factors of $\text{sgn}({\epsilon-\gamma/2})$ and therefore have opposite signs for photoelectrons and photoholes.

C.5 Dynamic solution to $\mathcal{O}[(u/v)^{0}]$

In order to obtain the time dependence of the photocurrent after a sudden illumination, we add the solution of the homogeneous Boltzmann equation to our steady state solution, Eq. (54b). Under the assumption of thermal equilibrium at $t=0$ we obtain.

[TABLE]

Inserting this into the definition of the current leads to $\tau_{\rm CPGE}\rightarrow T(t)$ with

[TABLE]

We assume that $\zeta_{0}=+1$ and find, after some algebra, an expression by means of $\overline{\Gamma}=\frac{\Gamma_{+}+\Gamma_{-}}{2}$ and $\Delta\Gamma=\Gamma_{+}-\Gamma_{-}$ to be

[TABLE]

It is worth to notice, that \Big{[}1-\frac{\Delta\Gamma(1-\kappa)+\overline{\Gamma}\kappa}{\sqrt{(\Delta\Gamma)^{2}(1-\kappa)(1-2\kappa)+\overline{\Gamma}^{2}\kappa^{2}}}\Big{]} changes sign at $\Delta\Gamma=0$. This is the origin of the non-monotonic behavior reported in Fig. 3.

Appendix D Intraband response

Here we present details on the contribution from the Fermi surface. We concentrate on the limit $u/v=0$. The side jump contribution modifies the collision integral as follows

[TABLE]

Here, $\tilde{f}_{\zeta}(\epsilon,\hat{p})=f_{\zeta}(\textbf{{p}})$ and the sidejump is

[TABLE]

In this appendix, we investigate how the driving term proportional to $\dot{\textbf{{p}}}$ excites particle-hole pairs at the Fermi surface. We concentrate on the circular PGE. The according contribution to the distribution function is

[TABLE]

Here, $\underline{M}=\underline{\Gamma}_{v}/(2\pi n_{\rm imp})$, $i\underline{\Omega}=i\omega+\underline{\Gamma}_{v}$ and

[TABLE]

accounts for the side jump contribution. When evaluated by means of the definition of the current, which however now contains additional contributions:

[TABLE]

In the present case, the side jump accumulation is zero

[TABLE]

We obtain in linear response

[TABLE]

The second order response contains two contributions. First, the Berry curvature term leads to

[TABLE]

Clearly, in the clean limit contributions from opposite cones cancel up. Furthermore, there is a contribution from the usual group velocity combined with the side jump term

[TABLE]

In the second line, we introduced the dimensionless matrices

[TABLE]

Keep in mind that the prime ′ denotes derivative with respect to energy. In the limit $\kappa\rightarrow 0$ (no intervalley scattering) we have $\frac{v^{3}}{\gamma}[(\mathbf{1}+\kappa\underline{R})^{-1}\underline{S}]^{\prime}\underline{\nu}\rightarrow\text{diag}(1,-1)/[2\pi^{2}]$ so that

[TABLE]

Remarkably, the contribution from the side jump exactly cancels up the intrinsic contribution. Generally, we can state the the Fermi surface contributions are $\Gamma_{\zeta}/\omega\ll 1$ times smaller than the contributions from photocarriers.

Appendix E Microscopic Model

As the simplest model of the Weyl semimetal (WSM), we review the two-band multilayer Chern insulator Hamiltonian [Shapourian and Hughes, 2016; de Juan et al., 2016]. The disorder-free single-particle Hamiltonian on the cubic lattice reads

[TABLE]

where $\mathbf{r}=(x\,\,y\,\,z)$ is the coordinate of a site on the cubic lattice, $({\bf a}_{1}\,\,{\bf a}_{2}\,\,{\bf a}_{3})=(\hat{\mathbf{x}}\,\,\hat{\mathbf{y}}\,\,\hat{\mathbf{z}})$ the set of unit vectors, $C_{\mathbf{r}}^{\dagger}\equiv[c_{\uparrow}^{\dagger}(\mathbf{r})\,\,c_{\downarrow}^{\dagger}(\mathbf{r})]$ ($C_{\mathbf{r}}$) the creation (annihilation) operator of an electron spinor, and ${{\bm{\sigma}}}=(\sigma_{1},\,\sigma_{2},\,\sigma_{3})$ the Pauli matrices in the spin space. The Hamiltonian (102) describes a stack of Chern insulator layers [Eq. (102a)] coupled through tunneling [Eq. (102b)] and all the model parameters $t,\gamma,M$, and $\mu$ are real.

In reciprocal space the Hamiltonian takes the form

[TABLE]

with $\psi({\bf k})$ the two-spinor

[TABLE]

and the Hamiltonian matrix

[TABLE]

The energy dispersion is

[TABLE]

For $|M/t-2|<1$, a pair of Weyl valleys form about the momenta and energies

[TABLE]

where $\zeta=\pm 1$ indicates the valleys. For $|M/t-2|=1$ there is only one band-touching point and for $|M/t-2|>1$ band gap opens, for which the physics is less interesting.

In the Weyl semimetal regime $|M/t-2|<1$, we expand the Hamiltonian (105) about the Weyl point ${\bf K}_{\zeta}$ to obtain a low-energy continuous model. Substituting ${\bf k}={\bf K}_{\zeta}+{\bf p}$ into Eq. (105), for small momentum deviation $|{\bf p}|\ll K$ and up to linear order in $p_{x,y,z}$ we have

[TABLE]

and the corresponding low-energy two-spinor (104) takes the valley index $\zeta$ (chirality)

[TABLE]

Introducing the four-spinor

[TABLE]

and the Pauli matrices ${\bm{\kappa}}=(\kappa_{1},\,\kappa_{2},\,\kappa_{3})$ for valley space, we obtain the Weyl Hamiltonian

[TABLE]

Especially, at $M/t=2$, where $K=\pi/2$ and $E=\gamma$, the Weyl cones become isotropic and the Hamiltonian (E) takes the form of Eq. (1a) at $u/v=0$.