Topological Andreev Rectification

TL;DR

This paper develops a theory for topological Andreev rectification in Josephson junctions with 2DEGs, showing that a robust, quantized conductance signal can reveal the Fermi sea's topology even in realistic device geometries.

Contribution

It extends previous narrow junction models to realistic geometries, demonstrating the robustness of the topological conductance signal under certain scattering conditions.

Findings

Rectified conductance remains quantized with suppressed large momentum scattering.

Wider junctions and point contacts still exhibit robust conductance signals.

Simulations suggest experimental observation is feasible in ballistic InAs and graphene devices.

Abstract

We develop the theory of an Andreev junction, which provides a method to probe the intrinsic topology of the Fermi sea of a two-dimensional electron gas (2DEG). An Andreev junction is a Josephson $\pi$ junction proximitizing a ballistic 2DEG, and exhibits low-energy Andreev bound states that propagate $\textit{along}$ the junction. It has been shown that measuring the nonlocal Landauer conductance due to these Andreev modes in a narrow linear junction leads to a topological Andreev rectification (TAR) effect characterized by a quantized conductance that is sensitive to the Euler characteristic $\chi_F$ of the 2DEG Fermi sea. Here we expand on that analysis and consider more realistic device geometries that go beyond the narrow linear junction and fully adiabatic limits considered earlier. Wider junctions exhibit additional Andreev modes that contribute to the transport and degrade the quantization of the conductance. Nonetheless, we show that an appropriately defined $\textit{rectified conductance}$ remains robustly quantized provided large momentum scattering is suppressed. We verify and demonstrate these predictions by performing extensive numerical simulations of realistic device geometries. We introduce a simple model system that demonstrates the robustness of the rectified conductance for wide linear junctions as well as point contacts, even when the nonlocal conductance is not quantized. Motivated by recent experimental advances, we model devices in specific materials, including InAs quantum wells, as well as monolayer and bilayer graphene. These studies indicate that for sufficiently ballistic samples observation of the TAR effect should be within experimental reach.

Tables

Table 1. Table 1: Parameter values (in meV) used for simulating an Andreev junction in a hybrid InAs/Al system.

$t$	${\lambda }_R$	${\Delta }_0$	${\mu }_N$	${\mu }_S$
$15.0$	$0.75$	$0.15$	$10.0$	$10.0$, $10.3$

Table 2. Table 2: Parameter values for simulating an Andreev junction in hybrid graphene–superconductor systems.

	$t$	$t_{⟂}$	${\gamma }_3$	${\Delta }_0$
MLG	$2.8\text{eV}$	n.a.	n.a.	$0.2\text{meV}$ (Al)
				$1\text{meV}$ (Nb or MoRe)
BLG	$2.8\text{eV}$	$0.1t$	$0.1t$	$0.0025t$

Equations

$$G_{21}(V)\equiv {\frac{d{I}_2}{d{V}_1}|}_{V_1=V}=\{\begin{array}{ll}{c}_e{e}^2/h& \text{for}V<0,\\ c_h{e}^2/h& \text{for}V>0.\end{array}$$

$$c_e-c_h={\chi }_F,$$

$$\delta G_{21}(V)\equiv G_{21}(V)-G_{21}(-V).$$

$${\chi }_F=\sum _l(-1{)}^l{b}_l,$$

$${\chi }_F=n_0-n_1+n_2,$$

$${\chi }_F=\sum _k(2-2g_k-b_k),$$

$${\chi }_F=n_e-n_h,$$

$${\chi }_F=c_e-c_h,$$

$$\mathcal{H}=-i{\tau }_z{\sigma }_z{v}_y{\partial }_y+{\Delta }_1(y){\tau }_x+{\Delta }_2(y){\tau }_y,$$

$$\begin{array}{rl}{H}_{\text{BdG}}(k_x)=& \left[\theta \left(y-W/2\right)-\theta \left(-y-W/2\right)\right]{\Delta }_0{\tau }_y\\ & +\left(-\frac{1}{2m^{\ast }}{\partial }_y^2+v_x\delta k_x\right){\tau }_z.\end{array}$$

$$\epsilon (k_x)=\pm \left[\frac{v_x\delta k_x}{2}+\mathrm{sgn}(m^{\ast })\sqrt{{\left(\frac{v_x\delta k_x}{2}\right)}^2+\frac{{\Delta }_0^2}{2}}\right].$$

$$G_{21}\equiv \frac{d{I}_2}{d{V}_1}\stackrel{T\to 0}{=}\frac{e^2}{h}{\left[T_{21}^{ee}(\epsilon )-T_{21}^{he}(\epsilon )\right]}_{\epsilon =e{V}_1}.$$

$$T_{21}^{ee}(\epsilon >0)=c_e\text{and}{T}_{21}^{ee}(\epsilon <0)=c_h,$$

$$G_{21}=\frac{e^2}{h}\left[c_e\theta (e{V}_1)+c_h\theta (-e{V}_1)\right].$$

$$\delta G_{21}(V_1)=\mathrm{sgn}(e{V}_1)\frac{e^2}{h}{\chi }_F,$$

$$\left(\begin{array}{c}{b}_{1L}^{\text{lead}}\\ {\varphi }_R^{\text{ABS}}\end{array}\right)=S_1\left(\begin{array}{c}{a}_{1R}^{\text{lead}}\\ {\varphi }_L^{\text{ABS}}\end{array}\right),\left(\begin{array}{c}{b}_{2R}^{\text{lead}}\\ {\varphi }_L^{\text{ABS}}\end{array}\right)=S_2\left(\begin{array}{c}{a}_{1L}^{\text{lead}}\\ {\varphi }_R^{\text{ABS}}\end{array}\right),$$

$$S_{1}=\left(\begin{array}[]{c|c}r_{1L}&t_{1L}\\ \hline\cr t_{1R}&r_{1R}\end{array}\right),\qquad S_{2}=\left(\begin{array}[]{c|c}r_{2R}&t_{2R}\\ \hline\cr t_{2L}&r_{2L}\end{array}\right),$$

$$\left(\begin{array}{c}{b}_{1L}^{\text{lead}}\\ {\varphi }_R^{\text{ABS}}\\ b_{1L}^{\overline{\mathrm{lead}}}\end{array}\right)=S_1\left(\begin{array}{c}{a}_{1R}^{\text{lead}}\\ {\varphi }_L^{\text{ABS}}\\ a_{1R}^{\overline{\mathrm{lead}}}\end{array}\right),$$

$$S_{1}=\left(\begin{array}[]{c|c|c}r_{1L}&t_{1L}&A\\ \hline\cr t_{1R}&r_{1R}&\bm{0}_{N\times(M_{1}-N)}\\ \hline\cr B&\bm{0}_{(M_{1}-N)\times N}&\rho_{1L}\end{array}\right).$$

$$t_{1L}^{†}A=t_{1R}{B}^{†}=0.$$

$$S_{1}=\left(\begin{array}[]{c|c|c}r_{1L}&t_{1L}&0\\ \hline\cr t_{1R}&r_{1R}&0\\ \hline\cr 0&0&\rho_{1L}\end{array}\right).$$

$$s=s_{1}\odot s_{2}=\left(\begin{array}[]{cc|cc}r_{L}^{he}&t_{L}^{hh}&0&0\\ t_{R}^{ee}&r_{R}^{eh}&0&0\\ \hline\cr 0&0&r_{L}^{eh}&t_{L}^{ee}\\ 0&0&t_{R}^{hh}&r_{R}^{he}\\ \end{array}\right),$$

$$\left(\begin{array}{cc}{t}_R^{ee}& t_R^{eh}\\ t_R^{he}& t_R^{hh}\end{array}\right)$$

$$\left(\begin{array}{cc}{r}_L^{ee}& r_L^{eh}\\ r_L^{he}& t_L^{hh}\end{array}\right)$$

$$G_{21}(+V_1)$$

$$G_{21}(-V_1)$$

$${\left(r_L^{he}\right)}^{†}{r}_L^{he}+{\left(t_R^{ee}\right)}^{†}{t}_R^{ee}$$

$$r_L^{he}{\left(r_L^{he}\right)}^{†}+t_L^{hh}{\left(t_L^{hh}\right)}^{†}$$

$$\mathrm{Tr}\left[{\left(t_R^{ee}\right)}^{†}{t}_R^{ee}\right]-\mathrm{Tr}\left[{\left(t_L^{hh}\right)}^{†}{t}_L^{hh}\right]=N_e-N_h,$$

$$\delta G_{21}(V_1)=\frac{e^2}{h}{\left(N_e-N_h\right)}_{\epsilon =e{V}_1}.$$

$$\displaystyle\begin{split}\hat{H}_{0}&=\sum_{\sigma=\uparrow,\downarrow}\sum_{\bm{n}}\Bigg{[}\left(2t_{x}+2t_{y}-\mu\right)\hat{c}_{\bm{n},\sigma}^{\dagger}\hat{c}_{\bm{n},\sigma}\\ &-\left(t_{x}\hat{c}_{\bm{n}+\bm{e}_{x},\sigma}^{\dagger}\hat{c}_{\bm{n},\sigma}+t_{y}\hat{c}_{\bm{n}+\bm{e}_{y},\sigma}^{\dagger}\hat{c}_{\bm{n},\sigma}+\mathrm{h.c.}\right)\Bigg{]},\end{split}$$

$${\hat{H}}_{\Delta }$$

Full text

I Introduction

A powerful method for characterizing quantum many-body phases of matter is to identify quantized response functions that probe topological features of the phase. This type of analysis was initiated by the integer quantized Hall effect, which probes the Chern number characterizing the topology of a gapped two-dimensional (2D) electronic phase [1, 2]. A related type of quantized response occurs in a one-dimensional (1D) ballistic metal, where the Landauer conductance exhibits steps that are quantized in units of $e^{2}/h$ reflecting the number of occupied one-dimensional subbands [3, 4, 5]. This integer can be interpreted as the Euler characteristic of the Fermi sea, $\chi_{F}$, a topological invariant that counts the number of components of the filled Fermi sea of a one-dimensional metal. Though the quantized Landauer conductance is insensitive to weak electron-electron interactions, this type of quantized response is less robust than the quantized Hall effect because it requires the conducting channels to be perfectly transmitted and the reflected channels to be perfectly reflected. This will be approximately the case provided transport is sufficiently ballistic, the contacts are sufficiently reflectionless and the channel is sufficiently long. Despite these conditions, quantized Landauer transport has been observed in a variety of systems, including relatively short quantum point contacts [6, 7, 8, 9].

The feasibility of quantized Landauer transport in one dimension motivated the search for probes of Fermi sea topology in higher-dimensional metals. Frequency-dependent and time-domain nonlinear responses [10, 11, 12], as well as equal-time density correlations [13], have been shown to probe $\chi_{F}$, though unlike the Landauer conductance in one dimension, quantization of these quantities is only precise in the absence of interactions. It was argued that a measure of multipartite entanglement in an interacting Fermi liquid probes $\chi_{F}$ [13], but that is challenging to probe experimentally.

Recently, two of the authors introduced a method to probe $\chi_{F}$ of a two-dimensional electron gas (2DEG) that is robust in the presence of interactions [14]. The proposal involves proximitizing the 2DEG with a grounded superconductor to form a linear111This type of device architecture has also been called planar in the literature [15, 16, 17, 18, 19, 20]. Here we use the term linear to highlight the one-dimensional channel in the 2DEG between the superconductors. We will also consider a point-contact geometry, where the 1D channel is reduced to a point. SNS (S: superconductor; N: normal metal) Josephson junction with phase difference $\pi$. Then $\chi_{F}$ is probed by measuring the Landauer transport along the junction between source and drain contacts to the 2DEG at either end of the superconducting junction, as shown in Fig. 1. The current in the drain ($I_{2}$) in response to a voltage applied to the source ($V_{1}$) is carried by Andreev bound states (ABSs) that are extended along the junction. We will refer to this setup as an Andreev junction, so as to contrast with the more conventional Josephson junction where the current of interest flows between the superconductors.

It is instructive to compare transport in this Andreev junction setup to Landauer transport in an ordinary point contact or quantum wire. In the case of ordinary Landauer transport, the 1D channel is formed by pinching off a 2DEG using an insulating energy gap. For the Andreev junction, the 1D channel is formed by pinching off with a superconducting energy gap. For a wide junction, relative to the superconducting coherence length, there will be many channels of ABSs within the junction. For a narrow junction, provided the phase difference across the junction is $\pi$, there will remain at least a single pair of ABSs that disperse as a function of the momentum along the junction.

In Ref. [14], it was found that the source to drain Landauer conductance of an adiabatic and narrow Andreev junction is given by

[TABLE]

Here $c_{e}$($c_{h}$) are the number of propagating Andreev modes in the 1D channel that become electrons (holes) upon adiabatic evolution into the leads. For a narrow linear junction, these, in turn, are determined by the number of electronlike (holelike) critical points on the Fermi surface, where the Fermi velocity is parallel to the channel, and the Fermi surface is locally convex (concave). A novel feature of this result is that the conductance is different at positive and negative bias voltage. This leads to a rectification effect, in which a low-frequency AC voltage in the source, $V_{\text{AC}}=V_{0}\cos(\omega t)$, leads to a DC current flowing into the drain, $I_{\text{DC}}=-(e^{2}|V_{0}|/\pi h)(c_{e}-c_{h})$. While the integers $c_{e}$ and $c_{h}$ depend on the specific geometry of the Fermi surface, as well as its orientation relative to the 1D channel, their difference depends only on its topology:

[TABLE]

where $\chi_{F}$ is the Euler characteristic of the 2D Fermi sea. In general, $\chi_{F}$ for any 2D Fermi sea can be expressed as the sum over all disconnected components of the Fermi surface, where an electronlike (holelike) Fermi surface contributes $+1$ ($-1$) and an open Fermi surface contributes [math]. We will henceforth refer to this effect as topological Andreev rectification (TAR). It motivates us to define the rectified conductance

[TABLE]

It follows from Eq. (1) and Eq. (2) that $\delta G_{21}(V)$ probes the intrinsic topology of the 2D Fermi sea.

The Andreev junction resembles similar Josephson junction devices that have been studied in an effort to engineer proximity-induced topological superconductivity [21, 22, 23, 15, 16, 17]. Such systems based on InAs or HgTe quantum wells [23, 15, 16, 17] are among the promising venues for studying the TAR effect. Nonlocal transport measurements have proven to be a useful diagnostic for probing topological superconductivity in Josephson junction setups [18, 19, 20], as well as in proximitized quantum wires [24, 25]. In that context, a rectification effect has been studied [26], which results from the specific structure of the Andreev levels in a proximitized 1D nanowire. Rectification in Andreev interferometers that require broken time-reversal symmetry has also been studied [27]. However, these are all distinct from the TAR effect, which probes the intrinsic topology of the Fermi sea of a 2DEG.

In this paper, we will expand on the analysis in Ref. [14], and in an effort to hasten experimental demonstration of the TAR effect, we will consider more realistic device geometries that go beyond the narrow junction and fully adiabatic limits considered in Ref. [14]. We will begin with a review of the narrow junction limit in Section II.1. Wider junctions will be considered in Section II.2. In the latter case there exist additional Andreev modes bound to the SNS junction. These modes will be critical to our analysis because they can contribute to the nonlocal conductance, and when there is scattering among them, the quantization is degraded. We will identify a physically accessible regime in which interchannel scattering among modes with nearly the same momentum, that are associated with the same Fermi surface critical point, degrades the quantization of the nonlocal conductance $G_{21}$. However, we will show by analyzing the multichannel Landauer transmission problem that, provided large momentum scattering between Fermi surface critical points is absent, the rectified conductance $\delta G_{21}$ remains quantized. Thus, the rectified conductance provides a more robust marker for the intrinsic topology of the Fermi sea.

In order to verify these predictions, we have performed extensive simulations of realistic device geometries using the Kwant package [28, 29]. In Section III we introduce a simple toy model that allows us to consider a wide variety of Fermi sea topologies and real-space junction geometries that access several regimes of interest, including narrow junctions, wide junctions, as well as a point-contact geometry. We find that in all of these cases, the TAR effect can be observed provided the transport is sufficiently ballistic. In Section IV, we model junctions in specific real materials that are motivated by recent experimental studies of Josephson junctions in InAs quantum wells, as well as in monolayer and bilayer graphene. In these simulations we assess the feasibility of measuring the TAR effect using experimental parameters.

II Topological Andreev rectification

In this section we develop the analytic framework for topological Andreev rectification. We will begin in Section II.1 by reviewing the analysis in Ref. [14] of the fully adiabatic narrow-junction limit. In Section II.2 we introduce a multichannel scattering analysis that allows for the analysis of wider junctions that are not perfectly adiabatic. We will show that provided large momentum scattering is suppressed, the rectified conductance $\delta G_{21}$ remains a quantized probe of the Fermi sea topology, even when the quantization of the nonlocal conductance $G_{21}$ is degraded.

II.1 Review

It was established in Ref. [14] that the rectified nonlocal conductance of a narrow and adiabatic Andreev junction is quantized and probes the Euler characteristic of the 2D Fermi sea. Here we review that analysis.

II.1.1 Euler characteristic

The Euler characteristic is an integer topological invariant that characterizes any topological space. Here we review several equivalent formulations of the Euler characteristic, $\chi_{F}$, of a two-dimensional Fermi sea. For concreteness, consider the hypothetical Fermi sea with $\chi_{F}=-1$ shown in Fig. 2.

a. In terms of Betti numbers: The Euler characteristic may be written as [30, 31]

[TABLE]

where the Betti numbers $b_{\ell}$ count the number of topologically distinct $\ell$-cycles. In Fig. 2, the Fermi sea has a single connected component, hence $b_{0}=1$. The two holes in the Fermi sea imply two independent $1$-cycles, hence $b_{1}=2$, and $b_{\ell}=0$ for $\ell>1$.

b. In terms of critical points of $E({\bm{k}})$: According to Morse theory, $\chi_{F}$ can be expressed in terms of the critical points of a Morse function [32, 33]. The electronic dispersion $E({\bm{k}})$ is a natural Morse function, and the critical points where $\nabla_{\bm{k}}E({\bm{k}})=0$ are characterized by a Morse index $\eta=(0,1,2)$ for a (minimum, saddle, maximum). Then,

[TABLE]

where $n_{\eta}$ is the number of critical points inside the Fermi sea with index $\eta$. Figure 2 identifies one minimum $(\eta=0)$ and two saddles ($\eta=1$). This formulation makes it clear that a topological Lifshitz transition occurs when a Fermi sea critical point passes through $E_{F}$ [34, 35].

c. In terms of genus and boundary components: In two dimensions $\chi_{F}$ can be expressed as a sum over all connected components $k$ of the Fermi sea,

[TABLE]

where $g_{k}$ is the genus of component $k$ and $b_{k}$ is the number of boundary components. In Fig. 2, there is a single component with $g=0$ and $b=3$, so $\chi_{F}=-1$. A non-trivial genus can arise because the Brillouin zone is a torus. For example, the complement of the Fermi sea in Fig. 2 has three components. Each of the two holes has $g=0$ and $b=1$, while the outer component wraps the Brillouin zone and thus has $g=1$ and $b=1$, leading to $\chi_{F}=+1$ altogether. In general, $\chi_{F}$ is odd under the exchange of electrons and holes for an even-dimensional Fermi sea.

d. In terms of the Fermi surface: In two dimensions every component of the Fermi surface has the topology of a circle, but they can be distinguished by whether the circle encloses electrons, encloses holes, or is open (i.e., encloses neither electrons nor holes). Then,

[TABLE]

where $n_{e(h)}$ counts the number of electronlike (holelike) Fermi surfaces, whereas open Fermi surfaces contribute zero. In Fig. 2, we have $n_{e}=1$ and $n_{h}=2$.

e. In terms of Fermi surface critical points: If we specify an arbitrary unit vector $\hat{\zeta}$, then we may identify critical points on the Fermi surface to be the points where the velocity $\bm{v}(\bm{k})=\hbar^{-1}\nabla_{\bm{k}}E({\bm{k}})\parallel\hat{\zeta}$. Fermi surface critical points can be distinguished by the sign of the curvature of the Fermi surface at that point. We call the critical point convex, or electronlike, if $\partial^{2}E/\partial k_{\perp}^{2}>0$ (where $k_{\perp}$ is the component of ${\bm{k}}$ perpendicular to $\hat{\zeta}\,$). In the vicinity of such a critical point the Fermi surface resembles an electron pocket. Alternatively, the critical point is concave, or holelike, if $\partial^{2}E/\partial k_{\perp}^{2}<0$, and the Fermi surface locally resembles a hole pocket. In terms of these data,

[TABLE]

where $c_{e(h)}$ gives the number of electronlike (holelike) critical points. In Fig. 2, the electronlike and holelike critical points for $\hat{\zeta}=+\hat{x}$ are, respectively, shown as the red and blue dots, with $c_{e}=2$ and $c_{h}=3$.

For the purposes of this paper this last formulation provides the most natural connection to the Andreev junction. We will show that the low-energy transport in the junction is dominated by the Fermi surface critical points where the velocity is parallel to the junction (i.e., with $\hat{\zeta}$ parallel to the NS interfaces). It is worth emphasizing that the integers $c_{e}$ and $c_{h}$ are not by themselves topological invariants. They depend on the specific geometry of the Fermi surface as well as on the chosen direction for defining the critical points. However, the difference $c_{e}-c_{h}$ depends only on the intrinsic topology of the Fermi sea. For example, if instead in Fig. 2 we had chosen $\hat{\zeta}=+\hat{y}$, then we would have $c_{e}=1$ and $c_{h}=2$.

II.1.2 Andreev zero modes at a $\pi$ junction

Fermi surface critical points, as defined above, are probed by Andreev bound states (ABSs). In this work, as well as in Ref. [14], we focus on the ABSs hosted by an SNS $\pi$ junction (i.e., with a superconducting phase difference of $\pi$, see Fig. 1). The reason for focusing on $\pi$ junctions is that in this case the presence of subgap ABSs is a robust topological feature. To see this, consider an infinitely-long channel oriented in the $\hat{x}$ direction, such that the electronic states can be indexed by the momentum $k_{x}$. The Andreev states bound to the junction are found by solving the Bogoliubov-de Gennes (BdG) equation for motion in the $\hat{y}$ direction for each value of $k_{x}$. For values of $k_{x}$ that are away from the Fermi surface critical points (defined for $\hat{\zeta}=\pm\hat{x}$) the electronic dispersion along the $k_{y}$ direction can be approximated in the vicinity of a Fermi surface point $k_{F}^{y}(k_{x})$ by linearizing: $E(k_{x},k_{y})\approx E_{F}+v_{y}(k_{x})\left[k_{y}-k_{F}^{y}(k_{x})\right]$ (here and throughout this section we set $\hbar=1$). When there is a single pair of Fermi surface points at $\bm{k}=(k_{x},\pm k_{F}^{y})$ with opposite velocity $\pm v_{y}$, the resulting BdG equation takes the form of a 1D Dirac equation with a spatially-varying mass term due to the pairing potential $\Delta=\Delta_{1}+i\Delta_{2}$. The corresponding BdG Hamiltonian is

[TABLE]

where $\sigma_{z}=\pm 1$ distinguishes the Fermi points, and $\bm{\tau}$ are Pauli matrices in the Nambu space. When the phase difference across the junction is $\pi$, we may set $\Delta_{1}=0$, and $\Delta_{2}(y)$ changes sign across the junction. Equation (9) then exhibits a pair of degenerate Jackiw-Rebbi zero modes indexed by $\sigma_{z}$ [36]. For values of $k_{x}$ where there are no Fermi surface points, then clearly there will be no zero modes. If $k_{x}$ is such that there are $2n$ Fermi surface points, then there will be $2n$ zero modes. Thus, the count of zero modes of Eq. (9) as a function of $k_{x}$ contains geometric information about the Fermi sea.

The exactness of the zero-mode solutions of Eq. (9) relies on the absence of any normal reflection at the interface, which guarantees that $[\mathcal{H},\sigma_{z}]=0$. In this work, we focus on junctions that are perfectly transparent, so that the Hamiltonian contains no scattering terms that connect the Fermi points. However, even for a perfect junction, the exactness of the zero modes is an artifact of the linear approximation in Eq. (9), which treats the $\sigma_{z}=\pm 1$ bands as completely independent. In reality, states at different Fermi points can derive from the same band, and a more complete description can be formulated in terms of a single-band model described by a second (or higher) order differential equation. As shown below, this leads to a splitting of the zero modes by an amount of order $|\Delta|^{2}/\mu(k_{x})$, where $\mu(k_{x})$ is the Fermi energy relative to the band edge at a given $k_{x}$.

When $k_{x}$ approaches a Fermi surface critical point, the linearized approximation necessarily breaks down because $\mu(k_{x})\rightarrow 0$. This is consistent with the fact that when $k_{x}$ crosses a critical point the zero mode count will change. In the vicinity of the critical points the energy of the nominal zero modes will disperse as a function of $k_{x}$. Our working assumption in this paper is that $|\Delta|\ll E_{F}$, where $E_{F}$ is the appropriate energy scale associated with the Fermi sea. Deep inside the Fermi sea (away from the surface critical points), the splitting of the zero modes is thus of order $|\Delta|^{2}/E_{F}$, which is much smaller than $|\Delta|$. In this case, the zero modes will exhibit significant dispersion as a function of $k_{x}$ only near the critical points. Thus, the dispersing Andreev modes are in correspondence with the critical points, and can be analyzed by focusing on the vicinity of each critical point. This analysis was introduced in Ref. [14], and is reviewed below.

II.1.3 Dispersive Andreev modes associated with Fermi surface critical points

Close to a Fermi surface critical point at which $v_{y}(k_{x})=\partial E/\partial k_{y}=0$, the energy can be linearized as a function of $k_{x}$, but it is necessary to expand to second order in $k_{y}$. This leads to the following model for an Andreev junction of width $W$ [14],

[TABLE]

Here $\delta k_{x}\equiv k_{x}-k_{F}$ is the $k_{x}$ momentum relative to the critical point at $\bm{k}=(k_{F},0)$, $\tau_{y}$ and $\tau_{z}$ are Pauli matrices acting in the electron-hole space, $\Delta_{0}$ is the magnitude of the proximity-induced pairing gap, and $m^{*}\gtrless 0$ is the effective mass at the convex/concave critical point. There are two symmetries useful for finding ABSs: (1) Chiral symmetry $\tau_{x}$, which relates states of energy $\pm\varepsilon$ for a fixed $\delta k_{x}$. (2) Mirror symmetry $\tau_{z}\mathcal{M}_{y}$ (where $\mathcal{M}_{y}$ takes $y\mapsto-y$), which labels each ABS by a mirror eigenvalue $\pm 1$. Clearly, electron-hole partners (related by $\tau_{x}$) acquire opposite mirror eigenvalues.

For $W=0$ exactly one pair of ABSs is found with dispersion [14]

[TABLE]

Deep inside the Fermi pocket where $v_{x}\delta k_{x}/m^{*}<0$, the dispersion of the ABS is flattened with $\varepsilon\approx\pm\Delta_{0}^{2}/(2v_{x}\delta k_{x})$. For $\delta k_{x}\sim k_{F}$, i.e., the momentum scale associated with the size of the Fermi sea, we have $v_{x}\delta k_{x}\sim E_{F}$ which recovers the $\Delta_{0}^{2}/E_{F}$ splitting mentioned above, and is consistent with the expectation of having approximate zero modes for $\Delta_{0}\ll E_{F}$. On the other hand, far outside the Fermi pocket, where $v_{x}\delta k_{x}/m^{*}>0$, the ABS energy approaches the dispersion of the metal with $\varepsilon\approx\pm v_{x}\delta k_{x}$, so the zero mode is absent. Right at the critical point $\varepsilon=\pm\Delta_{0}/\sqrt{2}$. While the exact form of the dispersion depends on microscopic details (e.g., $v_{x}$ and $W$), the presence of such a dispersive ABS for each Fermi surface critical point is a robust topological feature of the narrow $\pi$ junction. Figure 3(a) shows the spectrum of the ABS near a convex critical point for a finite-width $\pi$ junction with $W\ll|2m^{*}\Delta_{0}|^{-1/2}$, which is well approximated by (11).

Just as there are two types of Fermi surface critical points (convex or concave), there are also two types of dispersive ABSs, which are distinguished by whether they are electronlike or holelike at positive energy. We define the electron-hole character of an ABS excitation as $\langle\tau_{z}\rangle$ for the BdG mode with $\varepsilon>0$. This is equal to $+1$ ($-1$) for a normal electron (hole). Note that the electron-hole character is related to the slope of the ABS dispersion by the Hellman-Feynman theorem: $d\varepsilon/dk_{x}=v_{x}\langle\tau_{z}\rangle$. Thus, the electron-hole character, indicated by the color of the curves in Fig. 3, is related to the slope of the curves. As we have established that the ABS dispersion $\varepsilon\rightarrow\pm v_{x}\delta k_{x}$ far outside the Fermi pocket, we find $\langle\tau_{z}\rangle\rightarrow\pm 1$. A convex (electronlike) critical point [as shown in Fig. 3(a)] is associated to an electronlike ABS excitation with a positive electron-hole character. A concave (holelike) critical point will have a dispersion that resembles Fig. 3(a) rotated by $180^{\circ}$, and will be associated with a holelike ABS excitation with a negative electron-hole character.

The electron-hole character of an ABS excitation, in turn, determines the fate of the mode as it propagates into a normal lead. In Ref. [14], a lead was modeled in two ways. The simplest is to imagine adiabatically decreasing $\Delta_{0}$ (while keeping $W$ fixed). In that case, the energy of the single ABS excitation relative to $\Delta_{0}$ grows, so that in Fig. 3(a) $\varepsilon/\Delta_{0}$ moves away from zero. An initially low-energy Andreev state, with a mixed electron-hole character, will thus evolve into a state with a definite electron-hole character $\langle\tau_{z}\rangle=\pm 1$. A somewhat more realistic model for an adiabatic contact is to consider increasing $W$ sufficiently slow (keeping $\Delta_{0}$ fixed). As $W$ increases, there will be additional Andreev modes bound to the junction, as shown in Fig. 3(b). In the adiabatic approximation, a state will follow the same band as it smoothly evolves upon increasing $W$, and it can again be seen that deep inside the lead an ABS excitation will acquire a definite electron or hole character, which is determined by whether the critical point is convex or concave. This, in turn, determines the current carried by a quasiparticle populating the Andreev mode when it propagates into a lead. It also specifies the population of the Andreev mode, which depends on the chemical potential of the lead it is incident from.

We thus conclude that each electronlike (holelike) critical point is associated with a propagating ABS that provides a conducting channel connecting the source and the drain for electrons (holes). We next show that these modes lead to a quantized Landauer conductance that depends on the sign of the voltage bias.

II.1.4 Quantized transport via Andreev states

We now apply the relation between the convex or concave nature of a Fermi surface critical point and the electron or hole nature of the associated dispersive ABS to transport along a narrow Andreev junction as depicted in Fig. 1. We are interested in the nonlocal differential conductance defined as (see Appendix A),

[TABLE]

Here $T^{ee(he)}_{21}(\varepsilon)$ describes the transmission of an electron incident from lead N1 at energy $\varepsilon$ to lead N2, where it emerges as an electron (hole). For $\varepsilon<0$ this could equally well be interpreted as $T^{hh(eh)}_{21}(-\varepsilon)$, i.e., the transmission of a hole at energy $-\varepsilon$ in N1 to a hole (electron) in N2. The equivalence of these interpretations is guaranteed by the electron-hole symmetry of the BdG theory. When the junction length $L$ is much longer than the superconducting coherence length $\xi=v_{x}/\pi\Delta_{0}$, crossed Andreev reflection is suppressed and hence $T^{he}_{21}=0$. Moreover, the dispersive ABSs are the only transmission channels between leads N1 and N2. Since there are $c_{e}$ ($c_{h}$) number of electronlike (holelike) ABSs, and assuming reflectionless contacts and ballistic transport, we find

[TABLE]

which implies

[TABLE]

Notice that the quantization in $G_{21}$ is sensitive to the shape of the Fermi sea as well as the orientation of the junction. On the other hand, the rectified nonlocal conductance, as defined in Eq. (3), obeys

[TABLE]

which is a quantized probe of the intrinsic topology of the Fermi sea. We refer to the quantization in $\delta G_{21}$ as the topological Andreev rectification (TAR) effect.

Equation (14) is the central result of Ref. [14], which is derived for a narrow Andreev junction hosting one ABS per Fermi surface critical point, and under the assumption that the ABSs moving out of the junction would adiabatically evolve into definite electrons or holes inside the normal lead. In this case, Eq. (15) is simply a corollary. However, as we argue below, for a wider junction, which can have multiple ABSs, there is a regime in which the quantization of $G_{21}$ is degraded, while the quantization of $\delta G_{21}$ remains robust.

The reason for the existence of this regime can be understood by examining the Andreev mode dispersion near a critical point for a wider junction, as shown in Fig. 3(b). We first note that the ABSs of a wide Andreev junction tend to bundle in pairs with opposite mirror eigenvalues. This is related to the oscillations in the spectrum which originate from anticrossings between the same $\tau_{z}\mathcal{M}_{y}$ mirror sector and symmetry-protected crossings between different mirror sectors. In particular, at zero energy, different sectors have to cross because $\tau_{z}\mathcal{M}_{y}$ anticommutes with chiral symmetry $\tau_{x}$ [14]. The zero-energy crossings are determined by $k_{y}^{2}/2m^{*}+v_{x}\delta k_{x}=0$ where $k_{y}$ is obtained from the boundary conditions. In the limit $\Delta_{0}/v_{x}\delta k_{x}\gg 1$, the wave function vanishes at $y=\pm W/2$. This gives the condition $k_{y}=n\pi/W$ ($n=1,2,\ldots$) and hence $v_{x}\delta k_{x}/\Delta_{0}=-\operatorname{sgn}(m^{*})(n\pi/W\sqrt{|2m^{*}\Delta_{0}|}\,)^{2}$, which holds approximately in Fig. 3(b). We further note that the Andreev level spacing is smaller when $k_{x}$ is close to the critical point. In general, the Andreev level spacing is of order $v_{y}(k_{x})/W$. In Eq. (10), $v_{y}=k_{F}^{y}(k_{x})/m^{*}\sim\left|v_{x}\delta k_{x}/m^{*}\right|^{1/2}$. Deep inside the Fermi pocket, $\delta k_{x}\sim k_{F}$ and $v_{y}\sim v_{F}$, so the deep Andreev level spacing is of order $v_{F}/W$. However, $v_{y}\rightarrow 0$ when $\delta k_{x}$ is small. This leads to a level spacing of order $1/m^{*}W^{2}$ close to the critical point. We will focus on the regime in which the energy of the relevant dispersive Andreev modes, which is set by the bias voltage $V_{1}$, is of order (or less than) the deep level spacing $v_{F}/W$. In this case there can still be several Andreev bound states for $k_{x}$ close to the critical point. These modes have nearly the same momentum, which makes the adiabatic approximation less reliable for them. Moreover, even when the width of the junction increases adiabatically slow, the Fermi energy will pass through turning points where these modes come and go in a complicated way. This inevitably leads to a breakdown in the adiabatic approximation, resulting in a degradation of the quantization of the nonlocal conductance $G_{21}$. In the following, we will show that despite this breakdown in the quantization of $G_{21}$, the rectified conductance $\delta G_{21}$ remains robustly quantized provided large momentum scattering between the Fermi surface critical points is negligible.

II.2 Beyond the narrow-junction limit

In this section we consider a junction of width $W\gtrsim|2m^{*}\Delta_{0}|^{-1/2}$, which can have multiple ABSs associated with each Fermi surface critical point. We analyze the coupling of these modes to the leads within a multichannel Landauer-Büttiker formalism. Remarkably, we find that the quantization of the rectified nonlocal conductance $\delta G_{21}$ persists even for a wide junction as long as scattering between ABSs across the Fermi surface is suppressed, and either time-reversal symmetry $\mathcal{T}$ or mirror symmetry $\mathcal{M}_{x}$ ($x\mapsto-x$) along the transport direction is preserved. We further find that $\delta G_{21}$ is related to the number of occupied ABSs deep in the Fermi pockets. Hence, at sufficiently small bias voltage the plateau in $\delta G_{21}$ remains a good measure of the Fermi sea topology $\chi_{F}$.

We start by dividing the total system into the two leads and the junction, where the latter is defined as the region with constant width $W$. As illustrated in Fig. 4, scattering from the two leads to the junction is governed by two scattering matrices $S_{1}$ and $S_{2}$. We have

[TABLE]

with $a$ and $b$ the incoming and outgoing amplitudes of propagating modes in the leads, and the subscripts indicating the lead and the propagation direction. The amplitudes of right- or left-propagating ABSs are given by $\phi_{R/L}^{\text{ABS}}$. Here we have absorbed the phases accumulated during propagation into the $S$ matrices,

[TABLE]

where the subscript $R/L$ denotes transmission ($t$) or reflection ($r$) to the right/left.

Let us consider scattering between lead N1 and the junction in closer detail. The transmission matrix from lead N1 to the junction is the $N\times M_{1}$ matrix $t_{1R}$ where $N$ ($M_{1}$) is the number of right movers in the junction (lead N1). Since $M_{1}\gg N$, the matrix $t_{1R}^{\dagger}t_{1R}$ has at least $M_{1}-N$ zero eigenvalues by the rank-nullity theorem, which also holds for $t_{1L}t_{1L}^{\dagger}$. We now choose a basis for the incoming (outgoing) lead modes that diagonalizes $t_{1R}^{\dagger}t_{1R}$ (and $t_{1L}t_{1L}^{\dagger}$). In this basis, we can write

[TABLE]

where the lead modes with zero transmission eigenvalues are denoted as $\mathrm{\overline{lead}}$ and

[TABLE]

Here we have used the same notation even though the submatrices are generally different after the basis transformation. In the new notation $r_{1L}$, $r_{1R}$, $t_{1L}$, and $t_{1R}$ are $N\times N$ matrices222In the event that the original $t_{1R}^{\dagger}t_{1R}$ (or $t_{1L}t_{1L}^{\dagger}$) has more than $M_{1}-N$ zero modes, then there will exist zero modes of $t_{1R}t_{1R}^{\dagger}$ (or $t_{1L}^{\dagger}t_{1L}$), implying modes in the junction that are perfectly reflected and decoupled from the lead. In that case the new $t_{1L}$ and $t_{1R}$ matrices are not square, which complicates our analysis. To resolve that, we add a small perturbation that lifts the extra zero modes, making $t_{1L}$ and $t_{1R}$ square. That perturbation can then be set to zero in the decoupled scattering problem. , while $\rho_{1L}$ is square with dimension $M_{1}-N$. The latter contains both normal and Andreev reflection coefficients of the zero transmission lead subspace. Moreover, unitarity of the scattering matrix $S_{1}^{\dagger}S_{1}=S_{1}S_{1}^{\dagger}=\mathds{1}$ implies

[TABLE]

By construction, $t_{1L}^{\dagger}$ and $t_{1R}$ have trivial nullspaces. It follows that $A$ and $B$ are zero. We obtain

[TABLE]

Because we can perform the same procedure for the lead N2, the scattering problem is decoupled into three independent problems: two purely reflection problems of dimension $M_{i}-N$ ($i=1,2$ ) for the zero transmission lead modes and an $N$-dimensional problem involving both the junction and the lead modes with nonzero transmission eigenvalues. Since we are interested in the latter, we can henceforth focus on the smaller subsystem with $N\times N$ scattering matrices $s_{1}$ and $s_{2}$, as illustrated in Fig. 5.

Up to this point, the analysis has been completely general. We now restrict our attention to the case where all ABSs (at a given energy) have momentum $k_{x}\approx\pm k_{F}$. This is generally satisfied for energies $\varepsilon$ away from the Andreev levels deep inside the Fermi pocket, as shown in Fig. 3(b). These ABSs are associated to a Fermi surface critical point at $k_{x}=+k_{F}$, and given our definition for the critical point (with $\hat{\zeta}=+\hat{x}$), right movers at $k_{x}\approx+k_{F}$ are always electronlike while left movers are holelike. Then either $\mathcal{T}$ or $\mathcal{M}_{x}$ imply that the electronlike (holelike) ABSs at $k_{x}\approx-k_{F}$ are left (right) movers. Furthermore, in the adiabatic limit where the width of the junction varies slowly with respect to $k_{F}^{-1}$, scattering between states at $+k_{F}$ and $-k_{F}$ is negligible. Note that there are still normal and Andreev reflections in the zero transmission lead subspace, but these do not affect the transmission between leads N1 and N2. Therefore the total $S$ matrix of the reduced system takes on a block-diagonal form in the adiabatic limit,

[TABLE]

where $\odot$ stands for the composition of the reduced scattering matrices, see Fig. 5. The upper-left (lower-right) block is designated to the $+k_{F}$ ($-k_{F}$) modes, the subscript $L/R$ indicates the outgoing direction, and the superscript $\alpha\beta$ indicates the conversion of a $\beta$-like mode into an $\alpha$-like mode. Here the reflection and transmission matrices in Eq. (22) are related to those in $S_{1}$ and $S_{2}$. For example,

[TABLE]

Notice that in Eq. (22) we have assumed that the amplitudes $t_{R}^{eh}$ and $t_{R}^{he}$ in Eq. (23a), as well as $r_{L}^{ee}$ and $r_{L}^{hh}$ in Eq. (23b), all vanish due to the suppressed scattering across the Fermi surface.

The nonlocal differential conductance at zero temperature is now given by

[TABLE]

where we used the electron-hole symmetry of the BdG Hamiltonian in the last line. Unitarity of $s$ implies

[TABLE]

with $N=N_{e}+N_{h}$, and $N_{e}$ ($N_{h}$) the number of right (left) movers at $+k_{F}$ (see Fig. 5). Taking the trace and subtracting yields

[TABLE]

which are all evaluated at energy $\varepsilon=eV_{1}$. Finally, when either $\mathcal{T}$ or $\mathcal{M}_{x}$ is preserved, we may interchange $R$ and $L$ in Eq. (24b), so that the rectified conductance becomes

[TABLE]

The validity of this result is independent of the number of ABSs in the junction and depends only on the assumption that scattering across the Fermi surface is strongly suppressed for the nonzero transmission subspace, as well as the presence of $\mathcal{T}$ or $\mathcal{M}_{x}$ 333When both $\mathcal{T}$ and $\mathcal{M}_{x}$ are present, we have $\delta G_{21}(V_{1})=-\delta G_{11}(V_{1})\equiv G_{11}(-V_{1})-G_{11}(V_{1})$, so the Fermi sea topology can also be probed via the rectified local conductance in that case. See Appendix A for more discussion.. We thus see that the rectified conductance $\delta G_{21}$ is still quantized for a wide Andreev junction, provided that the applied bias does not excite any Andreev states deep inside the Fermi pocket. In particular, there always exists a window of small bias $eV_{1}>0$ such that $N_{e}-N_{h}=c_{e}-c_{h}=\chi_{F}$. Upon increasing the bias above the deep Andreev level spacing $v_{F}/W$, one obtains additional quantized plateaus in $\delta G_{21}$ at larger values.

In the narrow-junction limit ($W\ll|2m^{*}\Delta_{0}|^{-1/2}$) the number of positive energy electronlike (holelike) ABSs at $+k_{F}$ is given by the number of convex (concave) critical points of the Fermi surface $c_{e}$ ($c_{h}$) defined relative to the transport direction. As this is independent of the subgap energy we obtain $\delta G_{21}\rightarrow\operatorname{sgn}(eV_{1})(e^{2}/h)\chi_{F}$ for a narrow junction in accordance with Eq. (15).

III Toy model

To demonstrate and verify our results, we perform transport simulations with the Kwant package [28, 29]. We first consider a toy model on a square lattice (with lattice constant $a$) with a single $s$ orbital per site, as illustrated in Fig. 6(a). The corresponding setup for the transport calculation is shown in Fig. 6(b), where a small system is depicted for clarity. The scattering region is connected to two normal leads on the left (N1) and right (N2), and to two proximitized superconducting leads (S) on the top ($\Delta=\Delta_{0}>0$ ) and bottom ($\Delta=-\Delta_{0}$).

We consider a model with only nearest-neighbor hopping, with amplitudes $t_{x}$ and $t_{y}$ along the $x$ and $y$ directions, respectively. The mean-field Hamiltonian is then given by $\hat{H}=\hat{H}_{0}+\hat{H}_{\Delta}$, with

[TABLE]

where $\hat{c}_{\bm{n},\sigma}^{\dagger}$ ($\hat{c}_{\bm{n},\sigma}$) creates (destroys) an electron with spin $\sigma$ on site $\bm{n}=(n_{x},n_{y})\in\mathds{Z}^{2}$, $\Delta_{\bm{n}}$ is the superconducting pair potential at site $\bm{n}$, and $\mu$ is the chemical potential. The normal Hamiltonian $\hat{H}_{0}$ has dispersion

[TABLE]

where the Fermi sea is defined by $E(\bm{k})\leq 0$. In the following, we consider anisotropic hopping amplitudes with $t_{y}/t_{x}=0.5$, and $t_{x}>0$. In this way, we can access three types of Fermi sea topology:

[TABLE]

Representative cases are illustrated in Fig. 6(c). The count of Fermi surface critical points is $(c_{e},c_{h})=(1,0)$, $(1,1)$, $(0,1)$, respective to $\chi_{F}=+1$, [math], $-1$. Since the normal metal considered here possesses a spin-degenerate Fermi sea, both $\chi_{F}$ and $c_{e,h}$ in this section are counted per spin.

As explained in the previous section, information on the Fermi sea topology is contained in the nonlocal conductance, which can be computed as (see Appendix A)

[TABLE]

where $f_{1}=f_{0}(\varepsilon-eV_{1})$ is the Fermi distribution in lead N1. The total transmission functions $T_{ij}^{\alpha\beta}(\varepsilon)$ for a charge carrier of type $\beta$ to be transmitted from lead $\text{N}j$ to lead $\text{N}i$ as a charge carrier of type $\alpha$ ($\alpha,\beta=e,h$ and $i,j=1,2$) are calculated in the BdG formalism using Kwant. In the rest of the work, we focus on transport properties at temperature $T=0$. Quantized plateaus in $G_{21}(V_{1})$ and $\delta G_{21}(V_{1})$ for subgap bias $|eV_{1}|<\Delta_{0}$ shall remain at finite temperatures provided $k_{B}T\ll\Delta_{0}$.

We now move on to discuss the numerical results for three kinds of junction geometries: (A) narrow and long Andreev junctions; (B) wide and long Andreev junctions; and (C) Andreev point contacts, where the latter are defined in Fig. 7. Here a long junction has by definition a length $L\gg\xi$, where $\xi$ is the superconducting coherence length. We take $\Delta_{0}=0.1t_{x}$ and estimate

[TABLE]

A narrow junction is defined as hosting a single dispersive ABS per Fermi surface critical point, see Fig. 3(a), while for a wide junction there are multiple ABSs, see Fig. 3(b).

In the following, we take a scattering region of dimensions $L_{x}\times L_{y}$ with $L_{x}=L_{y}=100a$. For the cases (A) and (B), the linear junction has length $L=L_{x}$.

III.1 Narrow Andreev junction

For a long junction, both crossed Andreev reflection and electron tunneling are suppressed. The only transmission channel from lead N1 to lead N2 is through the ABSs which disperse along the junction. Here we consider a normal region of width $W=6a\lesssim\xi$. As we have demonstrated in Sec. II.1, each Fermi surface critical point is associated to one pair of ABSs for a narrow Andreev junction, see Fig. 3(a). This is the regime where we expect the theoretical results, given by Eq. (14) and Eq. (15), to be applicable.

The numerical results confirm our predictions. As shown in Figs. 8(a)–(c), the nonlocal conductance $G_{21}$ is quantized (in units of $2e^{2}/h$) to the integral value $c_{e}$ (or $c_{h}$) for subgap bias voltage $eV_{1}>0$ (or $eV_{1}<0$). The quantized conductance is contributed by the dispersive ABSs close to the Fermi surface critical points (with $k_{x}\approx\pm k_{F}$), while the non-universal zero-bias peak originates from the zero-energy ABSs deep inside the Fermi sea (with $k_{x}\approx 0$). The quantized plateaus at finite bias are universal in the sense that they depend only on the counting of Fermi surface critical points, but not on the detailed shape of the Fermi sea (such as its size). When $\chi_{F}=c_{e}-c_{h}\neq 0$, the nonlocal conductance is asymmetric in the bias voltage $V_{1}$. This asymmetry is quantified by the rectified conductance $\delta G_{21}$, as defined in Eq. (3), and is color-plotted as a function of $\mu/t_{x}$ and $eV_{1}/\Delta_{0}$ in Fig. 8(d). In almost the entire subgap regime, $\delta G_{21}$ is found to be quantized to $\chi_{F}$ (in units of $2e^{2}/h$). In this case, both $G_{21}$ and $\delta G_{21}$ serve as good markers for the Fermi sea topology. In contrast, the features near zero bias depend on microscopic details, such as the effect of normal reflections at the NS interfaces, which split the nominal zero modes away from zero energy, similar to Fig. 3(a). This energy splitting oscillates as a function of the chemical potential with an envelope that scales as $\Delta_{0}^{2}/\mu$ for $0<\mu/t_{x}<1$, as can be seen close to zero bias in Fig. 8(d). We have verified that the width of the zero-bias peak matches the ABS splitting for an infinitely-long junction with the same $W$ and $\mu$.

III.2 Wide Andreev junction

Next, we simulate a wide normal region with $W=40a\gg\xi$. In this case, there are always multiple ABSs at each Fermi surface critical point. The BdG spectrum near a convex critical point then resembles Fig. 3(b). As explained in Sec. II.2, backscattering among counterpropagating ABSs (which share opposite electron-hole character) is present even in the adiabatic limit as these modes are grouped around $k_{x}\approx+k_{F}$ (or $k_{x}\approx-k_{F}$). Consequently, we do not expect $G_{21}$ to be quantized. Even when these scattering events are negligible, the number of ABSs available for transmission would depend on microscopic details, such as the value of the chemical potential (for fixed $\chi_{F}$) and the bias voltage in the lead. Nevertheless, as long as scattering across the Fermi surface is negligible, we have argued in Sec. II.2 that $\delta G_{21}$ remains quantized.

These predictions are verified by the numerical results, which are summarized in Fig. 9. As shown in (a)–(c), $G_{21}$ is no longer quantized to $c_{e}$ or $c_{h}$. Moreover, there is a rather strong dependence on both $\mu$ and $eV_{1}$. Hence $G_{21}$ per se no longer indicates the Fermi sea topology. However, the rectification effect remains robust and quantized, depending only on the topology of the Fermi sea. This is illustrated in Figs. 9(d)–(f), where we show $\delta G_{21}$ corresponding to (a)–(c), respectively. We further show $\delta G_{21}$ as a function of $\mu$ and $eV_{1}$ in Fig. 9(g). For a given Fermi sea topology, irrespective of the precise value of $\mu$, we find that $\delta G_{21}$ always displays a plateau which attains a quantized value of $\chi_{F}$ (in units of $2e^{2}/h$). This confirms that $\delta G_{21}$ remains a good marker for the Fermi sea topology for transport along wide junctions as long as scattering between states at $k_{x}\approx+k_{F}$ and states at $k_{x}\approx-k_{F}$ is suppressed. In our simulation, the interface between the normal lead and the Andreev junction involves a sharp jump in the pair potential on the lattice scale. However, ordinary backscattering (i.e., leaving electron-hole character unchanged) from $+k$ to $-k$ is mediated by the pair potential as a second-order process and thus suppressed by a factor on the order of $(\Delta_{0}/t_{x})^{2}$.

The wide Andreev junction hosts additional transport features which are absent in a narrow junction. Particularly, for $\chi_{F}=\pm 1$, as shown in Fig. 9, there is an extended region in which $\delta G_{21}$ attains a quantized value of $\pm 3$ (in units of $2e^{2}/h$). This is consistent with our prediction from the $S$-matrix analysis, see Eq. (27), where $\delta G_{21}(V_{1})$ is quantized to an integer that represents the difference between the number of right-moving electronlike ABSs ($N_{e}$) and the number of left-moving holelike ABSs ($N_{h}$), at energy $\varepsilon=eV_{1}$. The plateau transitions happen when $N_{e}-N_{h}$ changes its integer value, which correspond to the energies where an extra ABS deep inside the Fermi pocket ($k_{x}\approx 0$) becomes occupied/depleted. This is illustrated in Fig. 10 for an Andreev junction which hosts five ABSs deep inside the Fermi sea. Correspondingly, the rectified conductance $\delta G_{21}$, shown in Fig. 10(b), exhibits plateaus at $\pm 1$, $\pm 3$, and $\pm 5$ for subgap bias $eV_{1}\gtrless 0$, respectively. Note that plateaus appear only at odd integer multiples of $2e^{2}/h$. This is because the condition for quantization of $\delta G_{21}$ (i.e., the absence of scattering across the Fermi surface) is only satisfied after two extra ABSs are occupied/depleted. If the phase difference across the junction is [math] instead of $\pi$, so that the topological zero modes are absent, we find similar plateaus, but at even integer multiples of $2e^{2}/h$.

For $\chi_{F}=0$, which corresponds to $2<\mu/t_{x}<4$, the Fermi sea is surrounded by open Fermi surfaces as depicted in the center panel of Fig. 6(c). Hence there are no ABSs deep inside the Fermi sea, which implies $N_{e}=N_{h}$ for all subgap biases such that no plateau transitions can be found in Fig. 9(e).

III.3 Andreev point contact

As a final example of our toy model, we consider a point-contact geometry where the normal regions on the two sides of the superconducting constriction are coupled at a point, as depicted in Fig. 7. In contrast to the long Andreev junctions considered before, neither electron tunneling nor crossed Andreev reflection is suppressed for the Andreev point contact. Hence, we do not expect $G_{21}$ to be universal nor quantized. The numerical results for $G_{21}$, as shown in Fig. 11 (a)–(c), indeed exhibit a strong dependence on $\mu$ (for fixed $\chi_{F}$) and $eV_{1}$. Nevertheless, an asymmetry with respect to the bias is still present for $\chi_{F}\neq 0$, and the direction of rectification clearly depends on $\operatorname{sgn}(\chi_{F})$.

The theoretical rationale behind the topological rectification effect assumes that the ABSs are the only transmission modes across the constriction, and that they adiabatically evolve into definite electrons/holes inside the normal leads. These assumptions no longer stand for a point contact. Thus, we do not a priori expect any plateau in $\delta G_{21}$ for the Andreev point contact. Remarkably, the numerical results shown in Figs. 11(d)–(f) for $\delta G_{21}$, which correspond to the same parameters used in Figs. 11(a)–(c), as well as the color plot in Fig. 11(g), display a robust quantization of $\delta G_{21}$ in terms of $\chi_{F}$ for a wide range of parameter values. Hence for a small but finite bias voltage $eV_{1}$, $\delta G_{21}$ remains a hallmark of the Fermi sea topology for an Andreev point contact.

IV Materials

In this section, we discuss numerical results geared towards realistic material platforms. Motivated by recent experimental studies on Josephson junctions [16, 17, 18, 19, 20, 37, 38, 39], we focus on a two-dimensional electron gas (2DEG) as realized in InAs quantum wells (with strong spin-orbit coupling) and graphene systems, and examine the feasibility of observing quantized transport in Andreev junctions based on these systems.

IV.1 InAs quantum well

In recent years, Josephson junctions in hybrid systems based on InAs have been studied extensively with the aim to realize topological superconductivity and Majorana zero modes [23, 15, 17, 16, 18]. Moreover, nonlocal conductance measurements have recently been performed on hybrid InAs/Al devices as a probe for topological phase transitions [19, 20]. In the following analysis, we demonstrate that essentially the same setup can be used to probe TAR pertinent to the Fermi sea topology of the Rashba-split 2DEG in an InAs quantum well.

We employ the following tight-binding model on a square lattice (with lattice constant $a=10\;\text{nm}$) to simulate an Andreev junction in an InAs/Al system: $\hat{H}=\hat{H}_{0}+\hat{H}_{\Delta}$ with

[TABLE]

where $\hat{c}_{\bm{n},\sigma}^{\dagger}$ ($\hat{c}_{\bm{n},\sigma}$) creates (destroys) an electron with spin $\sigma$ on site $\bm{n}=(n_{x},n_{y})\in\mathds{Z}^{2}$, $\mu_{\bm{n}}$ and $\Delta_{\bm{n}}$ are respectively the chemical potential and pair potential on site $\bm{n}$. The effective hopping amplitude is $t=\hbar^{2}/(2m^{*}a^{2})=15\;\text{meV}$, where an effective mass $m^{*}=0.025m_{e}$ ($m_{e}$ being the bare electron mass) has been adopted. The Rashba spin-orbit coupling strength is chosen to be $\lambda_{R}=0.75\;\text{meV}$. In the superconducting region proximitized by Al leads [cf. Fig. 6(b) and Fig. 7], the magnitude of the pairing gap is chosen to be $\Delta_{0}=0.15\;\text{meV}$. These are typical parameter values adopted from recent experimental studies [18, 19, 20].

The chemical potential in the normal region (including normal leads) is fixed at $\mu_{N}=10\;\text{meV}$, which gives rise to two electronlike Fermi seas as depicted in Fig. 12(a). In practice, the chemical potential can be tuned by gating within a range of $\sim 1-100\;\text{meV}$. For most calculations presented below, we set the chemical potential in the superconducting regions $\mu_{S}$ equal to $\mu_{N}$. Given $\mu\sim 10\;\text{meV}$, the superconducting coherence length $\xi\sim 500\;\text{nm}$, which matches the experimental values reported in Refs. [18, 19, 20]. The parameter values adopted in our Kwant simulation are summarized in Table 1.

We first simulate a narrow Andreev junction where the normal region has width $W=100\;\text{nm}$. This is a realistic value based on recent experiments [18, 19, 20]. The setup follows Fig. 6(b). The junction is attached to normal leads of width $L_{y}=1\;\mathrm{\mu m}$, and we consider junctions with lengths varying from $2\;\mathrm{\mu m}$ to $6\;\mathrm{\mu m}$. The corresponding numerical results are summarized in Fig. 13. While sufficiently long junctions ($L_{x}\gtrsim 6\;\mathrm{\mu m}$) are required for observing quantized plateaus in $G_{21}$ that reflect the geometry of the Fermi sea $(c_{e}=2,\;c_{h}=0)$, a shorter junction ($L_{x}\gtrsim 3\;\mathrm{\mu m}$) already exhibits quantization in the rectified conductance $\delta G_{21}$ which encodes the topology of the Fermi sea ($\chi_{F}=2$).

In a real device, the superconducting region is most likely at a different chemical potential than the normal region, as the Al contacts locally dope the 2DEG. Such an electrostatic gradient decreases the transparency of the junction by enhancing ordinary reflections. This results in an increased energy splitting of the nominal Andreev zero modes. Thus, a quantized response is anticipated only for biases $|eV_{1}|$ above this energy scale. To model this effect, we consider a mismatch in the chemical potential across the NS interface that is sharp on the lattice scale (with $\mu_{N}=10\;\mathrm{meV}$ and $\mu_{S}=10.3\;\mathrm{meV}$). Our simulations shown in Fig. 13 demonstrate that the rectification and the concomitant signatures of the Fermi sea topology in the nonlocal conductance are stable for a step mismatch $\mu_{S}-\mu_{N}\sim\Delta_{0}$. While in practice the mismatch may be larger than $\Delta_{0}$, a value of the order of $\Delta_{0}$ for a step mismatch serves as a proxy for a more realistic potential profile that varies slowly on the scale of the Fermi wavelength $\lambda_{F}$. In this case, the transparency of the NS interface is only slightly reduced and the TAR effect is expected to remain.

Next, we consider an Andreev point contact using the same setup as shown in Fig. 7. Our numerical results are shown in Fig. 14. While there is no quantization in the nonlocal conductance $G_{21}$ for the point contact, as expected, the rectification effect remains observable. Moreover, for systems that are large compared to $\xi$, the rectified nonlocal conductance $\delta G_{21}$ is quantized and reflects the Fermi sea topology. This conclusion is unchanged when there is a mismatch in the chemical potential across the NS interface that is sharp on the lattice scale and of the order of $\Delta_{0}$.

We note that all our numerical results were obtained for disorder-free systems. For the state-of-the-art devices reported in Refs. [18, 19, 20], an electron mean free path in the range of $\sim 200-600\;\text{nm}$ was realized (i.e., comparable to the superconducting coherence length). We anticipate that with improved fabrication techniques, the mean free path can be further increased to realize a long and ballistic Andreev junction, in which both quantized $\delta G_{21}$ and $G_{21}$ can be measured.

IV.2 Graphene

Another interesting platform for realizing our proposal is graphene. Josephson junctions in graphene have been extensively studied both theoretically and experimentally [40, 41, 42, 43, 44, 45, 46, 47, 48, 37, 49, 38, 50, 39]. As such, it seems reasonable that nonlocal conductance measurements—to demonstrate the TAR effect and probe the Fermi sea topology—should be achievable in graphene systems.

We first consider monolayer graphene (MLG) and simulate transport in Andreev junctions in the ballistic regime for reasonable device geometries using a scaled honeycomb lattice, and taking realistic values for the proximitized superconducting gap $\Delta_{0}$. Secondly, we consider Bernal bilayer graphene (BLG) with an interlayer bias. This system has a highly-tunable Fermi sea as a function of the chemical potential and interlayer bias, with $|\chi_{F}|$ ranging from 0 to 3 (per spin and valley). It is therefore an interesting candidate for investigating TAR.

In all our graphene simulations, we take a scattering region of dimensions $L_{x}\times L_{y}$ with $L_{x}=L_{y}$, and consider a linear junction of length $L=L_{x}$ and width $W$, similar to what is shown in Fig. 6(b) for the square lattice. Normal (superconducting) leads are attached to the left and right (top and bottom) of the scattering region.

IV.2.1 Monolayer graphene

To demonstrate TAR in monolayer graphene (MLG) we consider the simplest lattice model for the charge carriers close to charge neutrality [51]. To this end, we consider a lattice model $\hat{H}=\hat{H}_{0}+\hat{H}_{\Delta}$ with

[TABLE]

where $\bm{R}=n_{1}\bm{a}_{1}+n_{2}\bm{a}_{2}$ are lattice vectors with $n_{1/2}\in\mathds{Z}$ and $\bm{a}_{1/2}=a(\sqrt{3}/2,\pm 1/2)$ with $a\approx 0.25\;\text{nm}$, see Fig. 15(a). Here $\hat{a}_{\bm{R},\sigma}^{\dagger}$ and $\hat{b}_{\bm{R},\sigma}^{\dagger}$ ($\hat{a}_{\bm{R},\sigma}$ and $\hat{b}_{\bm{R},\sigma}$) create (destroy) an electron with spin $\sigma=\uparrow,\downarrow$ on sublattice A and B, respectively, at sites $\bm{R}$ and $\bm{R}+\bm{\tau}$ with $\bm{\tau}=a\hat{x}/\sqrt{3}$. The nearest-neighbor hopping amplitude is given by $t\approx 2.8$ eV and $\Delta_{\bm{r}}$ gives the proximity-induced pair potential at position $\bm{r}$ with $\Delta_{\bm{r}}=\Delta_{0}\left[\theta(y-W/2)-\theta(-y-W/2)\right]$.

In the low-density regime, the energy bands are given by a pair of spin-degenerate Dirac cones at the two distinct zone corners (valleys) $K$ and $K^{\prime}$, as shown in Fig. 15(b). In this case, the Fermi sea topology (per spin and valley) is given by $\chi_{F}=\operatorname{sgn}(\mu/t)$ where the count of Fermi surface critical points, per spin and valley, is $(c_{e},c_{h})=(1,0)$ for electron doping and $(c_{e},c_{h})=(0,1)$ for hole doping. This is illustrated in Fig. 15(c).

Using $\Delta_{0}=0.2\;\text{meV}$ as the proximity-induced gap for Al contacts, we find $\xi/a=\sqrt{3}t/2\pi\Delta_{0}\sim 4\times 10^{3}$. To make the computation feasible for an intermediately long ($L_{x}\gtrsim\xi$) and narrow ($W\lesssim\xi$) Andreev junction, and since we are only considering the regime with a linear dispersion, we consider a scaled honeycomb lattice [52]. This is equivalent to solving the Dirac equation on a hexagonal grid scaled by a factor $s$, where

[TABLE]

such that $\hbar v_{F}=\sqrt{3}\,ta/2$ is invariant. This is justified as long as $s\ll|t/\mu|$, or equivalently $\lambda_{F}\gg sa$, such that we remain in the linear regime. Since we want to simulate a large system (comparable to $\xi\sim 1\;\mathrm{\mu m}$) we need a sufficiently large scaling factor $s$ which puts an upper bound on $|\mu|$. To this end, we consider carrier densities $|n|=k_{F}^{2}/\pi<10^{11}$ cm*-2* which, assuming a linear dispersion, yields $|\mu|<30$ meV and thus $s\ll 100$. In practice, it is more favorable to consider a higher carrier density with $|\mu|\gtrsim 100\;\text{meV}$, as a spatial variation of the chemical potential $\delta\mu\sim 5\;\text{meV}$ is present even for high-quality graphene devices on hBN substrates [53]. When $\mu$ is too small, the presence of electron-hole puddles across the device, which exhibit different Fermi sea topology, can preclude low-energy Andreev states even in a $\pi$ junction [37].

In the following, we take $s=12$, which gives a scaled lattice constant $sa$ of about $3$ nm and $\xi/sa\sim 300$. We then numerically implement an Andreev junction of length $L_{x}=1000\>sa\approx 3\;\mathrm{\mu m}$ and width $W=100\>sa\approx 300$ nm in Kwant [28, 29]. These dimensions are comparable to those of state-of-the-art devices [45, 37, 49, 38, 39]. An overview of the parameters used for our Kwant simulation is given in Table 2.

The nonlocal conductance $G_{21}$ and the rectified conductance $\delta G_{21}$ are shown in Fig. 16 for two values of the chemical potential corresponding to $\chi_{F}=\pm 1$ (per spin and valley), and with $\Delta_{0}=0.2\;\text{meV}$ (Al contacts) and $\Delta_{0}=1\;\text{meV}$ (Nb or MoRe contacts). We first notice that the qualitative features of the nonlocal response are insensitive to whether transport is along the armchair or zigzag direction. This is anticipated, as the Fermi sea has full rotational symmetry in the Dirac regime. For Al contacts [Fig. 16(a) and Fig. 16(c)] we have $L_{x}/\xi=3$ and the nonlocal conductance $G_{21}$ is not quantized to the count of Fermi surface critical points. This can be attributed to the presence of electron-tunneling across the rather short Andreev junction. Nevertheless, $\delta G_{21}$ is found to be nearly quantized to $\pm 1$ (in units of $4e^{2}/h$) allowing TAR to be observed in this setting. We also considered Nb or MoRe contacts [Fig. 16(b) and Fig. 16(d)] with $L_{x}/\xi=15$, which is in the long-junction regime. In this case, we find almost perfect quantization in $G_{21}$ which captures the Fermi sea geometry in accordance to Eq. (1). For realistic devices with short mean free paths compared to $L_{x}$ and $\xi$, the quantization in $G_{21}$ and $\delta G_{21}$ would both be degraded, but when the chemical potential is tuned through charge neutrality one still expects to observe a sudden jump in $\delta G_{21}$ which reflects the change in the Fermi sea topology. We thus believe that it is possible to observe signatures of TAR in current state-of-the-art hybrid superconductor-graphene devices.

IV.2.2 Bernal bilayer graphene

Bernal-stacked bilayer graphene (BLG) is a potentially more interesting platform for probing Fermi sea topology using TAR. Upon the application of an interlayer bias, by varying top and bottom gate voltages, this system hosts a much richer landscape of different Fermi sea topologies as compared to MLG. One can thus potentially tune experimentally between different values of $\chi_{F}$ by varying the interlayer bias as well as the chemical potential. Moreover, SNS junctions have been realized experimentally in BLG proximitized with Al [47] ($\Delta_{0}\approx 0.2\;\text{meV}$) as well as NbSe2 thin films [46] ($\Delta_{0}\approx 1\;\text{meV}$).

To investigate transport along an Andreev junction in BLG, we use the following lattice model, $\hat{H}=\sum_{l=1,2}\hat{H}_{l}+\hat{H}_{\perp}+\hat{H}_{\Delta}$ with

[TABLE]

where $\mu_{1/2}=\mu\pm U/2$ with $U$ the interlayer bias and $\mu$ the chemical potential. Here $\hat{a}_{\bm{R},l,\sigma}^{\dagger}$ and $\hat{b}_{\bm{R},l,\sigma}^{\dagger}$ ($\hat{a}_{\bm{R},l,\sigma}$ and $\hat{b}_{\bm{R},l,\sigma}$) create (destroy) an electron with spin $\sigma=\uparrow,\downarrow$ on layer $l=1,2$ on sublattice A and B, respectively. The in-plane position of sublattices A and B are, respectively, $\bm{R}$ and $\bm{R}+\bm{\tau}$ for layer $1$ and $\bm{R}-\bm{\tau}$ and $\bm{R}$ for layer 2. The most important hoppings [54] are given by the intralayer nearest-neighbor hopping amplitude $t$ and the interlayer hopping between eclipsing sites on different layers $t_{\perp}\approx 0.1t$. We also include second-nearest-neighbor interlayer hopping between different sublattices $\gamma_{3}\lesssim 0.1t$. The latter gives rise to trigonal warping [54] and allows for a richer variety of Fermi sea topology. The lattice and hoppings are illustrated in Fig. 17(a). For simplicity, we take the same value for the proximity-induced pair potential in both layers: $\Delta_{\bm{r}}=\Delta_{0}\left[\theta(y-W/2)-\theta(-y-W/2)\right]$. The parameters used for the Kwant simulation are given in Table 2.

The interlayer bias $U$ breaks inversion symmetry and opens a band gap $E_{g}$ at the $K$ and $K^{\prime}$ points. In the absence of $\gamma_{3}$ hopping, the Fermi sea at each valley is given by an annulus for $E_{g}/2<|\mu|<U/2$ [54]. Upon turning on $\gamma_{3}$, a small energy window appears starting from the band edge where the annular Fermi sea fractures into three electronlike ($\mu>0$) or holelike ($\mu<0$) Fermi pockets. We henceforth focus on $\mu>0$. Owing to a chiral symmetry that relates the conduction and valence band, the results for $\mu<0$ are obtained by reversing the sign of $\delta G_{21}$. Figure 17(b) shows the Fermi sea topology as a function of $\mu$ and $U$ where representative examples are illustrated in Fig. 17(c). For $\mu>U/2$, the Fermi sea (per valley and spin) has the topology of either a disk, where $\chi_{F}=+1$, or a pair of concentric disks when the second conduction band is also occupied, where $\chi_{F}=+2$.

To make the transport calculation feasible, instead of scaling the lattice (which gets tricky due to $\gamma_{3}$), we consider a larger superconducting gap $\Delta_{0}=0.0025t\approx 7\;\text{meV}$. The superconducting coherence length is then $\xi/a\sim\hbar v_{F}/\pi\Delta_{0}\sim 100$. While this value for $\Delta_{0}$ is about seven times larger than the proximity-induced gap for Nb or MoRe contacts, our aim is to present a proof of principle for TAR in BLG, which exhibits a rich landscape of Fermi sea topology. Note that our results would be unchanged if we could perform simulations for a larger system and a realistic $\Delta_{0}$. In fact, assuming a clean system, the quantization of $\delta G_{21}$ would only improve since $\Delta_{0}$ limits the energy resolution of TAR. We consider here a long and narrow Andreev junction for transport along the armchair direction. Here we take $L_{x}=600a$ and $W=4a$. Our numerical results for the rectified conductance $\delta G_{21}$ are summarized in Fig. 18. We see that Fig. 18(a) matches well to Fig. 17(b), confirming that $\delta G_{21}$ provides an excellent probe of the Fermi sea topology. The quantized plateaus shown in Fig. 18(b) further confirm that TAR is in principle observable over a wide range of bias voltages.

Finally, we comment on a potential issue related to the experimental implementation of the Andreev junction in BLG. The application of an interlayer bias across the entire sample (including the regions coupled to the superconductor) by means of an additional top gate can be challenging in practice due to screening by the superconductor. However, it is expected that screening is less important for few-layer van der Waals superconductors such as NbSe2 [46]. Given the appropriate device fabrication, we believe that the rich landscape of Fermi sea topology in BLG can indeed be extracted experimentally in a ballistic Andreev junction by measuring $\delta G_{21}$.

V Summary and Outlook

In this work, we presented both theoretical and numerical studies to consolidate a recent proposal in Ref. [14] which relates the topology of a two-dimensional Fermi sea to the ballistic transport of Andreev bound states (ABSs) in a superconducting $\pi$ junction. To highlight the significance of transport along the junction (i.e., between the two normal leads connected to the two ends of the junction) we termed the proposed setup an Andreev junction. Reference  [14] considered the limits of adiabatic and ballistic transport in a narrow Andreev junction of width $W\ll\xi$. In this case, the junction hosts the minimal number of ABSs and the nonlocal conductance $G_{21}(V_{1})$ is quantized to $c_{e}$ or $c_{h}$ (in units of $e^{2}/h$) for subgap bias $eV_{1}\gtrless 0$. Here $c_{e}$ ($c_{h}$) counts the number of convex (concave) critical points on the Fermi surface. While $c_{e/h}$ depends on the geometry of the Fermi sea, as well as the orientation of the Andreev junction relative to the Fermi sea, their difference gives the Euler characteristic $\chi_{F}=c_{e}-c_{h}$, which is an intrinsic topological characterization of the Fermi sea. This motivated us to introduce the rectified nonlocal conductance $\delta G_{21}(V_{1})\equiv G_{21}(V_{1})-G_{21}(-V_{1})$, which is expected to be quantized to $\chi_{F}$. We referred to this effect as topological Andreev rectification (TAR).

To evaluate the feasibility of observing TAR, we considered more general device geometries. We first provided a scattering-matrix analysis which also incorporates wide Andreev junctions (i.e., $W\gtrsim\xi$) hosting multiple dispersive ABSs around each Fermi surface critical point. The assumption of adiabaticity is then relaxed to allow for backscatterings among counterpropagating ABSs associated to the same Fermi surface critical point, while scattering between distinct Fermi surface critical points is still assumed to be suppressed. We established that, while $G_{21}$ is no longer quantized for wide junctions, $\delta G_{21}$ remains quantized to $\chi_{F}$ for small but finite bias, whereas additional quantized plateaus in $\delta G_{21}$ are observed at larger bias. Our theoretical analysis is further supported by numerical calculations using Kwant [29]. We first simulated nonlocal transport in Andreev junctions with a toy model on a square lattice. Our simulations for a narrow linear junction, as well as for a wide linear junction, confirm the predicted TAR effect. Furthermore, TAR is also demonstrated for an Andreev junction in a point-contact geometry. While such Andreev junctions cannot be studied analytically, our numerical results clearly establish TAR as a robust phenomenon insensitive to both real-space and momentum-space microscopic details.

Motivated by recent experimental progress on transport and spectroscopic measurements in Josephson junctions, we also simulated the TAR effect in realistic materials. Specifically, we considered a 2DEG in an InAs quantum well and graphene, for which superconducting junctions have already been fabricated (for InAs, see Refs. [17, 18, 19, 20]; for graphene, see Refs. [42, 43, 44, 45, 46, 47, 48, 37, 49, 38, 50, 39]). Using experimentally relevant parameters, and assuming a clean system, we examined the criteria for TAR to be observed in these systems. For an InAs quantum well proximitized by Al contacts ($\Delta_{0}\approx 0.2\;\text{meV}$) TAR is seen for a junction length $L\gtrsim 3~{}\mu$m when $\mu\sim 10\;\text{meV}$ (counted from the band bottom). We find that the rectified conductance $\delta G_{21}$ is nearly quantized to $\chi_{F}=2$ (in units of $e^{2}/h$) over a wide range of subgap bias voltages. For monolayer graphene, we considered both Al and Nb (or MoRe) contacts, where the latter induces a larger gap $\Delta_{0}\approx 1\;\text{meV}$. In this case, we considered a device of dimensions $3\times 3\;\mu$m2 and observed TAR where $\delta G_{21}$ is now quantized to $\pm 4e^{2}/h$ for electron or hole doping, respectively. Quantization of $G_{21}$ to $c_{e}e^{2}/h$ (for $eV_{1}>0$) and $c_{h}e^{2}/h$ (for $eV_{1}<0$), which reflects the Fermi sea geometry, is only obtained with Nb (or MoRe) contacts, owing to a smaller coherence length that puts our setup in the long-junction regime. Finally, for Bernal bilayer graphene, which hosts a rich variety of Fermi sea topologies upon tuning the interlayer bias and chemical potential, we demonstrate once again that $\delta G_{21}$ serves as a good marker of the Fermi sea topology. Our numerical results thus indicate that both InAs quantum wells and graphene are promising candidate platforms for observing topological Andreev rectification.

While our simulations have mostly used experimentally relevant parameter values, we have largely ignored the role of disorder. Our transport simulations are thus in the ballistic regime and perfect quantization in both $G_{21}$ and $\delta G_{21}$ can be achieved for $L\gg\xi$, in accordance with theoretical predictions. However, the current experimental status suggests that the electron mean free path $\ell_{e}$ is of the order of $\xi$. We therefore believe that a major challenge for observing quantization in both $G_{21}$ and $\delta G_{21}$ lies in increasing the mobility in these platforms such that $\ell_{e}$ becomes at least one order longer than $\xi$. Importantly, however, achieving quantization in $\delta G_{21}$ (i.e., TAR) is much easier than quantization in $G_{21}$. As indicated by our simulations, TAR is observed for shorter junctions for which $G_{21}$ is not quantized. Furthermore, a linear junction geometry is not required, as TAR has also been consistently demonstrated for Andreev point contacts.

Finally, even if quantized topological Andreev rectification cannot be observed in existing devices due to disorder, as long as the variation in chemical potential across the sample is not so drastic to create multiple regions with distinct Fermi sea topology (e.g., for graphene we would require the absence of electron-hole puddles), the signatures of Andreev rectification would still remain for $\chi_{F}\neq 0$. In particular, $\delta G_{21}$ is expected to change drastically as the quantized value of $\chi_{F}$ jumps under a Lifshitz transition of the underlying metal in the Andreev junction. In short, we believe that the topological Andreev rectification effect is experimentally testable in the near future.

Acknowledgements.

We thank M. Claassen for computational resources. This research was funded in whole, or in part, by the Luxembourg National Research Fund (FNR) (project No. 16515716). Work by P.M.T. and C.L.K. was supported by a Simons Investigator Grant to C.L.K. from the Simons Foundation.

Appendix A Landauer-Büttiker formalism

We consider a four-terminal setup with two normal leads (leads N1 and N2) and two superconducting leads (leads S1 and S2) as illustrated in Fig. 1 of the main text. Here the bias at the normal leads is given by $\mu_{1}=eV_{1}$ and $\mu_{2}=eV_{2}$ while the superconducting leads are grounded with $\mu_{S}=0$. When transport is coherent across the system, the current that flows into the normal leads is given by [55, 56]

[TABLE]

with $e=-|e|$ the electron charge and $f_{0}(\varepsilon)=1/\left[\exp\left(\varepsilon/k_{B}T\right)+1\right]$ the Fermi-Dirac distribution function, where the excitation energy $\varepsilon$ is defined relative to the equilibrium chemical potential $\mu$. Here we suppressed the subscript $N$ in the total (including spin) transmission functions between the normal leads

[TABLE]

with $i,j=1,2$ and $\alpha,\beta=e,h$. The matrices $s_{ij}^{\alpha\beta}$ contain the amplitudes for scattering from $\beta$-type incident modes in lead N$j$ to $\alpha$-type outgoing modes in lead N$i$. In the subgap regime $|\varepsilon|<\Delta_{0}$, the total number of electron modes in the normal leads is given by $M_{i}^{e}=\sum_{j}\left(T_{ij}^{ee}+T_{ij}^{eh}\right)$. Note that electron-hole symmetry implies $T_{ij}^{ee}(\varepsilon)=T_{ij}^{hh}(-\varepsilon)$ and $T_{ij}^{eh}(\varepsilon)=T_{ij}^{he}(-\varepsilon)$, while reciprocity implies $T_{ij}^{\alpha\beta}(\varepsilon)=\overline{T}_{ji}^{\beta\alpha}(\varepsilon)$ where $\overline{T}$ is calculated for the time-reversed system. The local and nonlocal differential conductance are then, respectively,

[TABLE]

where $f_{1}(\varepsilon)=f_{0}(\varepsilon-\mu_{1})$ with $\mu_{i}=eV_{i}$ ($i=1,2$). In the zero-temperature limit, we obtain

[TABLE]

The rectificed nonlocal conductance is

[TABLE]

where electron-hole symmetry and reciprocity have been used. For a transport setup that is both time-reversal symmetric and mirror symmetric along the transport direction (i.e., under $\mathcal{M}_{x}:\;x\mapsto-x$ for the setup shown in Fig. 1 of the main text), we have $\overline{T}_{12}^{he}=T_{12}^{he}=T_{21}^{he}$, and hence contributions from crossed Andreev reflections are canceled. The quantization in $\delta G_{21}$, as established in Sec. II.2, still relies on having a long junction with $L\gg\xi$ so that $T^{ee}_{21}$ and $T^{hh}_{21}$ only contain contributions from the dispersive ABSs (and not from tunneling).

Let us also compute the rectified local conductance $\delta G_{11}$ here, which is not studied explicitly in the main text. Using $M^{e}_{1}=T^{ee}_{11}+T^{he}_{11}+T^{ee}_{21}+T^{he}_{21}$ and electron-hole symmetry, we have

[TABLE]

As before, reciprocity and time-reversal symmetry $\mathcal{T}$ give $T^{he}_{11}=\overline{T}^{eh}_{11}=T^{eh}_{11}$, such that the first two terms in Eq. (51) cancel each other. Hence, when $\mathcal{T}$ is preserved, we have

[TABLE]

When, in addition to $\mathcal{T}$ symmetry, $\mathcal{M}_{x}$ is preserved, we have $T^{eh}_{21}=T^{he}_{21}$ and hence $\delta G_{11}=-\delta G_{21}$. If $\mathcal{M}_{x}$ is not present, but the junction is long enough such that crossed Andreev reflection is suppressed, we also have $\delta G_{11}\approx-\delta G_{21}$. Moreover, when $\mathcal{T}$ symmetry is broken but $\mathcal{M}_{y}\mathcal{T}$ is conserved, using reciprocity, the relations given in Eq. (52) and Eq. (53) still hold.

In this work we have focused on the nonlocal conductance $G_{21}$ instead of the local conductance $G_{11}$, as not only is $\delta G_{21}$ quantized to reflect the Fermi sea topology, in the case of a narrow Andreev junction $G_{21}$ is itself also quantized to reflect the Fermi sea geometry.