Designer curved-space geometry for relativistic fermions in Weyl metamaterials

TL;DR

This paper introduces Weyl metamaterials where relativistic fermions experience a tunable 3D curved geometry and gauge fields, enabling new electronic device functionalities such as a 3D electron lens effect.

Contribution

It establishes a systematic way to design 3D geometries and gauge fields for relativistic carriers in Weyl semimetals based on symmetry-breaking configurations.

Findings

Derived explicit connection between geometry and symmetry-breaking fields.

Proposed fabrication methods using inhomogeneous magnetization.

Demonstrated potential for novel electronic devices like electron lenses.

Abstract

Weyl semimetals are recently discovered materials supporting emergent relativistic fermions in the vicinity of band-crossing points known as Weyl nodes. The positions of the nodes and the low-energy spectrum depend sensitively on the time-reversal (TR) and inversion (I) symmetry breaking in the system. We introduce the concept of Weyl metamaterials where the particles experience a 3d curved geometry and gauge fields emerging from smooth spatially varying TR and I breaking fields. The Weyl metamaterials can be fabricated from semimetal or insulator parent states where the geometry can be tuned, for example, through inhomogeneous magnetization. We derive an explicit connection between the effective geometry and the local symmetry-breaking configuration. This result opens the door for a systematic study of 3d designer geometries and gauge fields for relativistic carriers. The Weyl metamaterials provide a route to novel electronic devices as highlighted by a remarkable 3d electron lens effect.

Tables

Table 1. Table 1: The four groups of γ μ ​ ν subscript 𝛾 𝜇 𝜈 \gamma_{\mu\nu} and their transformation properties.

	TR	I
$𝐛=\left({\gamma }_{23}{\gamma }_{31}{\gamma }_{12}\right)$	$-1$	$+1$
$𝐩=\left({\gamma }_{14}{\gamma }_{24}{\gamma }_{34}\right)$	$+1$	$-1$
${𝐛}^{\prime }=\left({\gamma }_{15}{\gamma }_{25}{\gamma }_{35}\right)$	$-1$	$+1$
$\epsilon ={\gamma }_{45}$	$+1$	$-1$

Equations

$$H_0=n(\mathbf{k})\mathbb{I}+\mathbf{\kappa }(\mathbf{k})\cdot \mathbf{\gamma }+m(\mathbf{k}){\gamma }_4,$$

$$H=H_0+\mathbf{u}\cdot \mathbf{b}+\mathbf{w}\cdot \mathbf{p}+{\mathbf{u}}^{\prime }\cdot {\mathbf{b}}^{\prime }+f\epsilon ,$$

$$H^{\prime }=\left(\begin{array}{cc}{D}_0(\mathbf{k}){\sigma }^0+\mathbf{D}(\mathbf{k})\cdot \mathbf{\sigma }& 0\\ 0& d_0(\mathbf{k}){\sigma }^0+\mathbf{d}(\mathbf{k})\cdot \mathbf{\sigma },\end{array}\right),$$

$$\begin{array}{rl}{d}_1& =-\frac{(\hat{\mathbf{u}}\times \mathbf{\kappa }(\mathbf{k}))\cdot \hat{\mathbf{z}}}{\sqrt{1-{\hat{u}}_3^2}}\\ d_2& =-\frac{[\hat{\mathbf{u}}\times (\mathbf{\kappa }(\mathbf{k})\times \hat{\mathbf{u}})]\cdot \hat{\mathbf{z}}}{\sqrt{1-{\hat{u}}_3^2}}\\ d_3& =u-\sqrt{(\hat{\mathbf{u}}\cdot \mathbf{\kappa }(\mathbf{k}){)}^2+m^2},\end{array}$$

$$H_W^{\prime }=d_0({\mathbf{k}}_W){\sigma }^0+\text{tensor}e{}_{\nu }^j(\mathbf{r}){\sigma }^{\nu }(k_j-k_j^W).$$

$$\text{tensor}e{}_{\nu }^i(\mathbf{r})={\frac{\partial d_{\nu }}{\partial k_i}|}_{\mathbf{k}=\mathbf{k}W},$$

$$\mathfrak{g}=\left(\begin{array}{cc}-1& 0\\ 0& \mathbf{g}\end{array}\right),$$

$$g^{ij}=({\partial }_{k_i}\mathbf{\kappa })\cdot ({\partial }_{k_j}\mathbf{\kappa })-\frac{m^2}{u^2}({\partial }_{k_i}\mathbf{\kappa }\cdot \hat{\mathbf{u}})({\partial }_{k_j}\mathbf{\kappa }\cdot \hat{\mathbf{u}})$$

$$H_W$$

$$=\text{tensor}e{}_{\nu }^j{\sigma }^{\nu }\left(-i{\partial }_j-q{a}_j-\frac{i}{2}\text{tensor}\tilde{e}{}_j^{\mu }{\partial }_n\text{tensor}e{}_{\mu }^n\right).$$

$${\dot{r}}^i={\partial }_{k_i}\mathcal{E}-{\Omega }_{k_i{r}^l}{\dot{r}}^l-{\Omega }_{k_i{k}_l}{\dot{k}}_l$$

$$\dot{k_i}=-{\partial }_{r^i}\mathcal{E}+{\Omega }_{r^i{k}_l}{\dot{k}}_l+{\Omega }_{r^i{r}^l}{\dot{r}}^l-q(E_i+{\epsilon }_{ijk}{\dot{r}}^j{B}^k),$$

$$\Delta z=2\pi \frac{v_z}{v\omega }\sqrt{\frac{1-\frac{m^2}{u^2}}{\frac{m^2}{u^2}+\frac{u}{v_z{k}_z}\sqrt{1-\frac{m^2}{u^2}}}},$$

$$H=\mathbf{\kappa }(\mathbf{k})\cdot \mathbf{\gamma }+m(\mathbf{k}){\gamma }_4+\mathbf{u}\cdot \mathbf{b}.$$

$$\{b_i,b_j\}=2{\delta }_{ij},[b_i,b_j]=2i{\epsilon }_{ijk}{b}_k.$$

$$\mathbf{u}\cdot \mathbf{b}\to e^{i\frac{\phi }{2}\hat{\mathbf{n}}\cdot \mathbf{b}}(\mathbf{u}\cdot \mathbf{b})e^{-i\frac{\phi }{2}\hat{\mathbf{n}}\cdot \mathbf{b}}$$

$$\mathcal{U}(\mathbf{u})=e^{i{b}_2\frac{(\theta -\pi /2)}{2}}{e}^{i{b}_3\frac{\phi }{2}},$$

$$\begin{array}{r}\mathcal{U}(\mathbf{u})H{\mathcal{U}}^{†}(\mathbf{u})=(\mathbf{\kappa }(\mathbf{k})\cdot \hat{\mathbf{u}}){\gamma }_1+\frac{(\hat{\mathbf{u}}\times \mathbf{\kappa }(\mathbf{k}))\cdot \hat{\mathbf{z}}}{\sqrt{1-{\hat{u}}_3^2}}{\gamma }_2\\ +\frac{[\hat{\mathbf{u}}\times (\mathbf{\kappa }(\mathbf{k})\times \hat{\mathbf{u}})]\cdot \hat{\mathbf{z}}}{\sqrt{1-{\hat{u}}_3^2}}{\gamma }_3+m(\mathbf{k}){\gamma }_4+u{\gamma }_{23}.\end{array}$$

$$\begin{array}{rl}& \mathcal{RU}(\mathbf{u})H{\mathcal{U}}^{†}(\mathbf{u}){\mathcal{R}}^{†}=u{\gamma }_{23}+\frac{(\hat{\mathbf{u}}\times \mathbf{\kappa }(\mathbf{k}))\cdot \hat{\mathbf{z}}}{\sqrt{1-{\hat{u}}_3^2}}{\gamma }_2\\ & +\frac{[\hat{\mathbf{u}}\times (\mathbf{\kappa }(\mathbf{k})\times \hat{\mathbf{u}})]\cdot \hat{\mathbf{z}}}{\sqrt{1-{\hat{u}}_3^2}}{\gamma }_3+m(\mathbf{k})\sqrt{1+\frac{(\mathbf{\kappa }(\mathbf{k})\cdot \hat{\mathbf{u}}{)}^2}{m^2(\mathbf{k})}}{\gamma }_4,\end{array}$$

$$\begin{array}{rl}{\gamma }_1& ={\tau }_x,{\gamma }_2={\tau }_z{\sigma }_x\\ {\gamma }_3& ={\tau }_z{\sigma }_y,{\gamma }_4={\tau }_z{\sigma }_z.\end{array}$$

$$\begin{array}{rl}{d}_1^{\pm }& =\pm \frac{(\hat{\mathbf{u}}\times \mathbf{\kappa }(\mathbf{k}))\cdot \hat{\mathbf{z}}}{\sqrt{1-{\hat{u}}_3^2}}\\ d_2^{\pm }& =\pm \frac{[\hat{\mathbf{u}}\times (\mathbf{\kappa }(\mathbf{k})\times \hat{\mathbf{u}})]\cdot \hat{\mathbf{z}}}{\sqrt{1-{\hat{u}}_3^2}}\\ d_3^{\pm }& =u\pm m(\mathbf{k})\sqrt{1+\frac{(\mathbf{\kappa }(\mathbf{k})\cdot \hat{\mathbf{u}}{)}^2}{m^2(\mathbf{k})}}\end{array}$$

$$\begin{array}{rl}{d}_1^W& =\frac{(\hat{\mathbf{u}}\times \mathbf{\kappa }(\mathbf{k}))\cdot \hat{\mathbf{z}}}{\sqrt{1-{\hat{u}}_3^2}}\\ d_2^W& =\frac{[\hat{\mathbf{u}}\times (\mathbf{\kappa }(\mathbf{k})\times \hat{\mathbf{u}})]\cdot \hat{\mathbf{z}}}{\sqrt{1-{\hat{u}}_3^2}}\\ d_3^W& =u-\sqrt{m(\mathbf{k}{)}^2+(\mathbf{\kappa }(\mathbf{k})\cdot \hat{\mathbf{u}}{)}^2}\end{array},$$

$$\mathbf{\kappa }({\mathbf{k}}^W)=\sqrt{u^2-m^2({\mathbf{k}}^W{)}^2}\hat{\mathbf{u}}.$$

$$H^W=\text{tensor}e{}_i^j(k_j-k_j^W){\sigma }^i,$$

$$\text{tensor}e{}_1^j$$

$$\text{tensor}e{}_2^j$$

$$\text{tensor}e{}_3^j$$

$$\mp {\sqrt{1-\frac{m^2({\mathbf{k}}^W)}{u^2}}{\partial }_{q_j}(\mathbf{\kappa }(\mathbf{q})\cdot \hat{\mathbf{u}})|}_{\mathbf{q}={\mathbf{k}}^W},$$

$$\text{tensor}\tilde{e}{}_l^i=\frac{1}{2}\frac{{\epsilon }^{ijk}{\epsilon }_{lmn}\text{tensor}e{}_j^m\text{tensor}e{}_k^n}{{\epsilon }_{stu}\text{tensor}e{}_1^s\text{tensor}e{}_2^t\text{tensor}e{}_3^u},\text{tensor}\tilde{e}{}_0^0=1.$$

$$g^{ij}=({\partial }_{k_i}\mathbf{\kappa })\cdot ({\partial }_{k_j}\mathbf{\kappa })-\frac{m^2}{u^2}({\partial }_{k_i}\mathbf{\kappa }\cdot \hat{\mathbf{u}})({\partial }_{k_j}\mathbf{\kappa }\cdot \hat{\mathbf{u}})\pm \frac{m}{u}\sqrt{1-\frac{m^2}{u^2}}\left(({\partial }_{k_i}m)({\partial }_{k_j}\mathbf{\kappa }\cdot \hat{\mathbf{u}})+({\partial }_{k_j}m)({\partial }_{k_i}\mathbf{\kappa }\cdot \hat{\mathbf{u}})\right)+\frac{m^2}{u^2}({\partial }_{k_i}m)({\partial }_{k_j}m),$$

$$H=\mathbf{\kappa }(\mathbf{k})\cdot \mathbf{\gamma }+m(\mathbf{k}){\gamma }_4+\mathbf{w}\cdot \mathbf{p}.$$

$$\mathcal{U}(\mathbf{w})H{\mathcal{U}}^{†}(\mathbf{w})=\mathbf{K}\cdot \mathbf{\gamma }+m(\mathbf{k}){\gamma }_4+w{\gamma }_{14},$$

$$\begin{array}{rl}{e}^{i\frac{\nu }{2}{\gamma }_{23}}\mathcal{U}(\mathbf{w})H{\mathcal{U}}^{†}(\mathbf{w})e^{-i\frac{\nu }{2}{\gamma }_{23}}& =K_1{\gamma }_1+K_2\sqrt{1+\frac{K_3^2}{K_2^2}}{\gamma }_2\\ & +m(\mathbf{k}){\gamma }_4+w{\gamma }_{14}.\end{array}$$

Full text

I Introduction

Simulating relativistic phenomena in table-top systems has become a major theme in condensed matter physics volovik . Topological materials bernevig , exhibiting a myriad of connections to high-energy physics volovik , have been a significant inspiration for these developments. An important source of fascination has also been provided by semimetals which display emergent relativistic dynamics at low energies katsnelson ; hosur ; burkov1 ; wan ; lu ; xu ; liu ; liu1 ; neupane ; wang ; borisenko ; shekhar ; xu2 ; yang ; lv ; cano . This has given rise to a wide-spread interest in engineering artificial gauge fields in graphene katsnelson ; vozmediano ; amorim and 3d Dirac and Weyl semimetals cortijo1 ; cortijo2 ; grushin ; pikulin ; shapourian ; zabolotskiy . The phenomenology of general relativity and curved-space dynamics has also penetrated into condensed-matter research volovik . In addition to fundamental interest, curved-space physics may also have striking practical applications as electromagnetic metamaterials and transformation optics demonstrate pendry ; leonhardt1 ; leonhardt2 .

In this work, we study a class of structures which we call Weyl metamaterials. In these systems, the chiral Weyl fermions familiar from Weyl semimetals, depicted in Fig. 1a, are moving in an artificial 3d curved-space geometry. The great advantage of these systems is that, first of all, the effective geometry experienced by the particles can be widely tuned by external conditions. Since the band crossings in Weyl semimetals result from accidental degeneracy, the positions of the nodes and the low-energy spectrum can be tuned by a large amount by varying time-reversal (TR) and inversion (I) symmetry breaking perturbations, as depicted in Fig. 1b. In addition, a 3d curved geometry can be realized by atomically smooth TR or I breakings whose variation is significant only on the scale of the sample size, as illustrated in Fig. 1c. These two factors combined with the obvious richness offered by 3d space give Weyl metamaterials flexibility and enormous advantage over 2d materials katsnelson ; vozmediano ; amorim as a platform for relativistic curved-space physics. Weyl metamaterials can be fabricated from 3d semimetal or insulator parent states by generic inhomogeneous TR and I breaking couplings. Hence the concept reaches substantially beyond the strain-engineering through elastic deformations cortijo1 ; cortijo2 ; grushin ; pikulin ; shapourian . More generally, Weyl metamaterials pave the way for novel 3d electronic devices through curvature engineering as sketched in Fig. 1d. In this article, we provide an example of such a device where the Weyl metamaterial acts as a 3d lens for the itinerant charge carriers. A magnetization with a radial gradient in the plane perpendicular to the initial velocities gives rise to an effective curved space and emergent gauge field, which cause the trajectories of incoming charge carriers at different radii to converge at a focal point.

In Sec. II, starting from a generic four-band model, we derive a curved-space Weyl Hamiltonian $H=\tensor{e}{{}^{\mu}_{i}}(k_{\mu}-a_{\mu})\sigma^{i}$ describing the low-energy dynamics of Weyl metamaterials. The emergent gauge field $a_{\mu}$ and the frame fields $\tensor{e}{{}^{\mu}_{i}}$ which encode the metric tensor $g^{\mu\nu}$ are solved as functions of the TR and I breaking fields. These expressions provide the fundamental tool for systematic reverse-engineering of synthetic curved-space geometries and gauge fields in Weyl metamaterials. Our results are applicable to a wide class of different lattice and continuum models as illustrated by examples. The general theory is illustrated in Sec. III by proposing simple magnetic textures that give rise to remarkable electron lensing effects. In Sec. IV, we briefly discuss the experimental feasibility of inhomogeneous symmetry-breaking fields, and in the final section summarize our work and results, as well as propose directions for future research.

II Theory of Weyl metamaterials

II.1 Engineering Weyl metamaterials from general parent states

The central idea of our work is to start from a generic insulating or semimetal phase and introduce smooth spatially varying TR and I breaking (TRIB) fields that push the parent system into an inhomogeneous Weyl semimetal phase. The TR and I breaking terms are not assumed to be small in any sense, however their spatial dependence is assumed to be smooth. Besides an intrinsic Weyl semimetal, the parent state for a Weyl metamaterial could be a Dirac semimetal or TR invariant topological insulator bernevig ; zhang as it is well known that these systems generically undergo a transition to Weyl semimetals when the I murakami ; ojanen or TR symmetry is broken. A possible mechanism to break TR in topological insulators is magnetic doping, which has been demonstrated in thin films chang1 ; chang2 ; bestwick .

As the minimum models of Weyl semimetals have four energy bands, we write down the Clifford representation of the most general TR and I invariant four-band models and identify all TR and I breaking terms as was also done in Ref. burkov3 . Without any TRIB fields, the most general four-band Hamiltonian is

[TABLE]

where the parameters $n$, and $m$ are symmetric in the momentum $m({\bf k})=m(-{\bf k})$, while the kinetic term is antisymmetric $\boldsymbol{\kappa}({\bf k})=-\boldsymbol{\kappa}(-{\bf k})$. We have denoted the unit matrix by $\mathbb{I}$. Furthermore, $\gamma_{\mu}$ denotes the five $4\times 4$ gamma matrices satisfying the anticommutation relations $\{\gamma_{\mu},\gamma_{\nu}\}=2\delta_{\mu\nu}$. Since the Hamiltonian (1) respects TR and I symmetry, it follows that $\gamma_{1,2,3}$ are odd under both I and TR, while $\gamma_{4}$ must be even. From this, it is then evident that $\gamma_{5}=\gamma_{1}\gamma_{2}\gamma_{3}\gamma_{4}$ must also be odd under both operations. To write down the most general Hermitian $4\times 4$ Hamiltonian, we need to include ten additional matrices constructed using commutators of the gamma matrices, i.e. $\gamma_{ij}\equiv-i[\gamma_{i},\gamma_{j}]/2$. The symmetries of $\gamma_{ij}$ can be easily deduced from the symmetries of their constituent gamma matrices. As it turns out, we can group all ten matrices into four distinct sets: three sets consisting of three mutually anticommuting matrices each, and a fourth set consisting only of the matrix $\gamma_{45}$. The sets, as well as the transformation properties of each set under TR or I, are summarized in table 1 – in fact, we can see that each group breaks either TR or I symmetry.

The general $4\times 4$ Hamiltonian is expressed as

[TABLE]

where the functions ${\bf u}$, ${\bf w}$, ${\bf u}^{\prime}$ and $f$ characterize the TRIB fields and fix the position of the Weyl nodes and the low-energy spectrum. To extract the low-energy Weyl Hamiltonian, we must block diagonalize (2). To make analytical progress, we concentrate on the case ${\bf u}^{\prime}=0$, $f=0$. This restriction leaves room for substantial generality, allowing for a full three-component TR breaking field ${\bf u}$ and I breaking field ${\bf w}$. To illustrate the wide applicability of the general formulation, we give in App. B three concrete examples of Weyl metamaterials: one based on a topological insulator-ferromagnet layer structure burkov2 , a popular theoretical model describing magnetically-doped topological insulators vazifeh , and a 3d Dirac semimetal wang ; borisenko with broken inversion symmetry.

We note that the recently studied photonic Weyl systems lu ; xiao2 employ a two-mode approximation and, as such, falls outside the scope of the our formalism developed with electronic band structures in mind where the four-band Clifford representation of the Hamiltonian (2) is the fundamental starting point.

II.2 Weyl block reduction

Here we derive the local Weyl Hamiltonian from the original four-band model (2). The block reduction is achieved by the unitary transformation theory developed in App. A and is given by

[TABLE]

where Weyl nodes are now found within either the upper or lower block, while the other block is gapped. The low-energy physics in the vicinity of the Weyl nodes is contained in the local Weyl approximation $H^{\prime}_{W}=d_{0}({\bf k})\sigma^{0}+{\bf d}({\bf k})\cdot\boldsymbol{\sigma}$. Here, we have collected the Pauli matrices $\sigma^{1,2,3}$ into a vector $\boldsymbol{\sigma}$, and written the $2\times 2$ identity matrix as $\sigma^{0}$. For concreteness, let us relegate the case where ${\bf w}\neq 0$, ${\bf u}=0$, as well as the general case ${\bf u}\neq 0$, ${\bf w}\neq 0$ to App. A and for now consider the specific case $H=\boldsymbol{\kappa}({\bf k})\cdot\boldsymbol{\gamma}+m\gamma_{4}+{\bf u}\cdot{\bf b}$, which could describe a TR and I invariant insulator or semimetal in the presence of spatially varying magnetization ${\bf M}={\bf u}$. The local Weyl approximation yields $d_{0}=0$ and

[TABLE]

where $u=|{\bf u}|$ and $\hat{{\bf u}}={\bf u}/u$. This system has Weyl nodes ${\bf k}^{W}$ satisfying $\boldsymbol{\kappa}({\bf k}^{W})=\pm\sqrt{u^{2}-m^{2}}\hat{{\bf u}}$ whenever $m<|{\bf u}|$. The gapped block and the results for the general case ${\bf u}\neq 0$ and ${\bf w}\neq 0$ are given in the Appendix. In the local Weyl approximation $H_{W}^{\prime}$, the ${\bf d}$-vector depends not only on ${\bf k}$ but has an implicit position dependence through ${\bf u}$ which is slowly varying in space. The Weyl reduction (3), which is exact for position-independent TRIB fields, is an approximation when ${\bf u}$ and ${\bf w}$ depend on position. The unitary transformations that block diagonalizes the four-band model become spatially dependent through ${\bf u}({\bf r}),{\bf w}({\bf r})$. Regarding ${\bf k}$ in $H_{0}$ as a canonical conjugate operator of ${\bf r}$, the transformation also produces terms proportional to gradients (as well as higher-order spatial derivatives) of ${\bf u}({\bf r}),{\bf w}({\bf r})$. The gradient terms could, for example, give rise to off-diagonal elements in Eq. (3). These terms would give corrections to $H_{W}^{\prime}$ of the order of $\mathcal{O}\left[(\partial_{i}{\bf u})^{2},(\partial_{i}{\bf w})^{2}\right]$ which, along with other gradient corrections, are very small for slowly varying fields. This is confirmed by direct numerical comparisons between the local density of states (LDOS) of the exact four-band model and the local Weyl approximation $H_{W}^{\prime}$, and typically suggest an excellent agreement between the two. The LDOS, defined by $\nu(\mathbf{r},E)=\sum_{E_{n}}|\Psi_{n}(\mathbf{r})|^{2}\delta(E-E_{n})$, probes the local properties of the system and therefore is an ideal tool to investigate the local Weyl approximation. A short explanation on how the computation of the LDOS was performed can be found in App. F. We have calculated the LDOS for representative Weyl-metamaterial models with $H=\boldsymbol{\kappa}({\bf k})\cdot\boldsymbol{\gamma}+m\gamma_{4}+{\bf u}\cdot{\bf b}$, where $\kappa_{i}({\bf k})=t\sin{k_{i}}$ and $m$ and $t$ are constants and different TR-breaking textures ${\bf u}$. This is illustrated in Fig. 2 which establishes an excellent agreement between the exact result and the local Weyl approximation in the energy window $E_{g}$ determined by the energy gap of the gapped block. Away from the Weyl node but within $E_{g}$, the smoothed LDOS exhibits a quadratic envelope typical for Weyl semimetals. At the Weyl node, the LDOS typically becomes modified due to the emergent gauge field effect, as clarified in the next section.

II.3 Designer curved-space geometry for Weyl particles

To gain better insight into the low-energy physics, we expand $H_{W}^{\prime}$ around the local Weyl points. For each Weyl point, we get a low-energy Hamiltonian of the form

[TABLE]

Here we have employed the Einstein summation convention – Greek indices represent integers from 0 to 3, while Latin indices take integer values from 1 to 3, unless otherwise stated. We define the frame fields (or tetrads) as

[TABLE]

where ${\bf k}^{W}({\bf u},{\bf w})$ is the momentum corresponding to the local Weyl point around which we have expanded. Introducing a two-component spinor $\Psi$, the low-energy Schrödinger equation then takes the form $H_{W}^{\prime}\Psi=\mathcal{E}\Psi$ which when squared can be formally written in the form $g^{\mu\nu}p_{\mu}p_{\nu}\Psi=0$ with appropriately defined canonical momenta $p_{\mu}$. When the term proportional to the unit matrix vanishes $d_{0}=0$ in $H_{W}^{\prime}$, this equation describes a massless relativistic particle with dispersion $\mathcal{E}^{2}=g^{ij}p_{i}p_{j}$ moving in a curved space. Here $p_{i}=k_{i}-k^{W}_{i}$ and the contravariant metric $g^{\mu\nu}=\eta^{ab}\tensor{e}{{}^{\mu}_{a}}\tensor{e}{{}^{\nu}_{b}}$ is defined by the standard relation zee using tetrad fields $\tensor{e}{{}^{\mu}_{\nu}}$ which coincide with those defined in Eq. (6) with the addition $\tensor{e}{{}^{0}_{0}}=1$. The covariant metric $g_{\mu\nu}$ is defined as the inverse of $g^{\mu\nu}$ and comes in to play in the semiclassical dynamics through Christoffel symbols.

The frame fields $\tensor{e}{{}^{i}_{\nu}}$ for the model (4) can be straightforwardly evaluated and are given in App. A. We can also easily find the inverse (co-frame) elements $\tilde{e}$ that satisfy $\tensor{e}{{}^{\alpha}_{\nu}}\tensor{\tilde{e}}{{}_{\alpha}^{\mu}}=\delta_{\nu}^{\mu}$ and $\tensor{e}{{}^{\mu}_{\alpha}}\tensor{\tilde{e}}{{}_{\nu}^{\alpha}}=\delta^{\mu}_{\nu}$. The contravariant metric tensor is given by

[TABLE]

where the elements of the $3\times 3$ matrix $[\boldsymbol{g}]^{ij}=g^{ij}$ are

[TABLE]

evaluated at the local Weyl node ${\bf k}={\bf k}^{W}$. In the form (8) we have assumed that $m$ is constant in the vicinity of ${\bf k}^{W}$, the more general case where $m$ depends on ${\bf k}$ is also given in App. A. The relations (4), (6), and (8), along with their generalizations given in App. A, provide the fundamental relation between the effective geometry and the physical symmetry-breaking configuration and serve as a starting point in designing synthetic curved-space geometries for the Weyl quasiparticles. The formula Eq. (8) has a similar role in Weyl metamaterials as the relation between the metric and electromagnetic permeability and permitivity tensors in transformation optics pendry ; leonhardt1 ; leonhardt2 . It shows concretely how a Weyl metamaterial acts as a transformation media where the geometry is tuned by ${\bf u}$.

We postulate an effective long-wavelength quantum Hamiltonian by making the replacement ${\bf k}\to-i\nabla$ in Eq. (5). A direct substitution would leave $H_{W}^{\prime}$ non-Hermitian. This can be remedied by the standard prescription of symmetrizing the non-commuting product of gradient and the frame fields. Assuming for simplicity that $d_{0}=0$, we obtain an effective curved-space Weyl Hamiltonian

[TABLE]

Here we have defined an effective gauge field ${\bf a}={\bf k}^{W}({\bf r})/q$ – $q$ being the charge of the quasiparticle – which gets a small correction due to the symmetrization. Though the correction is parametrically small since it depends on the derivatives of the TRIB fields, we include it to preserve Hermiticity exactly. Symmetrization is a standard way to quantize a product of two noncommuting observables and the minimal procedure to ensure Hermiticity. The situation here is formally analogous to quantizing Hamiltonian with a spatially varying Rashba spin-orbit coupling. Numerical comparisons show that the effective Hamiltonian (II.3) produces the low-energy LDOS of the exact four-band model accurately, as seen in Figs. 2c and d. The linear approximation also illuminates the behavior of the exact model. The flat part of the LDOS in the vicinity of the Weyl node appears due to a chiral mode appearing in the gap of the bulk Weyl cones because of the Landau quantization induced by the local effective magnetic field $\tilde{{\bf B}}=\nabla\times{\bf k}^{W}/q$. The inset of Fig. 2c shows how a constant field of magnitude $|\tilde{{\bf B}}|$ give rise to a chiral mode representing the zeroth Landau level dispersing in the direction the field. The chiral mode has a 1d dispersion and a constant density of states. This interpretation is confirmed by the fact that the calculated LDOS is approximately constant in the energy window $2E_{\tilde{B}}=2v_{\tilde{B}}\sqrt{2q|\tilde{{\bf B}}|}$ which correspond to the gap between $-1$ and $1$ Landau levels where $v_{\tilde{B}}$ is the Fermi velocity along the magnetic field. In addition, as shown in Fig. 2d, for textures which the effective magnetic field vanishes, the LDOS displays a quadratic trend down to the Weyl node.

III Particle motion in Weyl metamaterials

III.1 Semiclassical dynamics

To study the quasiparticle dynamics, we employ the semiclassical theory xiao . We assume that the system is doped above the Weyl node and that the Landau quantization due to the real and effective gauge fields can be neglected. Then we can ignore interband dynamics and concentrate on a single positive energy band. The semiclassical wave-packet motion for inhomogeneous systems xiao is described by the equations of motion

[TABLE]

where $\mathcal{E}=(g^{ij}k_{i}k_{j})^{\frac{1}{2}}$ is the relativistic curved-space dispersion and $E$ and $B$ are electric and magnetic fields including the emergent field $\tilde{B}$. Quantum particles have orbital degrees of freedom encoded in the frame fields, leading to quantum-geometric effects given by the Berry forces on top of the classical effects. The novelty compared to the previously extensively studied flat case is that, in addition to the curved-space dispersion, we must consider a $6\times 6$ phase-space Berry-curvature tensor $\Omega$ instead of its $3\times 3$ momentum block. In the studied model, the Berry tensor has components $\Omega_{q_{1}q_{2}}=-\frac{\varepsilon^{ijk}}{2|d|^{3}}d_{i}\partial_{q_{1}}d_{j}\partial_{q_{2}}d_{k}$, where $d_{i}$ are the coefficients of the Pauli matrices in Eq. (3). In the linear approximation, the Berry curvatures can be expressed in terms of the frame fields and their derivatives by substituting $d_{j}=\tensor{e}{{}^{i}_{j}}k_{i}$. In the lowest non-trivial order, we can drop the $\Omega_{r^{i}r^{l}}$ term since it is the only term that involves two spatial derivatives and is parametrically small. Away from the Weyl nodes, the Berry curvatures vanish as $\Omega_{k_{i}k_{l}}\sim\frac{1}{|k|^{2}}$, $\Omega_{r^{i}k_{l}}\sim\frac{1}{|k|}$, so by increasing the doping, one recovers the classical geometric effects together with the emergent magnetic field. However, in general, both classical and quantum geometric effects need to be considered together. On the other hand, we have shown in App. C that in the absence of external electromagnetic fields $E=0,\ B=0$, the ballistic motion is governed by the geodesic equation for a charged particle $\ddot{r}^{l}+\tensor{\Gamma}{{}^{l}_{i}{}_{j}}\dot{r}^{i}\dot{r}^{j}=qg^{li}\varepsilon_{ijk}\dot{r}^{j}\tilde{B}^{k}/\mathcal{E}$ in the lowest nontrivial order in spatial derivatives. Here $\tensor{\Gamma}{{}^{l}_{i}{}_{j}}$ is the Christoffel symbol calculated from the metric $g_{ij}$. The geodesic motion is a direct manifestation of the curved geometry and indicates how the carrier motion can be tailored through geometry. The carrier dynamics in disordered solids also include scattering processes, which could be implemented by incorporating the semiclassical dynamics into a Boltzmann equation. However, this goes beyond the scope of the present work where we illustrate the designer geometries by solving Eq. (III.1).

III.2 3D Weyl electron lens

The negative index of refraction in optical metamaterials enable fabrication of exotic devices leonhardt2 such as the Veselago lens veselago which can also be realized for electrons in graphene cheianov . Inspired by these ideas, we propose a 3d Weyl electron lens as an application of Weyl metamaterials. This phenomenon can be realized, for example, in a two-node model $H=v_{i}k_{i}\gamma_{i}+m\gamma_{4}+{\bf u}\cdot{\bf b}$, which could represent a topological insulator-ferromagnet sandwich structure or magnetically-doped topological insulator with slowly rotating magnetization ${\bf u}=u\left[\cos\varphi\sin\theta(r),\sin\varphi\sin\theta(r),\cos\theta(r)\right]$ expressed in spherical coordinates. The texture is parametrized by the function $\theta(r)=\omega r$ and depicted in Fig. 3 a which illustrates the carrier trajectories in the setup. The carriers of one chirality entering the sample in the vicinity of $r=0$ will converge periodically in the vicinity of focal points. This motion results from a combined effect of the curved geometry and the effective magnetic field that can be straightforwardly calculated in our general framework. In Fig. 3 b we have solved the trajectories from the semiclassical equations corresponding to initial conditions where parallel beams enter the sample at $z=0$. Thus, the image of the object at infinity is reproduced at focal points. The situation in Fig. 3 c corresponds to the case where a point-like object is located at $z=0$ and imagined at the focal points. The carriers with the opposite chirality experience an opposite magnetic field and diverge away from $x=y=0$ axis as shown in Fig. 3 c. Thus the considered geometry will realize a chirality-selective 3d electron lens. This effect is similar to the chirality-dependent Hall effect discovered recently in inhomogeneous Weyl semimetals yang2 .

It can be shown analytically (see App. D) that for small-amplitude oscillations, assuming $\sqrt{k_{x}^{2}+k_{y}^{2}}/k_{z}\ll 1$, $v_{x}=v_{y}=v$, the distance between focal points $\Delta z$ is

[TABLE]

if we neglect the effects of the Berry curvature – the exact numerical solution shows that the motion is not completely planar due to quantum geometric effects. In the absence of the effective magnetic field, the deviation $\delta l$ from the plane when the particle traverses distance $\Delta z$ along the $z$ axis can be shown from the Berry curvatures together with Eqs. (III.1) to scale as $\delta l/\Delta z\sim\omega k_{\parallel}/k_{z}^{2}$, where $k_{\parallel}$ stands for the momentum within the plane, perpendicular to the $z$ axis. In the presence of the effective magnetic field, this scaling relation is less exact in its momentum dependence, but nevertheless provides the correct order of magnitude for the deviation. Thus, in the regime of interest it holds that $\delta l\ll\Delta z$ so the lens motion is nearly planar. The expression for $\Delta z$ shows how the length scale for the lensing effect is controlled by the rotation of the texture through $1/\omega$. An arbitrarily slowly varying texture would yield a similar lens effect with a longer separation of the focal points. The lens effect depends on the energy of the incoming particles only weakly through $k$ in the denominator of the square root.

Interestingly, a similar 3d lens effect can also be realized in unidirectional magnetic textures as discussed in App. E. We also discovered a 2d lens effect for a simpler magnetic texture that only rotates in plane. As smoking-gun observable evidence of the lens phenomenon one could consider the quantum mirage effect pointed out in Ref. cheianov : a point-like impurity will in general give rise to localized Friedel oscillations whose image due to the lens geometry can be observed away from the impurity, giving rise to a mirage perturbation. For example, a density variation due to an impurity located at the origin in Fig. 3 c would have its image appear at the convergence of the trajectories near $(0,0,11)$. The 3d Weyl lens illustrates that Weyl metamaterials posses richness beyond the 2d materials considered for curvature and gauge field engineering. Interesting future challenges include, for example, finding realizations for 3d electronic invisibility devices.

In both the numerical and analytical analyses of the 3D lens geometry, we have not considered surface effects that might occur when injecting wave packets through sharply defined boundaries. To study sharp interface effects, one should complement the present theory of smooth TRIB fields with additional boundary conditions. This will be addressed in future work. Nevertheless, the TRIB-field configuration will result in oscillatory behavior for one chirality of the incident particles once they are inside the metamaterial.

IV Discussion

Experimental realizations of the proposed Weyl metamaterial structures require controlled fabrication of systems with inhomogeneous TR or I breaking. Weyl-metamaterial structures can be fabricated from topological insulators and topological semimetals – both of which have an ever-growing number of experimental realizations. In App. B, we have outlined how four much-studied models where uniform TR and I breaking has been promoted to inhomogeneous can serve as platforms for Weyl metamaterials. Local breaking of inversion through strain is already considered feasible grushin ; pikulin and, considering the relentless experimental advances in the field, local magnetic manipulation will become accessible in the near future. For example, artificially constructed layered Weyl semimetals burkov2 are now becoming experimentally realizable hasan ; artificially grown layer structures allow a degree of manipulation of the magnetization in the magnetic layers. Another interesting option for realizing varying magnetization arises from the natural magnetization dynamics in magnetically doped topological semimetals. These systems are predicted to display textured ground states or excitations araki which could be employed in metamaterial structures.

Now we turn to the question of experimental requirements important in realizing functional Weyl metamaterial structures. Considering the relatively mild requirements for physical parameters and the rapid advance of experimental methods, we are optimistic that the first Weyl metamaterial structures can be fabricated in the near future. It is promising that the required rate of variation of the symmetry breaking fields do not need to exceed specific critical value. The characteristic scale of variation simply fixes the overall size of functional units. For example, reducing the rate of rotation of the magnetization in the 3d lens geometry increases the focal length proportionally. To be more precise and to give simple estimates to guide the experimental realizations, it is illuminating to consider a two-node example such as the Weyl lens Hamiltonian. In order to observe the effects of the curvature and the emergent gauge field, their length scales must be comparable or smaller than the length scale of the sample, i.e. for a sample with linear dimensions $L$, we have $|\partial_{i}{\bf u}|L/u\gtrsim 1$, where $u$ and $|\partial_{i}{\bf u}|$ are the characteristic magnitudes of the symmetry-breaking field ${\bf u}$ and its gradient in the sample volume. On the other hand, the Weyl description breaks down if the spatial variation of the magnetic texture is too large, enabling scattering between different Weyl nodes. For Weyl nodes separated by a distance $\Delta k^{W}$, the spatial variation must hence satisfy $|\partial_{i}{\bf u}|/(\Delta k^{W}u)\ll 1$. The two requirements can in general be satisfied in systems where $\Delta k^{W}\gg L^{-1}$, which is not a serious restriction for Weyl metamaterial structures. In the particular case of the Weyl lens, we can estimate that $|\partial_{i}{\bf u}|/u=\omega$, leading to constrains $\omega L\gtrsim 1$ and $\omega\ll v_{F}\sqrt{u^{2}-m^{2}}$. We see that the requirements can be satisfied in a sufficiently large sample as long as the system is in the Weyl semimetal phase ($u>m$). Of course, in realistic systems there could be other scales affecting the performance such as the momentum relaxation length that need to be taken into account. Nevertheless, there are no fundamental restrictions for the size of the gradients of the magnetic textures and inversion breaking other than being slow enough as to not allow for scattering between Weyl nodes. Therefore it seems plausible that simple metamaterial structures could be realized in the near future.

Furthermore, the discussion in the previous paragraph points to an important feature of Weyl semimetals not only restricted to the lens geometry. Namely, the curved-space effects in Weyl metamaterials follow almost exclusively from the gradients of the TRIB fields, while the role of the absolute strengths of the fields is only minor in this regard. Certainly, in order for the Weyl-metamaterial picture to be valid, the local field strengths must be such that a corresponding homogeneous sample would be in a Weyl semimetal phase. So any restrictions put on the field strengths in a Weyl metamaterial are the same restrictions as for any corresponding Weyl semimetal in the same sample. Put slightly differently: any physical realization of a Weyl semimetal in an electronic system should in principle also be tunable into a Weyl metamaterial.

V conclusion

In this work, we proposed Weyl metamaterials as a highly tunable platform for relativistic fermions in curved space. Analogous to optical metamaterials, we show that Weyl metamaterials offer a 3D electronic platform where local manipulations of TR and I breaking fields allow efficient control over the particle propagation. We laid the ground work for the theoretical description of these systems and established the explicit connection between the local TR and I breaking, and the effective geometry and gauge fields experienced by the carriers.

From the point of view of applications, Weyl metamaterials offer a novel route to 3D electronic devices. The functionality of the devices arise from the curved-space geometry and emergent gauge fields and they can be systematically designed with the help of our theory framework. In the present work, we provided a concrete example, namely the 3D Weyl electron lens. This structure, realized by relatively simple magnetic textures, focuses carriers depending on their chirality. This is just the first example of the new possibilities offered by Weyl metamaterials. An interesting venue for future work is to study what other exotic phenomena can be engineered in these systems. For example, a possibility to realize electronic cloaking devices from Weyl metamaterials is particularly intriguing.

From the fundamental point of view, our results can be employed in systematic reverse-engineering of 3D curved geometries and studying quantum effects birrell in designer geometries. One direction is then to study the curved-space modifications to the exotic response properties such as the chiral anomaly and its manifestations in Weyl semimetals liu2 ; huang ; xiong ; li ; zhang2 ; behrends . Also, fabrication of condensed-matter analogues of cosmological horizons could bring new fundamental phenomena within the reach of experimental studies. Considering the rapid developments in materials and realizations of topological semimetals, the different aspects of curved-space quantum physics will likely become experimentally accessible soon.

Acknowledgements.

The authors acknowledge the Academy of Finland for support.

Appendix A Weyl block reduction

A.1 TR breaking without I breaking

In this supplement, we derive the unitary transformations that reduce the $4\times 4$ Hamiltonians for different models into blocks of ($2\times 2$) Weyl Hamiltonians.

It turns out that it is beneficial to first consider the Hamiltonian in Eq. (2) of the main text with ${\bf w}=0$, i.e.

[TABLE]

As can be easily shown by inspection, the matrices in ${\bf b}$ have identical properties to those of the regular $2\times 2$ Pauli matrices, i.e.

[TABLE]

These matrices, like the Pauli matrices, thus exponentiates to the special unitary group SU(2), which is two-to-one homomorphic to the orthogonal group O(3). In practice, this means that the transformation

[TABLE]

corresponds to rotating the vector ${\bf u}$ by an angle $\varphi$ around the axis defined by the unit vector $\hat{{\bf n}}$. Furthermore, by observing that the components of the other vectors in table 1 in the main text along with $\boldsymbol{\gamma}$ can be written as ${\bf p}=-\gamma_{5}{\bf b}$, ${\bf b}^{\prime}=\gamma_{4}{\bf b}$, and $\boldsymbol{\gamma}=-\gamma_{45}{\bf b}$, respectively, we can conclude that the exponential map has the same effect (it commutes with $\gamma_{4}$ and $\gamma_{5}$) on terms containing these vectors as well.

The first step is to rotate ${\bf u}$ in Eq. (12) to point along the $x$ axis. This can for example be achieved by first rotating ${\bf u}$ into the $xz$ plane and then rotating it around the $y$ axis, aligning it with the $x$ axis. These two transformations are effected by

[TABLE]

which is expressed in terms of the polar ($\theta$) and azimuthal ($\varphi$) angles of ${\bf u}$. Applying this transformation to the whole Hamiltonian yields

[TABLE]

To block-diagonalize this, we need to eliminate either $\gamma_{1}$ or $\gamma_{4}$. This is achieved using $\mathcal{R}=e^{-i\frac{\rho}{2}\gamma_{14}}$ with $\tan\rho=(\boldsymbol{\kappa}({\bf k})\cdot\hat{{\bf u}})/m$. The final result is then

[TABLE]

which splits into blocks in the basis given by for example

[TABLE]

With this choice, the ${\bf d}$ vectors for the upper and lower block become

[TABLE]

Depending on the sign of $m({\bf k})$, the nodes will either be located in the upper or the lower block. The ${\bf d}$ vector corresponding to the Weyl block will regardless have the form

[TABLE]

where the potential minus sign on the first and second component have been taken care of using a unitary transformation $\sigma_{z}$. The above vector vanishes whenever all three components are zero, i.e.

[TABLE]

Furthermore, note that if ${\bf d}^{W}({\bf k}^{W})=0$, it follows immediately that ${\bf d}^{W}(-{\bf k}^{W})=0$, which is just a manifestation of the fact that Weyl nodes come in pairs.

To arrive at an effective curved-space Weyl Hamiltonian, we linearize ${\bf d}^{W}$ around the Weyl nodes. The Weyl Hamiltonian can then be written in terms of frame fields as

[TABLE]

where the frame fields $\tensor{e}{{}^{j}_{i}}$ can be found to be

[TABLE]

where the sign in the last component corresponds to either ${\bf k}^{W}$(upper sign) or $-{\bf k}^{W}$(lower sign). The co-frame fields $\tilde{e}$ are easy to evaluate whenever $\tensor{e}{{}^{j}_{0}}=\tensor{e}{{}^{0}_{j}}=0$, and are given by

[TABLE]

Finally, the metric tensor can be straightforwardly calculated to be

[TABLE]

where we have have omitted writing out the arguments. The whole expression must then be evaluated at the relevant Weyl node.

A.2 I breaking with preserved TR

In this section we show how to obtain similar results for a system with only inversion-symmetry breaking. The starting Hamiltonian is

[TABLE]

The last term can be rotated using $\mathcal{U}({\bf w})$ to give

[TABLE]

where ${\bf K}$ can be directly extracted from Eq. (16). On this we can then simply apply a rotation $\exp(i\nu\gamma_{23}/2)$ with $\tan\nu=K_{3}/K_{2}$ to get

[TABLE]

Employing, for example, the representation

[TABLE]

we see that this Hamiltonian is block diagonal with blocks given by

[TABLE]

It follows that the Weyl block is

[TABLE]

where, if $\operatorname{sgn}(K_{2}({\bf k}^{W}))<0$, we apply a transformation $\sigma_{z}$ to get rid of the minuses in front of $K_{1}=\hat{{\bf w}}\cdot\boldsymbol{\kappa}({\bf k})$ and $m$. The condition for the existence of the Weyl nodes ${\bf k}^{W}$ become

[TABLE]

In this case, the minimum number of nodes is four. We can now linearize this around the Weyl points, yielding

[TABLE]

We can now immediately read off the frame fields to be

[TABLE]

and the metric tensor

[TABLE]

where all terms should be evaluated at ${\bf q}={\bf k}^{W}$.

A.3 Block reduction with TR and I breaking

The problem becomes significantly more cumbersome when both ${\bf u}$ and ${\bf w}$ are non-vanishing. In this case, it is useful to write ${\bf w}={\bf w}_{\perp}+{\bf w}_{\parallel}$, where ${\bf w}_{\perp}\cdot{\bf u}=0$ and ${\bf w}_{\parallel}\parallel{\bf u}$. The goal is to find a transformation that removes all matrices belonging to ${\bf p}$ from the Hamiltonian, such that we are left with the special case discussed in the first part of this appendix.

Before tackling the full problem, first consider the case where ${\bf w}_{\perp}=0$. Then it follows that $\mathcal{U}({\bf u})$ also aligns ${\bf w}$ with the $x$ axis. The Hamiltonian is then

[TABLE]

where ${\bf K}$ can be directly extracted from Eq. (16). We can now remove $\gamma_{14}$ by a simple unitary transformation $\mathcal{F}=e^{i\delta\gamma_{1}/2}$, provided that $\tan\delta=w/m$:

[TABLE]

Evidently, $\mathcal{F}$ also introduces new ${\bf b}$ terms, but at this point, we can simply refer back to the first half of this appendix, as the problem is mathematically identical to the case ${\bf w}=0$.

Solving the full problem can thus be done if we can find a transformation that aligns the vector associated with ${\bf b}$ and the vector associated with ${\bf p}$. For the sake of clarity, we will not concern ourselves with that of $H_{0}$ until at the very end (as it turns out, all transformations merely rotates $\boldsymbol{\kappa}$), so our primary focus will now lie on

[TABLE]

As usual, we begin by applying $\mathcal{U}({\bf u})$, which gives us

[TABLE]

We can then make a rotation around $b_{1}$ with the unitary transformation $\mathcal{V}=e^{i\mu b_{1}/2}$ which for $\tan\mu=\tilde{w}_{\perp,3}/\tilde{w}_{\perp,2}$ eliminates the $p_{3}$ term:

[TABLE]

At this point, we make two transformations in succession; we first apply $\mathcal{V}_{b}=e^{i\alpha b_{3}/2}$ and then $\mathcal{V}_{p}=e^{i\beta p_{3}/2}$ which gives a seemingly unwieldy

[TABLE]

where $w_{\parallel}^{\prime}=\operatorname{sgn}({\bf w}\cdot{\bf u})w_{\parallel}$. We can now determine the angles $\alpha$ and $\beta$ by requiring that the terms proportional to $b_{2}$ and $p_{2}$ vanish. One possible solution is then

[TABLE]

From this then follows that

[TABLE]

to which we can apply our earlier methods.

The effects of these transformations on $H_{0}$ are then straightforwardly calculated. As usual, $\mathcal{U}({\bf u})$ gives us

[TABLE]

which upon application of $\mathcal{V}$ turns into

[TABLE]

With the identifications

[TABLE]

the Hamiltonian can be written as

[TABLE]

from where we can refer back to the case with ${\bf w}\parallel{\bf u}$, which combined with our earlier discussion on parallel ${\bf w}$ and ${\bf u}$ gives us

[TABLE]

on which we can apply the machinery developed for the case of ${\bf w}=0$ (Eq. (12)).

As a final comment on the general case, we note that the equations for the Weyl nodes now take the form

[TABLE]

where

[TABLE]

and,

[TABLE]

As two straightforward examples, we see that if ${\bf w}\parallel{\bf u}$, Eq. (58) gives us

[TABLE]

while for ${\bf w}\perp{\bf u}$ it yields

[TABLE]

Both of these cases reduce to the correct equation in the ${\bf w}=0$ limit.

Appendix B Realizations of Weyl metamaterials

To illustrate the utility of our theoretical framework and the physics of Weyl metamaterials, we consider three much-studied models as platforms for Weyl metamaterials. These include the topological insulator-ferromagnet layer structure [31], the popular toy model introduced by Vazifeh and Franz [29], and a Dirac semimetal model of Cd3As2 [9] with broken inversion symmetry. As detailed in Ref. [31], the topological insulator-ferromagnet heterostructure can be modelled by the Hamiltonian

[TABLE]

where $v_{F}$ is the Fermi velocity, ${\bf k}$ the momentum, $\hat{\Delta}(k_{z})=\Delta_{S}\tau^{x}+\Delta_{D}(\cos k_{z}\tau^{x}-\sin k_{z}\tau^{y})$, and $m$ is a time-reversal breaking magnetization. The parameters $\Delta_{S}$ and $\Delta_{D}$ characterize the tunnelling between different surfaces within and between neighboring layers of the heterostructure, respectively. While not unique, a suitable choice for the representation of the Dirac gamma matrices allows us to write Eq. (66) as

[TABLE]

where we have promoted the model from having magnetization only in the $z$ direction $mb_{3}\to{\bf m}\cdot{\bf b}$. Here we see that ${\bf m}$ corresponds to the TR-breaking field ${\bf u}$ in Table 1 of the main text and letting it vary as a function of position, we obtain an example of a TR-breaking metamaterial. By applying our general framework, it is straightforward to study the curved-space low-energy physics of the model.

We next consider the toy model introduced in [29]. The model represents a 3d (topological) insulator with a finite magnetization and has a Hamiltonian

[TABLE]

where we have slightly deviated from the notation in [29] for the TR and I breaking terms. Here, the mass term is $M_{{\bf k}}=\mu-2t\sum_{\alpha=1}^{3}\cos k_{\alpha}$, where $\mu$ and $t$ are the chemical potential and hopping amplitude, respectively. The parameters $\lambda,\lambda_{z}$ characterize the strength of the spin-orbit hopping and ${\bf u}$ can be identified with magnetization. Expressing this Hamiltonian via Dirac gamma matrices gives us

[TABLE]

Within the scope of our present formalism, we cannot account for the I-breaking $\varepsilon$ term and must hence restrict ourselves to $u_{0}=0$, but otherwise this is now in a form which we can apply our formalism to, starting from Eq. (12). For example, we immediately get that the equations for the Weyl nodes are given by

[TABLE]

where $\lambda_{1,2}=\lambda$ and $\lambda_{3}=\lambda_{z}$. The frame fields and the metric can be now extracted from our general formulas.

Our third example is a 3D Dirac semimetal [9] with broken I symmetry. The Hamiltonian describes the spectrum of Cd3As2 (or Na3Bi) near the $\Gamma$ point and is given by

[TABLE]

where $\varepsilon_{0}({\bf k})=C_{0}+C_{1}k_{z}^{2}+C_{2}(k_{x}^{2}+k_{y}^{2})$ and $M({\bf k})=M_{0}+M_{1}k_{z}^{2}+M_{2}(k_{x}^{2}+k_{y}^{2})$ [9]. The above Hamiltonian can be directly translated into our gamma-matrix representation via the substitution $\tau_{z}\sigma_{j}=\gamma_{j}$ for $j=1,2$ and $\tau_{z}\sigma_{z}=\gamma_{4}$. We can add an inversion-breaking field to this model by introducing a coupling ${\bf w}\cdot{\bf p}$, such that we get

[TABLE]

While the above form with the I-breaking vector parameter ${\bf w}$ is a bit abstract, it is possible to implement it by, for example, elastic deformations due to strain [22-25]. Starting from Eq. (29), we can then then derive the desired quantities, viz. the frame fields and the induced metric.

Appendix C Semiclassical motion and the geodesic equation

The geodesic equation can be derived from the semiclassical equations of motion in Eq. (10) of the main text. In the absence of electromagnetic fields and Berry curvatures, the equations reduce to

[TABLE]

which – given that $\mathcal{E}$ is a constant of motion – can by a rescaling of the time $t\to\varepsilon t^{\prime}$ be brought into the form

[TABLE]

We recognize these as the Hamilton equations of motion, which follow from extremizing the action

[TABLE]

with the Lagrangian given by

[TABLE]

This problem is equivalent to extremizing the length of a path in a curved space given by the metric $g^{ij}$, and it immediately follows that ${\bf r}$ obeys the geodesic equation

[TABLE]

with the connection $\Gamma$ given by

[TABLE]

Reintroducing the Berry curvatures and the emergent magnetic field changes this picture: although the energy remains a constant of motion up to first order in spatial derivatives, the geodesic equation now becomes that of a charged particle in an external magnetic field. We have that

[TABLE]

If we neglect terms that are second order or higher in spatial derivatives, we can remove $\Omega_{{\bf r}{\bf r}}$ and write for example $(\mathbb{I}-\Omega_{{\bf r}{\bf k}})^{-1}\approx\mathbb{I}+\Omega_{{\bf r}{\bf k}}$. The equations reduce to

[TABLE]

By rewriting the first equation in component form and differentiating once more with respect to time, we have

[TABLE]

which follows from using the fact that $\mathcal{E}=\sqrt{g_{ij}k_{i}k_{j}}$, and dropping terms with more than one spatial derivative. To show that this is equivalent to a geodesic equation to linear order, note that the left term on the right-hand side is already linear in spatial derivatives, and we can thus simply substitute ${\bf k}=\mathcal{E}\boldsymbol{g}^{-1}\dot{{\bf r}}+\mathcal{O}(\partial_{{\bf r}})$, which follows directly from Eq. (80) and once again dropping all higher-order spatial derivatives. Hence, Eq. (C) now becomes

[TABLE]

where

[TABLE]

is the Christoffel symbol.

As a final remark, we point out that the above derivation used the fact that

[TABLE]

where the last equality follows from the antisymmetry of $\varepsilon_{jmn}$.

Appendix D Analytical solution for 3D lens geometry

In this appendix we present an analytical solution for the oscillatory behavior of the trajectories in the limit of small deviations from $x^{2}+y^{2}=0$ and $k_{x}/k_{z},k_{y}/k_{z}\ll 1$. The spatial part of the metric for the 3D lens is given by

[TABLE]

where the omitted elements that can be deduced from $\boldsymbol{g}$ being symmetric. For small $\theta$, the metric reduces to

[TABLE]

From Eq. (86), it is straightforward to calculate its spatial derivatives which in our approximation reduce to

[TABLE]

and $\partial_{z}\boldsymbol{g}=0$.

It follows from writing down the equations of motions using Eqs. (10) in the main text using these approximations that the system of first order differential equations separate in $x$ and $y$, such that we have

[TABLE]

for both of them. Here, we have ignored the contributions from the the Berry curvatures (which would couple the motion in the $x$ and $y$ directions) as their effects are small on short enough time scales. We have additionally taken into account the effective magnetic field arising from a position-dependent Weyl node, $\tilde{{\bf B}}=\nabla\times{\bf k}^{W}/q$. Since this magnetic field acts with different signs on the opposite nodes, we have introduced the parameter $c=\pm 1$ to distinguish the Weyl nodes ($c=0$ would correspond to pure geodesic motion). In the derivation of this linear differential equation, we have only kept terms up to linear order in $k_{x}$ and $k_{y}$ and used the fact that $\dot{k_{z}}=0$. The solution to this differential equation can be formally written as $\Psi(t)=\exp(\boldsymbol{A}t)\Psi(t=0)$, where $\boldsymbol{A}$ is the matrix on the right-hand side. From this, it is evident that if $\boldsymbol{A}$ has imaginary eigenvalues, we will obtain oscillatory solutions. To simplify, we now assume that $v_{x}=v_{y}\equiv v$. The eigenvalues are then

[TABLE]

from which it follows that we have sinusoidal solutions typically only for one chirality $c=1$. Furthermore, since $\dot{z}\approx v_{z}^{2}(1-m^{2}/u^{2})k_{z}/\mathcal{E}$, we have that the length of one oscillation with period $T=2\pi/|\lambda|$ is

[TABLE]

The distance between two focal points is $\Delta z/2$ which depends on the energy of particles weakly through $k_{z}$ in the square root. The analytical predictions are confirmed by exact numerics. The Berry-curvature effects that are neglected here affect the long time dynamics, most clearly by resulting in a non-planar motion that is not evident for short times. The analytical results are applicable also for different textures when $\omega$ is replaced by a characteristic measure of the variation of $u_{z}$ component near the lens axis.

Appendix E Alternative lens realizations

In the main text, we considered an electron lens in a two-node model $H=v_{i}k_{i}\gamma_{i}+m\gamma_{4}+{\bf u}\cdot{\bf b}$ with a skyrmion-like magnetic texture. This is not necessary for the lens effect to occur as we can see in Fig.  4, which illustrates how a similar effect can be realized by a unidirectional texture ${\bf u}=\left[0,0,u(r)\right]$. If $u(r)$ decreases radially away from the lens axis as in Fig. 4 a, the carrier trajectories near the axis are qualitatively similar to the rotating texture. However, even simpler textures result in nontrivial lensing; the case where $u(r)$ varies only in the $x$ direction, as shown in Fig. 4 b, still leads to a 2d lens effect. It is encouraging from a practical point of view that no elaborate 3d textures are needed in realizing remarkable effects. In Figs. 4 c and d, we have plotted the particle trajectories corresponding to a texture like in b. We have used a Gaussian envelope although the precise form is not crucial for the qualitative features. The motion in the $x-z$ plane looks qualitatively similar to the 3d-lens motion but only in that plane – not for any arbitrary plane parallel to the $z$-axis as in a 3d lens.

Appendix F Numerical calculation of the local density of states

We employed the LDOS as a tool for comparing the exact four-band model and the linearized two-band Weyl block for inhomogeneous TR-breaking textures ${\bf u}({\bf r})$. This was done by considering a system that is inhomogeneous in the $z$-direction but homogeneous with periodic boundary conditions in the $x-y$ plane. The considered model has a Hamiltonian $H=\boldsymbol{\kappa}({\bf k})\cdot\boldsymbol{\gamma}+m\gamma_{4}+{\bf u}\cdot{\bf b}$, where $\kappa_{i}({\bf k})=t\sin{k_{i}}$ and $m$ and $t$ are constants. In the numerics, we employed the representation $\gamma_{i}=\tau_{z}\sigma_{i}$ for $i=1,2,3$ and $\gamma_{4}=\tau_{x}$. The LDOS was studied for rotating and linear textures ${\bf u}(z)$. In the studied geometry, $k_{x}$ and $k_{y}$ are good quantum numbers so the LDOS can be expressed as $\nu(z,E)=\sum_{E_{k_{x},k_{y},n}}|\Psi_{k_{x},k_{y},n}(z)|^{2}\delta(E-E_{k_{x},k_{y},n})$, where $n$ denotes the discrete quantum number in $z$ direction. We considered a system on lattice with $N$ lattice sites in the $z$ direction, replaced $t\sin{k_{z}}$ with a standard hopping matrix, and imposed hard-wall boundary conditions. In the continuum model, the momentum operator $\partial_{z}$ was replaced by a hopping matrix representing a discretized derivative. In order to evaluate the LDOS, we numerically calculated $E_{k_{x},k_{y},n}$ and $\Psi_{k_{x},k_{y},n}(z)$ for the 1d tight-binding models with a sufficiently dense $(k_{x},k_{y})$-lattice for the four-band model and the two-band models. To obtain continuous functions in energy, we replaced the energy delta functions by Lorentzians $\delta(E-E_{n})\to\frac{1}{(E-E_{n})^{2}+\eta^{2}}$ with broadening $\eta$.