Quantum non-local theory of topological Fermi arc plasmons in Weyl semimetals

TL;DR

This paper develops a quantum non-local theory for Fermi arc plasmons in Weyl semimetals, revealing how topological surface states interact with bulk electrons and include quantum effects like penetration and dissipation.

Contribution

It introduces a novel quantum-mechanical non-local framework for topological Fermi arc plasmons in Weyl semimetals, incorporating quantum penetration and dissipation effects.

Findings

Quantum effects of Fermi arc penetration into bulk are significant.

Dissipation arises mainly from bulk electron-hole pair decay processes.

The theory captures non-local electron-electron interactions in Weyl semimetals.

Abstract

The surface of a Weyl semimetal (WSM) displays Fermi arcs, i.e. disjoint segments of a two-dimensional Fermi contour. We present a quantum-mechanical non-local theory of chiral Fermi arc plasmons in WSMs with broken time-reversal symmetry. These are collective excitations constructed from topological Fermi arc and bulk electron states and arising from electron-electron interactions, which are treated in the realm of the random phase approximation. Our theory includes quantum effects associated with the penetration of the Fermi arc surface states into the bulk and dissipation, which is intrinsically non-local in nature and arises from decay processes mainly involving bulk electron-hole pair excitations.

Equations

$${\mathcal{H}}_m(\mathbf{k})=\frac{\hslash v}{2b}(k_x^2-m){\sigma }_x+\hslash v(k_y{\sigma }_z+k_z{\sigma }_y).$$

$$\begin{array}{r}{u}_{\mathbf{k},\lambda }=\frac{1}{\sqrt{2}}\left(\begin{array}{c}\sqrt{1+\lambda \cos ({\beta }_{\mathbf{k}})}{e}^{-i{\varphi }_{\mathbf{k}}}\\ \lambda \sqrt{1-\lambda \cos ({\beta }_{\mathbf{k}})}\end{array}\right),\end{array}$$

$$\begin{array}{r}{\Psi }_{{\mathbf{k}}_{\parallel }}^{\mathrm{A}}(\mathbf{r})=\sqrt{\frac{b^2-k_x^2}{bS}}\left(\begin{array}{c}1\\ 0\end{array}\right)\Theta (-z)e^{z/\ell }{e}^{i{\mathbf{k}}_{\parallel }{\mathbf{r}}_{\parallel }}\end{array}$$

$${\Psi }_{\lambda ,\mathbf{k}}^{\mathrm{B}}({\mathbf{r}}_{\Vert },z)=\frac{\Theta (-z)}{\sqrt{V}}(u_{\mathbf{k},\lambda }{e}^{i{k}_{z}z}+r_{\mathbf{k}}{u}_{\bar{\mathbf{k}},\lambda }{e}^{-i{k}_{z}z})e^{i{\mathbf{k}}_{\parallel }\cdot {\mathbf{r}}_{\parallel }},$$

$$n_1(z,{\mathbf{q}}_{\parallel },\omega )={\int }_{-\infty }^{0}d{z}^{\prime }{\tilde{\chi }}_{\text{nn}}(z,z^{\prime },{\mathbf{q}}_{\parallel },\omega )V_{\mathrm{sc}}(z^{\prime },{\mathbf{q}}_{\parallel },\omega ).$$

$$V_{\mathrm{sc}}(z,{\mathbf{q}}_{\parallel },\omega )$$

$$1=\frac{2\pi e^2}{q_{\parallel }}{\chi }_{\mathrm{eff}}({\mathbf{q}}_{\parallel },\omega ),$$

$${\chi }_{\mathrm{eff}}({\mathbf{q}}_{\parallel },\omega )$$

$${\tilde{\chi }}_{nn}(q_z,q_z^{\prime },{\mathbf{q}}_{\parallel },\omega )$$

$${\chi }_{nn}^{(0)}(q_z,q_z^{\prime },{\mathbf{q}}_{\parallel },\omega )$$

$${\chi }_{nn}^{(0)}={\chi }_{\mathrm{AA}}^{(0)}+{\chi }_{\mathrm{BB}}^{(0)}+{\chi }_{\mathrm{AB}}^{(0)}.$$

$$\mathrm{Re}[{\chi }_{\mathrm{AA}}^{\mathrm{eff}}]$$

$${\Omega }_{\theta }^{\mathrm{FA}}=\frac{{\alpha }_{\mathrm{ee}}vb}{\pi }\cos (\theta ),$$

$$\mathrm{Re}[{\chi }_{\mathrm{BB}}^{\mathrm{eff}}({\mathbf{q}}_{\parallel },\omega )]\to \frac{n{q}_{\parallel }}{m_{\mathrm{eff}}{\omega }^2},$$

$${\rm Re}[\chi^{\rm eff}_{\rm AB}({\bm{q}}_{\parallel},\omega)]\to-\frac{q^{2}_{\parallel}\big{[}1+\sin^{2}(\theta)\big{]}}{3\pi^{2}}\frac{1}{\sqrt{2E_{\rm F}|\hbar\omega|}}~{}.$$

$$\Omega_{\theta}=\frac{1}{2}\Big{[}\Omega_{\theta}^{\rm FA}+\sqrt{(\Omega^{\rm FA}_{\theta})^{2}+4\Omega^{2}_{\rm s}}\Big{]}$$

$$\delta \Omega ({\mathbf{q}}_{\parallel })={\alpha }_{\mathrm{ee}}v{q}_{\parallel }\mathcal{I}(\theta ).$$

$$\mathcal{I}(\theta )$$

$$\Omega ({\mathbf{q}}_{\parallel })={\Omega }_{\theta }^{\mathrm{FA}}+q_{\Vert }\left[v\cos (\theta )-|\sin (\theta )|\frac{{\Omega }_{\theta }^{\mathrm{FA}}}{2b}\right].$$

$$\Gamma ({\mathbf{q}}_{\parallel })=\frac{\mathrm{Im}[{\chi }_{\mathrm{eff}}({\mathbf{q}}_{\parallel },\Omega ({\mathbf{q}}_{\parallel })]}{{\partial }_{\omega }\mathrm{Re}[{\chi }_{\mathrm{eff}}({\mathbf{q}}_{\parallel },\omega )]{|}_{\omega =\Omega ({\mathbf{q}}_{\parallel })}}.$$

$$\begin{split}{\rm Im}[\chi^{\rm eff}_{\rm BB}({\bm{q}}_{\parallel},\omega)]\simeq\int\frac{dq_{z}}{2\pi}\frac{2q_{\parallel}^{2}}{\big{(}q_{\parallel}^{2}+q_{z}^{2}\big{)}^{2}}{\rm Im}[\chi^{({\rm h})}_{\rm BB}(\bm{q}_{\rm 3D},\omega)]~{}.\end{split}$$

$$\begin{array}{r}\mathrm{Im}[{\chi }_{\mathrm{BB}}^{\mathrm{eff}}({\mathbf{q}}_{\parallel },\omega )]\to -\frac{q_{\Vert }^2}{32{\pi }^2\hslash \omega }{\left(\frac{2E_{\mathrm{F}}}{\hslash \omega }\right)}^2.\end{array}$$

$$\begin{split}{\rm Im}[\chi^{\rm eff}_{\rm AB}({\bm{q}}_{\parallel},\omega)]\to-\frac{q_{\parallel}^{2}\big{[}1+\sin^{2}(\theta)\big{]}}{3\pi^{2}\hbar\omega}\left(\frac{2E_{\rm F}}{\hbar\omega}\right)^{-1/2}~{}.\end{split}$$

$$\begin{array}{r}\Gamma ({\mathbf{q}}_{\parallel })=\frac{{\Omega }_{\theta }}{32}\frac{q_{\parallel }}{k_{\mathrm{F}}}\frac{1}{\frac{\cos (\theta )}{2}\frac{b}{k_{\mathrm{F}}}+\frac{1}{6}\left(\frac{2E_{\mathrm{F}}}{\hslash {\Omega }_{\theta }}\right)}{\left(\frac{2E_{\mathrm{F}}}{\hslash {\Omega }_{\theta }}\right)}^2.\end{array}$$

$$v_{\mathrm{ext}}(z,{\mathbf{q}}_{\parallel },t)=V_0{e}^{-q_{\parallel }|z|}{e}^{-i{\omega }_{0}t},$$

$$V_{\mathrm{sc}}(0,{\mathbf{q}}_{\parallel },{\omega }_0)=\frac{V_0}{1-\frac{2\pi e^2}{q_{\parallel }}{\chi }_{\mathrm{eff}}({\mathbf{q}}_{\parallel },{\omega }_0)}.$$

$$V_{\rm sc}(0,{\bm{q}}_{\parallel},\omega_{0})\simeq\frac{V_{0}}{\displaystyle-\frac{2\pi e^{2}}{q_{\parallel}}\left(\frac{\partial\chi_{\rm eff}({\bm{q}}_{\parallel},\omega)}{\partial\omega}\right)_{\omega=\Omega({\bm{q}}_{\parallel})}\big{[}\omega_{0}-\Omega({\bm{q}}_{\parallel})+i\Gamma({\bm{q}}_{\parallel})\big{]}}~{}.$$

$$V_{\rm sc}(0,{\bm{r}}_{\parallel},\omega_{0})\simeq V_{0}\int\frac{d^{2}{\bm{q}}}{(2\pi)^{2}}\frac{e^{i{\bm{q}}_{\parallel}{\bm{r}}_{\parallel}}}{\displaystyle-\frac{2\pi e^{2}}{q_{\parallel}}\left(\frac{\partial\chi_{\rm eff}({\bm{q}}_{\parallel},\omega)}{\partial\omega}\right)_{\omega=\Omega({\bm{q}}_{\parallel})}\big{[}\omega_{0}-\Omega({\bm{q}}_{\parallel})+i\Gamma({\bm{q}}_{\parallel})\big{]}}~{}.$$

$$V_{\rm sc}(0,{\bm{r}}_{\parallel},\omega_{\rm TIP})\simeq V_{0}\int_{0}^{\infty}q_{\parallel}dq_{\parallel}\int_{0}^{2\pi}d\theta\frac{1}{(2\pi)^{2}}\frac{e^{iq_{\parallel}r_{\parallel}\cos(\theta-\theta^{\prime})}}{\displaystyle-\frac{2\pi e^{2}}{q_{\parallel}}\left(\frac{\partial\chi_{\rm eff}(q_{\parallel},\theta,\omega)}{\partial\omega}\right)_{\omega=\Omega({\bm{q}}_{\parallel})}\big{[}\omega_{0}-\Omega(q_{\parallel},\theta)+i\Gamma(q_{\parallel},\theta)\big{]}}~{},$$

$$p_1=\frac{{\omega }_0-{\Omega }_{\theta }}{{\Gamma }_{\theta }^2+[v\alpha \mathcal{I}(\theta ){]}^2}v\alpha \mathcal{I}(\theta )$$

Full text

I Introduction

Plasmons are self-sustained oscillations of the charge density that occur in metals and doped semiconductors Pines_and_Nozieres ; Platzman_and_Wolff ; Giuliani_and_Vignale . The coupling of these matter excitations to photons enables to squeeze electromagnetic radiation from the visible Maier07 to the Terahertz (THz) pablo_naturenano_2017 spectral range into nanoscale devices. However, plasmons suffer scattering from ever-present extrinsic mechanisms in real solid-state devices khurgin_naturenano_2015 , including phonons and disorder. Therefore, when confinement is significant, i.e. when the plasmon wavelength $\lambda_{\rm p}$ in the material is much smaller than the illumination wavelength $\lambda_{0}=2\pi c/\omega$, losses tend to be high (at room temperature) and hamper potential technological breakthroughs.

Substantial efforts have been recently made to increase the lifetime of these propagating modes at room temperature, without decreasing the associated compression ratio $\lambda_{\rm p}/\lambda_{0}$. For example, one can utilize high-quality graphene sheets encapsulated in hexagonal boron nitride geim_nature_2013 , where graphene plasmons scatter essentially only against the acoustic phonons of the two-dimensional (2D) carbon lattice woessner_naturemater_2015 , which are weakly coupled to the electronic degrees of freedom. Another possible pathway is to use plasmons in topologically-non-trivial materials song_pnas_2016 ; kumar_prb_2016 ; dipietro_naturenano_2013 ; autore_aom_2015 . In the particular case of crystals displaying broken time-reversal symmetry (BTRS), the existence of unidirectional propagating modes akin to the ultra-long-lived kumada_prl_2014 topological jin_naturecommun_2016 edge magnetoplasmons that occur in 2D electron systems in the quantum Hall regime Giuliani_and_Vignale is expected. Technologically, it would be extremely useful to use materials where BTRS occurs without the aid of an external magnetic field. Natural candidates among topological materials with BTRS are recently discovered Weyl semimetals (WSMs) Hosur_2013 ; Gorbar_2016 ; Hasan_2017 ; Yan_2017 ; Burkov_2017 ; Armitage2018 . These are semimetals with protected linear band crossings in the Brillouin zone, which act as power-law-decaying sources of Berry curvature Hosur_2013 ; Gorbar_2016 ; Hasan_2017 ; Yan_2017 ; Burkov_2017 ; Armitage2018 . Some of these compounds do display intrinsic BTRS Shekhar_arXiv_2016 and, at the same time, have intriguing topological surface states called “Fermi arcs” (FAs) Hosur_2013 ; Gorbar_2016 ; Hasan_2017 ; Yan_2017 ; Burkov_2017 ; Armitage2018 . These are disjoint segments of a 2D Fermi contour—see Fig. 1(b)—which have so far been imaged only with angle-resolved photoemission spectroscopy (ARPES) xu_science_2015 ; lv_prx_2015 . Provided that losses are not too strong and that the plasmon electric field leaks sufficiently outside the material, the plasmonic excitations of FA surface states can in principle be explored with spatial resolution using scanning-type near-field optical spectroscopy (for a recent review see e.g. Ref. basov_science_2016, ).

Theoretically, a few aspects of plasmons in WSMs have been studied for both cases of bulk pellegrino_prb_2015 and surface hofmann_prb_2016 modes. In particular, the authors of Ref. pellegrino_prb_2015, carried out a study of bulk plasmons in WSMs, by coupling the Maxwell equations as modified by the axion term wilczek_prl_1987 with a local approximation for the constitutive equation, $J_{\alpha}({\bm{q}},\omega)\approx\sum_{\beta}\sigma_{\alpha\beta}({\bm{0}},\omega)E_{\beta}({\bm{q}},\omega)$. The local conductivity tensor $\sigma_{\alpha\beta}({\bm{0}},\omega)$ was calculated pellegrino_prb_2015 by using semiclassical Boltzmann transport theory in the clean limit. The authors of Ref. hofmann_prb_2016, studied WSM surface plasmon polaritons by using the same approach and deep in the retarded limit, i.e. for surface plasmon wave numbers $|{\bm{q}}_{\|}|\sim\omega/c$. Finally, the authors of Ref. song_arxiv_2017, carried out a semiclassical study of WSM surface plasmons, which is intrinsically valid in the long-wavelength limit. Non-local corrections were heuristically added to the semiclassical equations of motion. An important open question on the plasmon lifetime remains. In order to quantify this, the quantum-mechanical coupling between Fermi arc surface states and bulk states needs to be evaluated.

In this work, we present a fully quantum-mechanical theory of WSM Fermi arc (FA) plasmons that goes beyond the state-of-the-art pellegrino_prb_2015 ; hofmann_prb_2016 ; song_arxiv_2017 . The present derivation focuses on the simplest microscopic model Hamiltonian of a (type-I) WSM with BTRS Hosur_2013 ; Gorbar_2016 ; Hasan_2017 ; Yan_2017 ; Burkov_2017 and is based on linear response theory Pines_and_Nozieres ; Platzman_and_Wolff ; Giuliani_and_Vignale and the random phase approximation (RPA) Pines_and_Nozieres ; Platzman_and_Wolff ; Giuliani_and_Vignale . We focus on the electrostatic regime, where the plasmon wave number $|{\bm{q}}_{\|}|$ is much larger than the photon one ($\omega/c$), enabling a great concentration of electromagnetic energy. We claim that quantum non-local effects are crucial to understand WSM FA plasmon physics. First of all, our theory captures analytically the first non-local correction to the FA plasmon dispersion in the limit $|{\bm{q}}_{\|}|\ll k_{\rm F}$, where $k_{\rm F}$ is the bulk Fermi wave number. Second, since the FA wavefunctions are in contact (i.e. strong spatial overlap) with a bulk of gapless excitations, FA plasmons are susceptible to Landau damping even at zero temperature and deep in the long-wavelength $|{\bm{q}}_{\|}|\ll k_{\rm F}$ limit. We quantify this intrinsic dissipation mechanism, which is dominated by processes whereby FA plasmons decay by emitting electron-hole pairs in the bulk. Our calculations on the FA plasmon lifetime pose strict bounds on the observability of certain angular portions of the highly-anisotropic FA plasmon dispersion. Indeed, since FA plasmon modes can be strongly overdamped along certain angular directions, it is essential to consider the dispersion relation as well as the angular dependence of the losses. Finally, we study the pattern of FA plasmon waves that can be measured by carrying out a scattering-type near-field optical experiment (s-SNOM) basov_science_2016 on the surface of a WSM with BTRS. In this class of experiments, light is focused on a metallized tip with the aim of launching propagating surface plasmons. We calculate the screened potential associated with FA plasmons, which clearly shows characteristic features due to their peculiar propagation dynamics.

Our Article is organized as following. In Sect. II we present the single-particle Hamiltonian and corresponding bulk and surface eigenvalues and eigenstates of a type-I semi-infinite WSM with BTRS. In Sect. III we lay down a general theory of surface plasmons hosted by electron systems occupying a semi-infinite space. In Sect. IV we present a formal expression of the linear density-density response function of a WSM based on the Lehmann representation. We also introduce analytical asymptotic formulae for the same quantity in the relevant limits with respect to the 2D conserved wave vector ${\bm{q}}_{\|}$ and frequency $\omega$. Sect. V is devoted to the FA plasmon dispersion, while intrinsic damping of these modes is discussed in Sect. VI. Finally, in Sect. VII we present a calculation of the spatial pattern of FA plasmon waves that would be seen in a s-SNOM experiment on the surface of a WSM with BTRS. A brief summary is reported in Sect. VIII. Five appendices report a number of very useful technical details.

II Single-particle physics of semi-infinite WSMs

In order to describe a three-dimensional (3D) WSM with BTRS and two Weyl nodes separated by a vector $2{\bm{b}}$ (see Fig. 1), we use the following family of Hamiltonians Murakami depending on the parameter $m$:

[TABLE]

In Eq. (1), $v$ is the Weyl fermion’s velocity and $\sigma_{i}$ are ordinary $2\times 2$ Pauli matrices. The first term on the right-hand side of Eq. (1) breaks time-reversal symmetry. For $m=b^{2}>0$, the low-energy spectrum of ${\cal H}_{m}({\bm{k}})$ contains two Weyl nodes separated by $2{\bm{b}}=2b\hat{\bm{k}}_{x}$ with $b>0$. The spectrum of ${\cal H}_{m}({\bm{k}})$ is $E_{\lambda}({\bm{k}})=\lambda\hbar v\sqrt{K^{2}_{x}+k^{2}_{y}+k^{2}_{z}}$, where $\lambda=\pm 1$ denotes conduction/valance band states and $K_{x}\equiv(k_{x}^{2}-b^{2})/(2b)$. For $k_{x}\approx\pm b$, $E_{\lambda}({\bm{k}})$ is linear and isotropic, with two Weyl cones located at ${\bm{k}}=\pm b\hat{{\bm{k}}}_{x}$. The corresponding eigenstates are $\Psi_{\lambda,{\bm{k}}}({\bm{r}})=V^{-1/2}u_{{\bm{k}},\lambda}e^{i{\bm{k}}\cdot{\bm{r}}}$, where $V=SL_{z}$ is the 3D electron system volume, $S$ the surface area, $L_{z}$ its thickness, and $u_{{\bm{k}},\lambda}$ is a spinor given by:

[TABLE]

with $\cos{(\beta_{\bm{k}})}=k_{y}/\sqrt{K^{2}_{x}+k^{2}_{y}+k^{2}_{z}}$, $\sin{(\beta_{\bm{k}})}=\sqrt{(K^{2}_{x}+k^{2}_{z})/(K^{2}_{x}+k^{2}_{y}+k^{2}_{z})}$, and $e^{i\phi_{\bm{k}}}=(K_{x}+ik_{z})/\sqrt{K^{2}_{x}+k^{2}_{z}}$.

In this work, we are interested in the surface states of a 3D WSM, which we obtain from the following procedure. We first note that, for $m<0$, ${\cal H}_{m}({\bm{k}})$ describes an insulator with a gap $E_{\rm g}=-\hbar vm/(2b)$. The vacuum can be therefore modelled by taking the $\lim_{m\rightarrow-\infty}{\cal H}_{m}({\bm{k}})$. Fermi arcs (FAs) are present only on a surface parallel to Rao ; Yan_2017 $2{\bm{b}}=2b\hat{\bm{k}}_{x}$. Orienting the crystal in such a way that $(\hat{\bm{k}}_{x},\hat{\bm{k}}_{y},\hat{\bm{k}}_{z})\|(\hat{\bm{x}},\hat{\bm{y}},\hat{\bm{z}})$, we locate the vacuum/WSM interface hosting the FAs at $z=0$. Surface states emerge by assuming that that the parameter $m$ in Eq. (1) changes sign with $z$, i.e.  $m(z)=b^{2}\Theta(-z)-\tilde{m}\Theta(z)$, where $\Theta(z)$ is the Heaviside step function and $\tilde{m}$ is a positive constant. This choice describes a WSM for $z<0$ and the vacuum for $z>0$, provided that one takes the limit $\tilde{m}\rightarrow+\infty$. Breaking translational invariance along the $\hat{\bm{z}}$ direction requires the change $k_{z}\to-i\partial_{z}$. FAs are described by states that are bound to the surface. Following Ref. Murakami, , we find that their wave functions are

[TABLE]

where ${\bm{k}}_{\parallel}=(k_{x},k_{y})$, $\ell=2b/({b}^{2}-k_{x}^{2})$, and the index “A” in Eq. (3) stands for “arc”. These states have chiral dispersion, $E_{\rm A}({\bm{k}}_{\|})=\hbar vk_{y}$, group velocity given by $v$, density of states given by ${\cal N}_{\rm A}=b/(2\pi^{2}\hbar v)$ (which is independent of the bulk doping), and exist only between the two Weyl nodes, i.e. for $-b\leq k_{x}\leq b$. Note that $\ell$, which physically represents the extension of the FA state into the WSM bulk $z<0$, goes to infinity for $k_{x}\to\pm b$. This means that for $k_{x}$ near the location of the Weyl nodes, leakage of the FA states into the bulk cannot be neglected. For the sake of simplicity, the FA states (3) will be denoted by the shorthand $\ket{\rm A}$.

The presence of the WSM/vacuum interface at $z=0$ affects also the propagating (i.e. bulk) states. Since $\bm{k}_{\|}$ is a good quantum number (because of translational invariance in the $\hat{\bm{x}}$-$\hat{\bm{y}}$ plane), we obtain propagating states in the presence of the interface by taking linear combinations of bulk plane-wave states with positive ($k_{z}>0$) and negative ($-k_{z}<0$) values of the $\hat{\bm{z}}$ component of the wave vector ${\bm{k}}$. Assuming specular reflection at the interface, we find

[TABLE]

where ${\bm{\bar{k}}}=({\bm{k}}_{\parallel},-k_{z})$ is the reflected wave vector and $r_{\bm{k}}$ is a reflection coefficient (with $|r_{\bm{k}}|=1$) to be determined by imposing suitable boundary conditions. Of course, $\Psi^{\rm B}_{\lambda,{\bm{k}}}({\bm{r}}_{\|},z)$ is an eigenstate of the half-space Hamiltonian $\lim_{\tilde{m}\to+\infty}{\cal H}_{m(z)}({\bm{k}}_{\|},k_{z}\to-i\partial_{z})$ with eigenvalue $E_{\lambda}({\bm{k}})$. As explained in Appendix A , the continuity of the wave functions at the interface requires $r_{\bm{k}}=-1$. Because the “reflected” states (4) originate from the WSM bulk states $\Psi_{\lambda,{\bm{k}}}({\bm{r}})$, we will denote them by the shorthand $\ket{\rm B}$.

III Surface plasmons of semi-infinite WSMs

Due to the presence of the interface at $z=0$, it is natural to lay down linear response theory Giuliani_and_Vignale by imposing translational invariance in the direction perpendicular to $\hat{\bm{z}}$ and that the wave vector ${\bm{q}}_{\parallel}$ parallel to the surface is a good quantum number. The linear density response $n_{1}$ induced by an external potential is therefore a function of $z$, ${\bm{q}}_{\parallel}$, and frequency $\omega$. It can be expressed in terms of the screened potential Giuliani_and_Vignale , $V_{\rm sc}(z,{\bm{q}}_{\parallel},\omega)$, and the proper density-density response function Giuliani_and_Vignale $\tilde{\chi}_{nn}(z,z^{\prime},{\bm{q}}_{\parallel},\omega)$ according to

[TABLE]

For a self-sustained oscillation occurring in the absence of an external potential, the screened potential is related to the induced density by

[TABLE]

where the integrals over $z^{\prime}$ and $z^{\prime\prime}$ span the half-space $z<0$ and $v(z,z^{\prime},{\bm{q}}_{\parallel})=2\pi e^{2}\exp{(-q_{\parallel}|z-z^{\prime}|)}/q_{\parallel}$ is the Fourier transform of the Coulomb potential ($q_{\parallel}=|{\bm{q}}_{\parallel}|$).

Plasmons are non-trivial solutions of the integral equation (6). In this work we are solely interested in plasmons whose associated electric field is bound to the WSM/vacuum interface. We therefore set $V_{\rm sc}(z,{\bm{q}}_{\parallel},\omega)={\bar{v}}_{\rm sc}({\bm{q}}_{\parallel},\omega)e^{-q_{\parallel}|z|}$, which describes a mode bound to the surface with a “localization” length scale equal to $1/q_{\parallel}$. Utilizing this Ansatz for $V_{\rm sc}(z,{\bm{q}}_{\parallel},\omega)$, we evaluate pitarke_rpp_2007 Eq. (6) for $z>0$. We find that this integral equation can be written as a standard algebraic 2D plasmon equation,

[TABLE]

provided that one introduces the effective 2D proper response function:

[TABLE]

Here,

[TABLE]

No approximation has yet been made on $\tilde{\chi}_{nn}(q_{z},q_{z}^{\prime},{\bm{q}}_{\parallel},\omega)$. The only key assumption we used to derive Eq. (7) was that $\tilde{\chi}_{nn}(z,z^{\prime},{\bm{q}}_{\parallel},\omega)\propto\Theta(-z)\Theta(-z^{\prime})$.

Screening due to the dielectric background can be taken care of by replacing $e^{2}\to e^{2}/\bar{\epsilon}$ in Eq. (7), where $\bar{\epsilon}\equiv(1+\epsilon_{\rm b})/2$ and $\epsilon_{\rm b}$ is the high-frequency WSM bulk dielectric constant. This prescription, which is valid only for surface modes, is demonstrated in Appendix B.

IV Non-interacting response function of a semi-infinite WSM

In the RPA Giuliani_and_Vignale , the exact proper response function $\tilde{\chi}_{nn}(q_{z},q_{z}^{\prime},{\bm{q}}_{\parallel},\omega)$ in Eq. (9) is replaced by the non-interacting response function ${\chi}^{(0)}_{nn}(q_{z},q_{z}^{\prime},{\bm{q}}_{\parallel},\omega)$, which can be evaluated with the aid of the Lehmann representation Giuliani_and_Vignale :

[TABLE]

Here, $\ket{n}$, $\ket{m}$ ($E_{n}$, $E_{m}$) denote the eigenstates (eigenenergies) of the non-interacting Hamiltonian, $\hbar\omega_{n,m}=E_{n}-E_{m}$ are the excitation energies, $f_{n}$ are the occupation numbers, and $\hat{n}_{{\bm{q}}_{\parallel},q_{z}}\equiv e^{-i{\bm{q}}_{\|}\cdot\hat{{\bm{r}}}_{\parallel}-iq_{z}\hat{\bm{z}}}$. As we have seen above, the states of our semi-infinite WSM comprise the FA states (3) and the reflected states (4). We can therefore naturally decompose $\chi^{(0)}_{nn}(q_{z},q_{z}^{\prime},{\bm{q}}_{\parallel},\omega)$ into the sum of three contributions:

[TABLE]

In the first contribution, ${\chi}^{(0)}_{\rm AA}$, the sum over $n$ and $m$ in Eq. (10) spans only the FA states (3), i.e. $\ket{n}=\ket{A}$, $\ket{m}=\ket{A}^{\prime}$. In the second contribution, ${\chi}^{(0)}_{\rm BB}$, the sum spans only the reflected states (4), i.e. $\ket{n}=\ket{B}$, $\ket{m}=\ket{B}^{\prime}$. Finally, the third contribution, ${\chi}^{(0)}_{\rm AB}$, takes into account the remaining cross processes in which $\ket{n}=\ket{A}$ and $\ket{m}=\ket{B}$ or $\ket{n}=\ket{B}$ and $\ket{m}=\ket{A}$. Because of Eq. (8), the same decomposition holds true for $\chi_{\rm eff}({\bm{q}}_{\parallel},\omega)$, i.e. $\chi_{\rm eff}=\chi^{\rm eff}_{\rm AA}+\chi^{\rm eff}_{\rm BB}+\chi^{\rm eff}_{\rm AB}$. The FA plasmon dispersion is controlled by ${\rm Re}[\chi_{\rm eff}({\bm{q}}_{\parallel},\omega)]$, while the plasmon damping rate depends on dissipation, i.e. ${\rm Im}[\chi_{\rm eff}({\bm{q}}_{\parallel},\omega)]$.

To make analytical progress, we concentrate our attention on FA plasmons in the long-wavelength limit, defined by the regime in which $q_{\parallel}$ is the smallest wave vector scale in the problem. Eq. (8) implies $q_{z}\sim q_{z}^{\prime}\sim q_{\parallel}$. We denote by the shorthand $q$ the small quantity $q_{z}\sim q_{z}^{\prime}\sim q_{\parallel}$.

In the limit $q\ll b,\omega/v$, we find that the leading order FA effective response function reduces to (Appendix C)

[TABLE]

Like the FA density-of-states ${\cal N}_{\rm A}$, ${\rm Re}[\chi^{\rm eff}_{\rm AA}]$ does not depend on the bulk carrier density $n$. Neglecting for a moment the contribution to Eq. (7) from the reflected states, ${\rm Re}[\chi^{\rm eff}_{\rm BB}]$, and the cross term, ${\rm Re}[\chi^{\rm eff}_{\rm AB}]$, we find a long-wavelength FA plasmon dispersion

[TABLE]

where $\theta$ is the angle between ${\bm{q}}_{\|}$ and the $\hat{\bm{y}}$ axis (i.e. the direction along which the FA states have finite group velocity), $\cos(\theta)=q_{y}/q_{\|}$, and $\alpha_{\rm ee}=e^{2}/(\hbar v\bar{\epsilon})$ is a dimensionless coupling constant describing the strength of electron-electron interactions. The strong directionality dependence of $\Omega^{\rm FA}_{\theta}$ is evident. Interestingly, the group velocity ${\bm{v}}^{\rm FA}_{\theta}=\hat{\theta}q^{-1}_{\|}\partial_{\theta}\Omega^{\rm FA}_{\theta}=-\alpha_{\rm ee}vb\sin(\theta)\hat{\theta}/(\pi q_{\|})$ of the mode in Eq. (13) vanishes when the plasmon wave vector ${\bm{q}}_{\|}$ is in the ${\hat{\bm{y}}}$ direction.

We now turn to evaluate the contribution ${\chi}^{(0)}_{\rm BB}(q_{z},q_{z}^{\prime},{\bm{q}}_{\parallel},\omega)$ due to the reflected states (4). For a Fermi wave number $k_{\rm F}\ll b$, we can linearize the spectrum $E_{\pm}({\bm{k}})$ around the points $k^{(j)}_{x}=jb$, where $j=\pm 1$ is the Weyl cone index. Again, we work in the long-wavelength limit, which in this case is defined by $q_{\parallel}\ll\omega/v,k_{\rm F}$. Noting that $q_{z},q_{z}^{\prime}\sim q_{\parallel}$, we can state that $q_{\parallel}$ is the smallest wave vector in the problem. Moreover, we assume to be deep in the single-particle optical absorption gap, i.e. we consider $\hbar\omega\ll 2E_{\rm F}$ where $E_{\rm F}=\hbar vk_{\rm F}$ is the Fermi energy. In this regime, inter-band contributions to all response functions involving the reflected states (4) can be neglected and we can safely use the high-frequency moment expansion Giuliani_and_Vignale ; torre_prb_2017 for the intra-band contribution to ${\rm Re}[{\chi}^{(0)}_{\rm BB}]$. As detailed in Appendix D, we find

[TABLE]

where $n=k^{3}_{\rm F}/(6\pi^{2})$ is the density of bulk electrons and $m_{\rm eff}=E_{\rm F}/v^{2}$ the effective mass. Corrections to Eq. (14) scale as $q^{3}_{\|}$, while terms $\propto q^{2}_{\|}$ are exactly zero. Notice that inserting only ${\rm Re}[\chi^{\rm eff}_{\rm BB}({\bm{q}}_{\parallel},\omega)]$ in Eq. (7), we recover the well-known Ritchie relation ritchie_pr_1957 $\Omega_{\rm s}=\Omega_{\rm b}\sqrt{\epsilon_{\rm b}/(2\bar{\epsilon})}=vk_{\rm F}\sqrt{\alpha_{\rm ee}/(3\pi)}$ for the surface plasmon frequency, $\Omega^{2}_{\rm b}=4\pi ne^{2}/(\epsilon_{\rm b}m_{\rm eff})$ being the square of the bulk plasmon frequency. In the undoped $k_{\rm F}=0$ limit, Eq. (14) yields ${\rm Re}[\chi^{\rm eff}_{\rm BB}({\bm{q}}_{\parallel},\omega)]=0$.

Finally, we need to evaluate the cross contributions ${\chi}^{(0)}_{\rm AB}(q_{z},q_{z}^{\prime},{\bm{q}}_{\parallel},\omega)$ and ${\rm Re}[\chi^{\rm eff}_{\rm AB}({\bm{q}}_{\parallel},\omega)]$. Up to leading order in the long-wavelength limit and, as above, deep in the single-particle optical absorption gap $\hbar\omega/(2E_{\rm F})\ll 1$, we find (Appendix E)

[TABLE]

Clearly Eq. (15) cannot be used in the undoped $k_{\rm F}\to 0$ limit. In this case an explicit calculation ( Appendix E) yields $\lim_{k_{\rm F}\to 0}{\rm Re}[\chi^{\rm eff}_{\rm AB}({\bm{q}}_{\parallel},\omega)]=0$. Note that $\chi^{\rm eff}_{\rm AB}({\bm{q}}_{\parallel},\omega)$ does not depend on $b$ to leading order for $q_{\|}\ll b$. In this limit, indeed, $\chi^{\rm eff}_{\rm AB}({\bm{q}}_{\parallel},\omega)$ describes electron-hole transitions near each Weyl node.

V FA plasmon dispersion

We now solve Eq. (7) in the RPA, i.e. by replacing ${\rm Re}[\chi_{\rm eff}({\bm{q}}_{\parallel},\omega)]$ with ${\rm Re}[\chi^{\rm eff}_{\rm AA}({\bm{q}}_{\parallel},\omega)+\chi^{\rm eff}_{\rm BB}({\bm{q}}_{\parallel},\omega)+\chi^{\rm eff}_{\rm AB}({\bm{q}}_{\parallel},\omega)]$. In the following, we focus only on positive frequencies, since modes with negative frequency can be obtained from the relation $\Omega({\bm{q}}_{\parallel})=-\Omega(-{\bm{q}}_{\parallel})$.

For doped WSMs ($k_{\rm F}\neq 0$) we find that the plasmon dispersion in the long-wavelength limit is given by $\Omega({\bm{q}}_{\parallel})=\Omega_{\theta}+\delta\Omega({\bm{q}}_{\parallel})$, where $\Omega_{\theta}$ is the value of the plasmon frequency for $q_{\parallel}\rightarrow 0$ and $\delta\Omega({\bm{q}}_{\parallel})$ is the first-order non-local correction in $q_{\parallel}$:

[TABLE]

and

[TABLE]

In Eq. (17) we have introduced the quantity

[TABLE]

The first line in the previous equation is due to ${\rm Re}[\chi^{\rm eff}_{\rm AA}({\bm{q}}_{\parallel},\omega)]$, while the second line is due to ${\rm Re}[\chi^{\rm eff}_{\rm AB}({\bm{q}}_{\parallel},\omega)]$. Eq. (16) is formally identical hofmann_prb_2016 to the dispersion relation of a surface magnetoplasmon wallis_prb_1974 in the case of a magnetic field oriented along ${\bm{b}}$ and with $\Omega^{\rm FA}_{\theta}\to\omega_{\rm c}\cos(\theta)$, $\omega_{\rm c}$ being the cyclotron frequency. This result reflects the fact that the topological FA states play the same role of edge states in the quantum Hall regime. The FA plasmon dispersion relation $\Omega({\bm{q}}_{\|})$ is illustrated in Figs. 2(a) and 3(a). In Fig. 2(a) we clearly note that, for every value of $q_{x}$, there is a large asymmetry in the dispersion $\Omega({\bm{q}}_{\|})$ as a function of $q_{y}$. This stems from the fact that, for $q_{y}>0$, $\Omega({\bm{q}}_{\|})$ is dominated by the FA states. On the other hand, for $q_{y}<0$, the main contribution to $\Omega({\bm{q}}_{\|})$ is given by $\Omega_{\rm s}$, which stems from the reflected states (4) and is weakly dispersive. The weak dispersive features seen for $q_{y}<0$ reflect an hybridization between the FA and bulk channels. This mixing of channels with very different symmetries (the former highly directional, while the latter is isotropic) leads to a pair of saddle points in the dispersion relation, which are seen from the contours in Fig. 2(a) at $\hbar\Omega({\bm{q}}_{\|})/E_{\rm F}\approx 0.20$. In Fig. 3(a) we show cuts of $\Omega({\bm{q}}_{\|})$, each one taken at a fixed value of $q_{\|}$, and plotted as a function of $-\pi\leq\theta\leq\pi$.

For typical experimental values Armitage2018 , $\epsilon_{\rm b}=10$, $E_{\rm F}=40~{}{\rm meV}$, $b=0.05~{}{\rm{\AA}}^{-1}$, and $v=c/1000$, Eq. (16) yields a FA plasmon energy $\Omega_{\theta}$ in the forward $\theta\approx 0$ direction on the order of $30~{}{\rm meV}$, corresponding to a frequency $\nu_{\theta}\equiv\Omega_{\theta}/(2\pi)$ of $\sim 7.3~{}{\rm THz}$.

In the undoped $k_{\rm F}\to 0$ limit the plasmon dispersion is

[TABLE]

We will momentarily show that in the undoped limit the FA plasmon is strongly damped, but for non-zero $k_{\rm F}$, this damping is suppressed.

VI Intrinsic lifetime of FA plasmons

So far, we have treated the plasmon as a solution of Eq. (6) occurring on the real frequency axis. For fundamental reasons related to causality Giuliani_and_Vignale , however, retarded response functions must have poles located below the real axis, i.e. at $\omega=\Omega({\bm{q}}_{\parallel})-i\Gamma({\bm{q}}_{\parallel})$ with $\Gamma({\bm{q}}_{\parallel})>0$. When the imaginary part is small, i.e. when $\Gamma({\bm{q}}_{\parallel})\ll\Omega({\bm{q}}_{\parallel})$, the plasmon is a well-defined collective excitation of the many-body system Giuliani_and_Vignale . The plasmon lifetime $\tau({\bm{q}}_{\|})$ is $\Gamma^{-1}({\bm{q}}_{\|})$. Expanding Eq. (7) in the small parameter $\Gamma({\bm{q}}_{\parallel})/\Omega({\bm{q}}_{\parallel})$ we obtain:

[TABLE]

We define the quality factor as $Q({\bm{q}})\equiv\Omega({\bm{q}}_{\parallel})/[2\Gamma({\bm{q}}_{\parallel})]$. We first calculate $\partial_{\omega}{\rm Re}[\chi_{\rm eff}({\bm{q}}_{\parallel},\omega)|_{\omega=\Omega({\bm{q}}_{\parallel})}]$ at leading order in $q_{\parallel}$. Only ${\chi}^{(0)}_{\rm AA}$ and ${\chi}^{(0)}_{\rm BB}$ contribute to this quantity at leading order in $q_{\parallel}$. We find $\partial_{\omega}{\rm Re}[\chi_{\rm eff}({\bm{q}}_{\parallel},\omega)|_{\omega=\Omega({\bm{q}}_{\parallel})}]\to-q_{\parallel}[\cos(\theta)b/2+vk^{2}_{\rm F}/(3\Omega_{\theta})]/(\hbar\Omega^{2}_{\theta}\pi^{2})$. We then calculate ${\rm Im}[\chi_{\rm eff}({\bm{q}}_{\parallel},\omega)]$. As detailed in Appendix C, ${\rm Im}[\chi^{\rm eff}_{\rm AA}({\bm{q}}_{\parallel},\omega)]\to-q_{\parallel}b\cos(\theta)\delta(\omega-E_{\rm A}({\bm{q}}_{\parallel})/\hbar)/(2\hbar\pi)$. This result can be easily understood as the absorption spectrum of a one-dimensional system Giuliani_and_Vignale , whose role in our case is effectively played by the FA states. For the calculation of ${\rm Im}[\chi^{\rm eff}_{\rm BB}({\bm{q}}_{\parallel},\omega)]$, we use the well-known approximate relation Garrido $\chi_{\rm BB}^{\rm eff}(z,z^{\prime},{\bm{q}}_{\parallel},\omega)\simeq\Theta(-z)\Theta(-z^{\prime})[\chi^{({\rm h})}_{\rm BB}(z+z^{\prime},{\bm{q}}_{\parallel},\omega)+\chi^{({\rm h})}_{\rm BB}(z-z^{\prime},{\bm{q}}_{\parallel},\omega)]$. Here, $\chi^{({\rm h})}_{\rm BB}(z,{\bm{q}}_{\parallel},\omega)=\int dq_{z}\chi^{({\rm h})}_{\rm BB}({\bm{q}}_{\rm 3D},\omega)\exp(iq_{z}z)/(2\pi)$ is the Fourier transform of the homogenous response function $\chi^{({\rm h})}_{\rm BB}({\bm{q}}_{\rm 3D},\omega)$ of a bulk WSM lv_ijmp_2013 ; zhou_prb_2015 and ${\bm{q}}_{\rm 3D}=(q_{z},{\bm{q}}_{\|})$. Straightforward algebraic manipulations yield

[TABLE]

Eq. (21) is crucial as it states that, for every $\omega$, the decay rate (20) of the FA plasmon depends on the spectral density of electron-hole pairs ${\rm Im}[\chi^{({\rm h})}_{\rm BB}(\bm{q}_{\rm 3D},\omega)]$ in the bulk integrated over all values of $q_{z}$. Physically, this equation expresses the fact that a surface plasmon can decay without conserving the $\hat{\bm{z}}$ component of the wave vector $\bm{q}_{\rm 3D}$, for the presence of the surface breaks translational invariance in this direction relaxing the three-dimensional momentum conservation that a bulk plasmon would obey. The latter is encoded in regions of the plane ${\bm{q}}_{\rm 3D}$-$\omega$ where ${\rm Im}[\chi^{({\rm h})}_{\rm BB}(\bm{q}_{\rm 3D},\omega)]$ is zero. Because of the convolution in Eq. (21), this ceases to be true for ${\rm Im}[\chi^{\rm eff}_{\rm BB}({\bm{q}}_{\parallel},\omega)]$.

In the limit $\omega\gg vq_{\parallel}$ and at the leading order in $2E_{\rm F}/(\hbar\omega)\gg 1$ we obtain

[TABLE]

In the limit $k_{\rm F}\to 0$ we find ${\rm Im}[\chi^{\rm eff}_{\rm BB}({\bm{q}}_{\parallel},\omega)]=-q_{\|}/(24\pi\hbar v)$. This implies that, in the undoped limit, the decay rate $\Gamma$ tends to a finite value in the long-wavelength $q_{\|}\to 0$ limit. The FA plasmon is therefore not a well defined excitation in the case of an undoped WSM, since it easily decays by emitting inter-band electron-hole pairs in the gapless bulk.

Finally, we need to calculate ${\rm Im}[\chi^{\rm eff}_{\rm AB}(\bm{q}_{\|},\omega)]$. Following the steps described in Appendix E, we find

[TABLE]

We emphasize that the contribution to the decay rate coming from Eq. (23) is suppressed with respect to that coming from Eq. (22) by a factor $[2E_{\rm F}/(\hbar\omega)]^{-5/2}\ll 1$. We will therefore neglect the contribution (23) for the calculation of $\Gamma$. Physically, Eq. (23) represents processes in which a FA plasmon decays by emitting a composite electron-hole pair, with one partner of the pair belonging to the FA manifold of states and the other partner to the bulk manifold. In the limit $k_{\rm F}\to 0$, we get ${\rm Im}[\chi^{\rm eff}_{\rm AB}(\bm{q}_{\|},\omega)]=-2q^{2}_{\|}[1+\sin^{2}(\theta)]/(3\pi^{2}\hbar\omega)$.

In summary, for a doped WSM we get the following compact expression for the decay rate of a FA plasmon in the long-wavelength limit and deep in the single-particle optical gap $\hbar\omega\ll 2E_{\rm F}$:

[TABLE]

The quantity $\Gamma({\bm{q}}_{\parallel})/\Omega({\bm{q}}_{\parallel})$ is plotted in Figs. 2(b) and 3(b). In writing Eq. (24) we have omitted the contribution $\propto\delta(\omega-E_{\rm A}({\bm{q}}_{\parallel})/\hbar)$ stemming from ${\rm Im}[\chi^{\rm eff}_{\rm AA}]$. This $\delta$-function contribution is not present in the range of values of ${\bm{q}}_{\|}$ we have used to make the plot in Fig. 2(b). From Fig. 2(b) we immediately see that the WSM surface plasmons are highly directional and weakly damped only in the direction of propagation of the single-particle FA states. From Figs. 2(b)-3(b) we clearly see that damping increases (quality factor decreases) very rapidly as a function of the angle $\theta$.

Before concluding, we would like to comment on one evident limitation of the WSM model in Eq. (1). As discussed above, in this model the FA states disperse only in the $\hat{\bm{k}}_{y}$ direction and their dispersion is strictly linear, $E_{\rm A}({\bm{k}}_{\|})=\hbar vk_{y}$. More complicated models—see, for example, Ref. tchoumakov_prb_2017, —allow of course for more general FA dispersion relations. For example, taking $E_{\rm A}({\bm{k}}_{\|})=\hbar vk_{y}+\hbar^{2}(k_{x}^{2}-b^{2})/(2m_{x})$, one can demonstrate that Eq. (16) is not modified by the parabolic correction. Also, one can show that the first non-local correction in Eq. (17) is modified only for $\theta\sim\pm\pi/2$, i.e. along the directions parallel or anti-parallel to the Weyl-node separation vector $2{\bm{b}}$. Along these directions, however, the FA plasmon is anyway strongly damped—see Fig. 3(b)—and therefore quantitative changes to its dispersion relation for $\theta\sim\pm\pi/2$ are uninteresting.

VII Response of FA plasmons to an illuminated metallized tip

In this Section we present an elementary theory that describes the potential that one would observe by carrying out a scattering-type near-field optical experiment (s-SNOM) basov_science_2016 on the surface of a WSM with BTRS. In this experiment basov_science_2016 , light is focused on a metallic tip with the aim of launching propagating surface plasmons. In Fig. 4(a) we sketch the setup of such an experiment.

We model the tip as an external potential oscillating in time and with an exponential profile in the $\hat{\bm{z}}$ direction:

[TABLE]

where $\omega_{0}$ is the frequency of the laser that illuminates the tip. The exponential dependence on $z$ dramatically simplifies the algebra. In the Fourier transform, we find $v_{\rm ext}(z,{\bm{q}}_{\parallel},\omega)\sim\delta(\omega-\omega_{0})$. Since we are interested in what happens at the surface of a WSM, we set $z=0$. In this case, the external potential reduces to a delta function of ${\bm{r_{\parallel}}}$: $v_{\rm ext}(0,{\bm{r}}_{\parallel},t)=V_{0}\delta({\bm{r}}_{\parallel})e^{-i\omega_{0}t}$.

Working in the RPA, we can easily write the relationship between $V_{\rm sc}(0,{\bm{q}}_{\parallel},\omega)$ and the external potential $v_{\rm ext}(0,{\bm{q}}_{\parallel},\omega)$, which involves the effective 2D response function $\chi_{\rm eff}({\bm{q}}_{\parallel},\omega_{0})$:

[TABLE]

Near resonance, i.e. when the laser frequency $\omega_{0}$ is close to the WSM FA plasmon frequency $\Omega({\bm{q}}_{\parallel})$, we can write

[TABLE]

Fourier-transforming the previous result to real space we find

[TABLE]

This integral can be easily evaluated in polar coordinates:

[TABLE]

where we have introduced the angle $\theta^{\prime}$ defined by ${{\bm{r}}_{\parallel}\cdot{\bm{\hat{y}}}}=r_{\parallel}\cos(\theta^{\prime})$.

Using the result $\partial_{\omega}{\rm Re}[\chi_{\rm eff}({\bm{q}}_{\parallel},\omega)|_{\omega=\Omega({\bm{q}}_{\parallel})}]\to-q_{\parallel}[\cos(\theta)b/2+vk^{2}_{\rm F}/(3\Omega_{\theta})]/(\hbar\Omega^{2}_{\theta}\pi^{2})$ reported above, and noting that the integrand is peaked where the denominator vanishes, we can extend the integration range over $q_{\parallel}$ from $[0,\infty]$ to $[-\infty,\infty]$. This allows us to use the power of the residue theorem.

The integrand has a simple pole with respect to $q_{\parallel}$. We define it as $\tilde{p}(\omega_{0},\theta)=p_{1}+ip_{2}$, which satisfies the equation $\omega_{0}=\Omega(\tilde{p},\theta)-i\Gamma(\tilde{p},\theta)$. Using the definitions of $\mathcal{I}(\theta)$ and $\Gamma_{\theta}$ given above, we have

[TABLE]

and

[TABLE]

After straightforward mathematical manipulations we find

[TABLE]

which needs to be evaluated numerically.

Fig. 4(b) displays the calculated screened potential for the following parameters: $b/k_{\rm F}=5$, $\alpha_{\rm ee}=0.8$, and $\hbar\omega/E_{\rm F}=1.4$. Here, it is possible to recognize features due to the peculiar propagation dynamics of FA plasmons. The illustrated SNOM pattern is highly anisotropic because of the highly unidirectional character of both FA dispersion and damping. The anisotropy of the dispersion relation—shown in Fig. 2(a)—manifests through the FA plasmon group velocity, which is maximal around the angular directions $\theta\approx\pm\pi/4$. The anisotropy of the plasmon dissipation—see Fig. 2(b)—shows up in the propagation direction of FA plasmons, which predominantly occurs only along the $y>0$ direction.

VIII Summary and conclusions

In summary, we have presented a quantum-mechanical non-local theory of surface plasmons in semi-infinite Weyl semimetals with broken time-reversal symmetry. We have been able to derive a simple analytical formula—see Eqs. (16)-(18)—for the surface plasmon dispersion relation in the electrostatic limit $|{\bm{q}}_{\|}|\gg\omega/c$, which takes into account exquisite quantum effects associated with the penetration of the Fermi arc surface states into the gapless bulk. We have also included non-local corrections, which were crucial to investigate in a quantitative manner the surface plasmon damping rate (24), as determined by decay processes involving the excitation of electron-hole pairs.

Our calculations show that the intrinsic damping of topological Fermi arc plasmons is small at small values of the in-plane wave vector ${\bm{q}}_{\|}$, mainly in a specific directions, i.e. the direction along which the single-particle Fermi arc states disperse. Scattering-type near-field optical spectroscopy can therefore be used as an alternative to ARPES to carry out spatially-resolved investigations of these intriguing chiral modes occurring in Weyl semimetals with broken time-reversal symmetry, in the absence of an external magnetic field.

Acknowledgements.

It is a great pleasure to thank Iacopo Torre and Fabio Taddei for very useful discussions. M.P. wishes to thank Fondazione Istituto Italiano di Tecnologia for financial support. F.H.L.K. acknowledges financial support from the Government of Catalonia through the SGR grant (2014-SGR-1535), and from the Spanish Ministry of Economy and Competitiveness, through the ”Severo Ochoa” Programme for Centres of Excellence in R&D (SEV-2015-0522), support by Fundacio Cellex Barcelona, CERCA Programme / Generalitat de Catalunya and the Mineco grants Ramón y Cajal (RYC-2012-12281) and Plan Nacional (FIS2013-47161-P and FIS2014-59639-JIN). Furthermore, the research leading to these results has received funding from the European Union’s Horizon 2020 research and innovation programme under grant agreement No. 696656 “GrapheneCore1” and the ERC starting grant (307806, CarbonLight). M.P. is extremely grateful for the financial support granted by ICFO during a visit in August 2016.

Appendix A Boundary conditions

In order to determine the correct boundary conditions, we solve the problem for $z>0$, evaluating the evanescent states at $z=0^{+}$, and then we take the $\lim_{\tilde{m}\to+\infty}{\cal H}_{m(z)}({\bm{k}}_{\|},k_{z}\to-i\partial_{z})$ in the semi-infinite WSM model Hamiltonian introduced in the main text. We find

[TABLE]

where $\Psi_{E,{\bm{k}}_{\parallel}}({\bm{r}}_{\parallel},z)$ is the evanescent solution of the problem for $z>0$. Imposing the continuity of the wave function at $z=0$ between the reflected state in Eq. (4) of the main text and the evanescent state (33) we find $r_{{\bm{k}}}=-1$. Note that the square modulus of the reflected state in Eq. (4) is correctly normalized only in the limit $V\rightarrow\infty$.

Appendix B Dielectric background screening

We here describe how to include the effect of a dielectric background in Eq. (7) of the main text.

We assume that dielectric screening stems from the response of the empty (occupied) states in the bands above (below) the conduction (valence) band, far away in energy from the Weyl crossing.

The total response $\chi_{\rm tot}$ of the system is the sum of two terms: i) the response of the states described by the Weyl Hamiltonian in Eq. (1) of the main text, denoted by the symbol $\chi$, and ii) the contribution due to the high-energy bands, denoted by the symbol $\chi_{\rm b}$, i.e.

[TABLE]

For simplicity, we start from the homogeneous (h) case and we focus on a single pair of bands, with one energetically higher than the conduction band and the other one energetically lower than the valence band. By using the Lehmann representation we can write

[TABLE]

where the collective labels “L” and “H” refer to states of a given electronic band with lower (higher) energy than the valence (conduction) band, respectively, and ${\cal M}_{{\rm L},{\rm H}}({\bm{q}}_{\rm 3D})$ represents a suitable matrix element. Because the band L (H) is deep below (well above) the Fermi level, we replace $f_{\rm L}\to 1$ ($f_{\rm H}\to 0$). Furthermore, we define $\Delta\equiv\min(\epsilon_{\rm H}-\epsilon_{\rm L})>0$. Since we are interested in the frequency range $\omega\ll\Delta/\hbar$ we can approximate Eq. (35) as

[TABLE]

In the long-wavelength limit, because of the orthogonality of the states belonging to different bands, we have $|\mathcal{M}_{{\rm L},{\rm H}}({\bm{q}}_{\rm 3D})|^{2}\propto{\bm{q}}_{\rm 3D}^{2}$. We can therefore write

[TABLE]

where $\chi_{0}$ is a positive constant, independent of ${\bm{q}}$ and $\omega$. By exploiting the generality of the above argument, we can express the global contribution to the response function from all electronic bands far away in energy from the Weyl crossing as $\chi^{(\rm h)}_{\rm b}({\bm{q}}_{\rm 3D},\omega)\simeq-\chi_{\rm b}{\bm{q}}_{\rm 3D}^{2}$, where $\chi_{\rm b}>0$.

In order to recover the plasmon dispersion relation in the homogeneous case, we have to solve the well-known RPA equation Giuliani_and_Vignale

[TABLE]

or, equivalently,

[TABLE]

where $\epsilon_{\rm b}=1+4\pi e^{2}\chi_{\rm b}$.

We now follow a similar path for the case of surface plasmon modes. We first replace the expression $\chi_{\rm b}(z,z^{\prime},{\bm{q}}_{\parallel},\omega)\simeq\Theta(-z)\Theta(-z^{\prime})[\chi^{(\rm h)}_{\rm b}(z+z^{\prime},{\bm{q}}_{\parallel})+\chi^{(\rm h)}_{\rm b}(z-z^{\prime},{\bm{q}}_{\parallel},\omega)]$ in Eq. (34). We then use the result derived earlier, i.e. $\chi^{(\rm h)}_{\rm b}({\bm{q}}_{\rm 3D},\omega)\simeq-\chi_{\rm b}{\bm{q}}_{\rm 3D}^{2}$. Carrying out straightforward algebraic manipulations we find Eq. (7) in the main text, i.e.

[TABLE]

where $\bar{\epsilon}=(1+\epsilon_{\rm b})/2$.

Appendix C The AA response function

We here provide more details on the calculation of the response function ${\chi}^{(0)}_{\rm AA}(q_{z},q_{z}^{\prime},{\bm{q}}_{\parallel},\omega)$. Using the FA eigenstates in Eq. (3) of the main text and the Lehmann representation Giuliani_and_Vignale —Eq. (10) of the main text—we find the following exact expression:

[TABLE]

where

[TABLE]

Because in our model the FA states do not disperse along the $\hat{\bm{x}}$ direction, we have ${\chi}^{(0)}_{\rm AA}(q_{z},q_{z}^{\prime},q_{x},q_{y}=0,\omega)=0$. This means that, in our model, an external field can perturb the FA states only if it carries a finite momentum along the $\hat{\bm{y}}$ direction. Expanding $\mathcal{L}(q_{z},q_{z}^{\prime},q_{x})$ in a power series of $q$ and carrying out the integral in Eq. (C), we obtain

[TABLE]

From Eq. (43) one easily obtains ${\rm Re}[{\chi}^{\rm eff}_{\rm AA}({\bm{q}}_{\|},\omega)]$, as in Eq. (12) of the main text. Similarly, we find

[TABLE]

from which ${\rm Im}[{\chi}^{\rm eff}_{\rm AA}({\bm{q}}_{\|},\omega)]$ follows.

Appendix D The BB response function

We here report some useful technical details on the calculation of the response function ${\chi}^{(0)}_{\rm BB}(q_{z},q_{z}^{\prime},{\bm{q}}_{\parallel},\omega)$ for the doped case.

Since we are interested in studying the long-wavelength limit, we can use the (high-frequency) moment expansion Giuliani_and_Vignale ; torre_prb_2017 , i.e.

[TABLE]

where $m_{\rm eff}=E_{\rm F}/v^{2}$ is the effective mass, $\hat{n}_{q_{z}-q^{\prime}_{z}}=\sum_{j}e^{i(q_{z}-q_{z}^{\prime})\hat{z}_{j}}$ is the Fourier Transform of the density operator Giuliani_and_Vignale , and $\braket{\ldots}_{\rm GS}$ is a shorthand for the following ground-state expectation value $\bra{\rm GS}_{\rm B}\ldots\ket{\rm GS}_{\rm B}$. Here, $\ket{{\rm GS}}_{\rm B}=\prod_{|{\bm{k}}|<k_{\rm F}}b_{+,{\bm{k}}}^{\dagger}\ket{0}$, where $b_{+1,{\bm{k}}}^{\dagger}$ creates an electron in a state labeled by the band index $\lambda=+1$ and wavevector ${\bm{k}}$. The quantity $\ket{0}$ represents the vacuum state with no electrons in the conduction band.

By using the eigenstates in Eq. (4) of the main text, we find that

[TABLE]

where

[TABLE]

and

[TABLE]

where $u_{{\bm{k}},\lambda=\pm 1}$ has been introduced in Eq. (2) of the main text. As in the main text, $\bar{\bm{k}}=({\bm{k}}_{\|},-k_{z})$. The index “scl” (“int”) stands for “semiclassical” (“interference”). Indeed, $\braket{\hat{n}_{q_{z}}}_{\rm GS}^{({\rm scl})}$ ($\braket{\hat{n}_{q_{z}}}_{\rm GS}^{({\rm int})}$) is independent of (dependent of) the relative phases, see Eq. (4) of the main text.

We start from the case $q_{z}=0$. We find

[TABLE]

where $N$ is the total number of electrons in the semi-infinite WSM, while the interference term $\braket{\hat{n}_{0}}_{\rm GS}^{\rm int}\propto 1/L_{z}$ is negligible in the thermodynamic limit.

For a generic value of $q_{z}$, we obtain

[TABLE]

where $n\equiv 2N/(L_{z}S)$. We then calculate the interference term $\braket{\hat{n}_{q_{z}}}^{({\rm int})}_{\rm GS}$. We are interested in the small $q_{z}$ limit. We find $\braket{\hat{n}_{q_{z}}}^{({\rm int})}_{\rm GS}-\braket{\hat{n}_{0}}^{({\rm int})}_{\rm GS}\propto q_{z}$. In the thermodynamic limit, this can induce in ${\rm Re}[{\chi}^{(0)}_{\rm BB}(q_{z},q_{z}^{\prime},{\bm{q}}_{\parallel},\omega)]$ terms of ${\cal O}(q^{3}_{z},q^{\prime 3}_{z})$, which are beyond the interest of our leading-order long-wavelength theory.

Finally, taking into account only the semiclassical term in Eq. (45), we find

[TABLE]

which implies Eq. (14) of the main text.

Appendix E The AB response function

Here, we detail the calculation of the response function $\chi^{(0)}_{\rm AB}$, always for a generic doping. We neglect inter-Weyl-node electron-electron scattering processes and we linearize the spectrum $E_{\lambda}({\bm{k}})$ around each Weyl node. In order to simplify the notation we set $\hbar=1$. Under these simplifying assumptions, the quantity we need to calculate is

[TABLE]

where $\eta=0^{+}$, $j=\pm 1$ is a Weyl node index, $\lambda=\pm$ is a conduction/valence band index, $k_{+1}=k_{\rm F}$, and $k_{-1}=\Lambda$. Here, $\Lambda$ is an ultraviolet cutoff. We will take the limit $\Lambda\to\infty$ momentarily. The matrix element $\mathcal{M}^{j}_{\lambda}(q_{x},q_{z},q_{z}^{\prime})$ reads as follows

[TABLE]

We denote by the shorthand $q$ the small quantity $q_{z}\sim q_{z}^{\prime}\sim q_{x}$, and we expand the response function up to the quadratic order in $q$. Carrying out the sum over the index $j$ we find

[TABLE]

It is now useful to define the following dimensionless quantities, ${\bm{p}}\equiv{{\bm{k}}}/{k_{\rm F}}$ and $\tilde{\omega}\equiv{\omega}/{E_{\rm F}}$, and use cylindrical coordinates, $p_{z}=p_{\perp}\cos(\phi)$ and $p_{y}=p_{\perp}\sin(\phi)$. We perform the integral over $\phi$ and, exploiting the identity \lim_{\Lambda\rightarrow\infty}\big{[}\Theta(1-p_{y})-\Theta(\Lambda/k_{\rm F}-p)\big{]}=-\Theta(p_{y}-1), we obtain

[TABLE]

where we have defined

[TABLE]

We first analyze the real part of the response function. We start by studying the auxiliary function $\mathcal{J}_{1}(\omega)$. Introducing a new set of polar coordinates, $p_{y}=p\sin(\theta)$ and $p_{\perp}=p\cos(\theta)$, we find

[TABLE]

where “${\rm PV}$” denotes the Cauchy principal value. The following mathematical identity will be used below:

[TABLE]

Introducing the auxiliary variable $t=1-\sin(\theta)$ we can write

[TABLE]

We now note that the integrand in the previous equation is peaked at $t\sim 0$ and $t\sim\tilde{\omega}$. Furthermore, for $\omega\ll 2E_{\rm F}$, both peaks collapse at $t\sim 0$. We can therefore approximate (61) as

[TABLE]

which can be evaluated analytically, yielding

[TABLE]

A further asymptotic expansion in the limit $\omega\ll 2E_{\rm F}$ yields

[TABLE]

Similarly, we find

[TABLE]

while, in the same limits, ${\rm Re}[\mathcal{J}_{3}(\omega\ll 2E_{\rm F})]$ is a regular function having a negligible contribution with respect to ${\rm Re}(\mathcal{J}_{1})$ and ${\rm Re}(\mathcal{J}_{2})$. Replacing Eqs. (64)-(65) in Eq. (E) we finally obtain

[TABLE]

which gives Eq. (14) of the main text.

In order to evaluate the imaginary part of the response function we exploit the Dirac identity {1}/\big{[}{x^{2}-(\omega+i0^{+})^{2}}\big{]}={1}/\big{(}{x^{2}-\omega^{2}}\big{)}+{i\pi}\delta({x^{2}-\omega^{2}}). For $\omega<2E_{\rm F}$ we obtain

[TABLE]

which, replaced in Eq. (E), result into

[TABLE]

The leading term for $\omega\ll E_{\rm F}$ gives Eq. (23) of the main text.

Following similar steps in the undoped case ($E_{\rm F}=0$), we find

[TABLE]