Thermoelectric detection of Andreev states in unconventional superconductors

We theoretically describe a thermoelectric effect that is entirely due to Andreev processes involving the formation of Cooper pairs through the coupling of electrons and holes. The Andreev thermoelectric effect can occur in ballistic ferromagnet-superconductor junctions with a dominant superconducting proximity effect on the ferromagnet, and it is very sensitive to surface states emerging in unconventional superconductors. We consider hybrid junctions in two and three dimensions to demonstrate that the thermoelectric current is always reversed in the presence of low-energy Andreev bound states at the superconductor surface. A microscopic analysis of the proximity-induced pairing reveals that the thermoelectric effect only arises if even and odd-frequency Cooper pairs coexist in mixed singlet and triplet states. Our results are an example of the richness of emergent phenomena in systems that combine magnetism and superconductivity, and they open a pathway for exploring exotic surface states in unconventional superconductors.


I. INTRODUCTION
Thermoelectric effects, where temperature gradients and electric voltages are converted into each other, open a promising way to reuse waste heat in electronic devices [1]. In mesoscopic conductors, a thermoelectric current can be generated by breaking the electron-hole symmetry around the chemical potential. However, the current is usually small at low temperatures. By contrast, large thermoelectric effects have recently been predicted [2][3][4] and measured [5,6] in hybrid junctions between superconductors and magnetic materials. The novel field of superconducting spintronics [7][8][9] explores such phenomena that arise due to the interplay between superconductivity and magnetism [10][11][12][13][14][15]. Interestingly, by coupling the spin degree of freedom to the heat (or charge) transport in ferromagnet-superconductor hybrid junctions, it may be possible to develop more efficient thermoelectric devices [16][17][18], such as coolers [3,4,19] and thermometers [20].
Here, we analyze the electrical current induced by a temperature gradient across hybrid junctions involving unconventional superconductors in two and three dimensions, see Fig. 1. The thermoelectric ef-fect has been measured in ferromagnet-superconductor tunnel junctions under an external magnetic field [5] or in proximity to a ferromagnetic insulator [6]. For this reason, we concentrate here on the transport properties of ballistic ferromagnet-insulatorsuperconductor (F-I-S) and ferromagnet-ferromagneticinsulator-superconductor (F-Fi-S) junctions described by scattering theory [22,25,50]. Building on earlier work on the spin and charge transport of F-I-S junctions with unconventional pairings [51][52][53][54][55][56][57][58][59][60], we show that the thermoelectric current is dominated by Andreev processes, where electrons and holes are converted into each other at the surface of the superconductor. Interestingly, the thermoelectric current changes sign in the presence of SABS: it is positive (it runs into the superconductor) when the conductance features a zero-bias peak and it is negative otherwise, as illustrated in Fig. 1(b,c).
Earlier predictions of the thermoelectric effect in ferromagnet-superconductor junctions require that the exchange field from the magnet penetrates the superconductor, and thereby spin-splits its density of states [21]. Such a magnetic (or inverse) proximity effect can be achieved in thin films of conventional superconductors [15,61], which could be challenging to implement with unconventional ones. Our results, by contrast, rely on the coupling between the ferromagnet's exchange field and the Cooper pairs leaking out of the superconductor. When the superconducting proximity effect takes place [62], quasiparticles can tunnel into the superconductor only when the pairing rotational symmetry is broken at the interface. However, this additional contribution to the thermoelectric current is always smaller than the Andreev one for pairings featuring zero-bias peaks, and hence it does not affect the direction of the current flow. An analysis of the proximity-induced pair amplitude at the junction interface, based on Green's function techniques, allows us to show that the Andreev thermoelectric current is only possible when even-and oddfrequency Cooper pairs coexist [63][64][65][66][67]. The thermoelectric effect is thus a signature of proximity-induced oddfrequency pairing, and it thereby becomes a sensitive tool to explore unconventional superconductivity. The rest of the paper is organized as follows. We introduce the theoretical methods in Section II and explain the Andreev thermoelectric effect and the current sign change in Section III. In Section IV, we use Green's function techniques to analyze the proximity-induced pairing on the ferromagnetic region. Section V discusses the quasiparticle contribution to the thermoelectric current. Finally, we present our conclusions in Section VI.

II. THEORETICAL FORMULATION
We consider ballistic junctions between a superconductor and a magnetic normal-state metal separated by an insulating barrier. Transport takes place in the xdirection, with the superconductor extending over the region x ≥ 0. In all cases below, we assume perfectly flat and clean interfaces. Transport in such junctions is described by the Bogoliubov-de Gennes (BdG) equationš with E being the excitation energy. Particle-hole symmetry imposes that every solution Ψ E of the BdG equations with excitation energy E is accompanied by another solution, Ψ −E , at energy −E. We then choose the basis in Nambu (particle-hole) and spin space, withĉ σ (k) being the annihilation operator for an electron with spin σ =↑, ↓ and momentum k. The general form of the BdG Hamiltonian in any region iš and energy (bottom) schematics of an Andreev reflection process. Assuming k ↑ > kF > k ↓ , Andreev reflection is only possible for angles of incidence below a critical value |θ| < θc, virtual Andreev reflections (VAR) take place for θc < |θ| < θc2, and all incident particles are backscattered for |θ| > θc2.
whereσ 0,1,2,3 are the usual Pauli matrices in spin space and µ and∆ the chemical and pair potentials, respectively. It is reasonable to assume that these quantities are constant in each region. For the sake of simplicity, we take the same value of the chemical potential on both sides of the junction. The pair potential is zero in the normal region, which is a good approximation if the superconducting coherence length is much larger than the Fermi wavelength. The superconducting phase φ only plays a role in junctions with several superconductors, hence, we omit it here. In the following we consider systems in either two or three dimensions with the single-particle Hamiltonian where m is the electron mass and U (x) = U 0 Θ(−x) is the exchange potential in the normal-state region with U 0 ≥ 0 and Θ(x) being the Heaviside step function. The potential V 0 describes the insulating barrier at the interface with U B accounting for the exchange potential of ferromagnetic insulators. Assuming translational invariance in the plane perpendicular to the direction of transport, r = (y, z), the momentum parallel to it, k = (k y , k z ), becomes a good quantum number. For a given incident wave, all related wave vectors lie in the same plane, see Fig. 2(a). We can then choose the coordinate system so that transport depends only on the azimuthal angle θ = sin −1 (|k |/k F ), where k F = √ 2mµ/ is the Fermi wave vector. Inserting Eq. (4) into Eq. (3), the BdG equations can be decoupled into two spin channels if the pair po-tential∆ is either a diagonal or an off-diagonal matrix, i.e., proportional toσ 0,3 orσ 1,2 , respectively. The spin-singlet state,∆(k) = ∆ 0 (k)(iσ 2 ) always fulfills this condition. The triplet state, which we parametrize aŝ ∆(k) = d(k) ·σ(iσ 2 ), withσ being a vector of Pauli matrices, fulfills the condition if d lies either along the z-axis or in the x-y plane. In the following, we only consider pair potentials that fulfill these conditions and that couple electrons with opposite spins. Therefore, the injected and Andreev reflected particles have opposite spins and, if U 0 = 0, Andreev retro-reflection does not take place [52]. Following a quasiclassical approximation, where the magnitude of the wave vectors is evaluated at the Fermi surface, translational symmetry at the interface entails that where k σ = k F √ 1 + ρX, with X = U 0 /µ being the polarization from the exchange field and ρ = 1 (−1) and σ =↓ (↑), when σ =↑ (↓). Consequently, assuming upspin injection and X > 0, total (normal) reflection takes place for θ > sin −1 (k F /k ↑ ) ≡ θ c2 , when k F < k ↑ , and the Andreev-reflected wave becomes evanescent for the an- The latter case corresponds to virtual Andreev reflections [52], which do not contribute to the current, see Fig. 2. For a given spin channel σ, the general solution for the decoupled BdG equations is given by where we have defined with m α σ (m α σ ) being the amplitude for an electron-like (hole-like) quasiparticle with direction of propagation α = ±. Here, θ + = θ and θ − = π − θ label right and left movers, respectively, and we have defined η α σ = ∆ α σ /|∆ α σ |, Γ α σ = |∆ α σ |/(E + Ω α σ ), and where k α e(h),σ = k F 1 + (−)Ω α σ /µ + (−)ρX and ∆ α σ = ∆ σ (θ α ) is the pair potential felt by an α-moving particle in spin channel σ. Figure 2(a) shows the scattering of an electron incident from the normal-state region into an Andreev reflected hole (AR), a normal reflected electron (NR), an electron-like quasiparticle transmitted to the superconductor (ELQ), or a hole-like transmitted quasiparticle (HLQ). The reflection amplitudes are obtained evaluating the wave function in Eq. (7) for this scattering process and imposing the boundary conditions witĥ When solving Eq. (9), we follow the Andreev approximation, discarding terms of order E/µ and ∆ 0 /µ, resulting in k α e(h),σ ≈ k F in the superconducting region (x > 0) and k α e(h),σ ≈ k σ(σ) for x < 0. Accordingly, for the scattering of an incident electron with spin σ the angles in Eq. (7) become θ = θ = θ S , for x > 0, and θ and θ = θ σ A for x < 0.
The angle-resolved conductance per spin reads [50,52] where a σ (E, k ) and b σ (E, k ) are the Andreev and normal reflection amplitudes, respectively, and we have defined λ σ = Re{kσ cos θ σ A /(k σ cos θ)} and θ σ A = sin −1 (k σ sin θ/kσ). By averaging over all incident angles, the normalized conductance becomes and the normal state transmissivity is with λ σ N = k σ cos θ/(k F cos θ S ) and P σ = (1 + ρX)/2. The angle averages are given by [68] for two and three dimensions, respectively. Finally, the electrical current is where e is the electron charge and δf (E, δT, identify the two main contributions to the current coming from Andreev processes (2I A ) and quasiparticle transfers (I qp ), respectively. If the junction is subjected to a temperature gradient δT in the absence of an applied voltage, V = 0, the difference of the Fermi distributions in each contact is an odd function of the energy E, i.e., δf (E, δT, 0) = −δf (−E, δT, 0). The thermoelectric current is then the integral of Eq. (12) multiplied by an odd function of the energy. Therefore, only the part of Eq. (12) that is odd in the energy will contribute to the integral. Using the anti-symmetry of the difference of Fermi functions, we then find , as in the case of conventional normal metal-insulator-superconductor (N-I-S) junctions. On the other hand, a thermoelectric current can be induced by an asymmetry in the conductance, as we explain in the following section.

III. ANDREEV THERMOELECTRIC EFFECT
To illustrate the effect of SABS on the thermoelectric current we start by considering a singlet (d xy + is)-wave pairing,∆ with ∆ 1 + ∆ 2 = 1. Henceforth, the pairing temperature dependence is included in the amplitude by the interpolation formula ∆ T /∆ 0 = tanh(1.74 T c /T − 1). This pair potential was proposed as a fragile state [49,69] that could explain the splitting of the zero-bias peak in high-T c superconductors [48]. We choose it because the parameter ∆ 1 = 1 − ∆ 2 controls the evolution from a gap profile (s-wave, ∆ 1 = 0) into a zero-energy state (d xywave, ∆ 1 = 1). The tunnel spectroscopy of the (d xy + is)wave state, Fig. 3(a), reveals that the zero-energy peak in the conductance from the d xy -wave pairing splits in the presence of a small s-wave component until the Andreev resonances reach the gap edges at E = ±∆ 0 , when ∆ 1 = 0. The splitting of the zero-bias peak becomes asymmetric with the energy for an F-I-S junction with X = 0. The energy-antisymmetric part of the resulting conductance, Fig. 3(b), presents a Fano-like sign change at the resonance positions. While well-separated resonances feature two Fano-like line shapes, one for positive and another for negative energies, there is only one Fano-like resonance when the SABS merge around zero energy (∆ 1 1). This different behavior for separated and merged resonances explains the sign change in the current, Fig. 3(c,d). For a merged resonance, the anti-symmetric contribution to the conductance around zero energy follows the same profile as the difference in Fermi distributions δf [indicated by a gray dotted line in Fig. 3(b)]. Consequently, the thermoelectric current is positive after integrating over the energy. By contrast, for separated resonances located at E = ±δE, with δE > 0, the anti-symmetric part of the conductance has the opposite sign to δf in the energy range |E| < δE. For experimentally relevant temperatures [5,6], T T c /2 and δT 0.3T c , this energy window provides the biggest contribution to the integral, and the resulting current is thus negative, see Fig. 3(c,d). The transition between these two regimes depends on the energy splitting between the SABS, with δE ∆ 0 /4 enough to avoid a temperature-induced sign change.
We stress here that the asymmetry in the conductance only occurs at subgap energies, where the normal and Andreev probabilities for each spin channel, B σ = |b σ | 2 and A σ = λ σ |a σ | 2 , respectively, cf. Eq. (11), fulfill B σ = 1 − A σ . Consequently, the conductance is only due to Andreev processes, i.e., σ = 2R N P ↑ A ↑ + P ↓ A ↓ k , and the resulting thermoelectric current in Eq. (15) is Andreev-dominated, that is, We have thus demonstrated the Andreev thermoelectric effect for an F-I-S junction with X = 0 and (d xy + is)wave pairing and how the current goes from positive to negative in the presence or absence, respectively, of lowenergy SABS. Within our approach, the necessary and sufficient condition for the Andreev thermoelectric effect is to break time-reversal and particle-hole symmetries per spin species so that the tunnel conductance becomes asymmetric with respect to the energy, cf. Eq. (16). This can be achieved for any pairing state in an F-I-S junction with X = 0, as a result of the different polarization of each spin channel shown in Eq. (12). Indeed, the Andreev thermoelectric effect can occur in well-known unconventional superconductors. The (d + is)-wave pairing of our example explains a specific state of high-T c superconductors, which, more generally, feature an anisotropic d-wave singlet pair potential [22,70]. Spin-triplet p-wave states have been suggested as an explanation of the properties of UPt 3 , Sr 2 RuO 4 and other Ru-based compounds [30][31][32]. Some of these materials can develop topological phases, with the emergence of chiral states [32]. Most of these cases can be described by the pairing potentials listed in Table I. Some threedimensional chiral superconductors feature more complicated potentials that we detail later. In this section, we assume that the pairing rotational symmetry is preserved at the junction interface. We also define the effective pairing potentials for the right-moving electron-and holelike quasiparticles in the superconducting region as ∆ + σ and ∆ − σ , respectively, see Fig. 2(a). Rotational symmetry guarantees that |∆ + σ | = |∆ − σ |, and transmitted quasiparticles with the same spin thus experience the same pairing.
The spectral differential conductance, Fig. 4, illustrates the effect of the exchange fields on the subgap Andreev states for the four types of junctions considered here. For a tunnel N-I-S junction, it reveals the dispersion of the superconducting pairing. For example, d x 2 −y 2 -wave features a nodal gap without SABS, while d xy -wave presents the same anisotropy but displays zero-energy states, cf. Fig. 4(a,b). The chiral p-and dwave states, Fig. 4(c) and Fig. 4(d), respectively, have SABS with linear dispersions, but only the chiral p-wave state crosses the zero energy point at zero momentum. In all cases, a ferromagnetic insulating barrier (N-Fi-S junction) causes a spin-splitting of the spectrum, very weak for the gapped d x 2 −y 2 -wave case, but pronounced if SABS are present. The splitting, however, does not cause any asymmetry. Consequently, there is no current flowing from a temperature bias in an N-Fi-S junction.
A spin-active barrier alone is not enough to generate a thermoelectric effect; breaking time-reversal and particlehole symmetries per spin species requires at least an F-I-S junction. When X = 0, the angle of incidence θ for each spin channel becomes different, see Fig. 2(b). The maps in Fig. 4 show the superposed contribution of both channels, Eq. (11), revealing a strong angular asymme-try. In F-Fi-S junctions, the spin-active barrier further splits the spin bands, which, in addition to the angular asymmetry, greatly affects the thermoelectric current. If the spin-active barrier favors the majority spin species, in our case when Z m > 0 if X > 0, the asymmetry at subgap energies is enhanced for all pairings, see Fig. 4. The resulting current is enhanced by an order of magnitude. Otherwise, when Z m < 0 for X > 0, the current flow is reversed, see Fig. 5. The resulting thermoelectric currents still have opposite signs for gapped pairings and those featuring SABS, but their magnitude is smaller than for the case with parallel exchange fields. Notice that the crossover point takes place at Z m < 0, since there is a finite current at Z m = 0, cf. Fig. 5(a). Fig. 5 shows how the sign of the thermoelectric current is positive (negative) for superconducting pairings featuring low-energy SABS (gap) in the conductance. To explain the direction of the current flow for the chiral por d-wave pairings, we need to take into account that, for a given angle θ, the SABS can be located at any subgap energy, see Fig. 4(c,d). After the angle-average, however, the biggest contributions come from the angles around normal incidence (θ ∼ 0). Therefore, the chiral p-wave pairing follows the behavior of a merged resonance, with a positive current, while the chiral d-wave case behaves like a separated one and thus provides a negative current. This different behavior of the chiral states is not always captured by the tunneling conductance, especially for 3d chiral superconductors [41].
The sign of the Andreev thermoelectric current thus pinpoints the presence or absence of low-energy Andreev states. If the SABS are exactly at zero energy, like for p xand d xy -wave cases, the positive current is usually one or two orders of magnitude stronger than the negative one for gapped pairings, see Fig. 4. We note here that the negative current for gapped superconductors is already much larger than the expected thermoelectric current for equivalent, non-superconducting junctions [2][3][4] [71].

IV. ODD-FREQUENCY PAIRING
We now show how a finite thermoelectric current is only possible if even and odd-frequency Cooper pairs coexist with similar amplitude. The Andreev reflection amplitude is related to the proximity-induced pairing, which can be described by the anomalous part of the Green's function associated to Eq. (3). Following Refs. [22,[72][73][74][75][76][77][78][79][80][81] (see all the details in the appendix), we can write the anomalous Green's function aŝ where we have identified the spin-symmetric singlet and triplet components f 0 and f 3 , respectively, as Here, a S,T are the singlet and triplet Andreev reflection amplitudes (a 1σ corresponds to an injected electron in spin channel σ) and we have approximated the wave vectors as k σ e,h ≈ k σ ± E/( v σ ), with k σ = k F cos 2 θ S + ρX and v σ = k σ /m. It simplifies the symmetry analysis to use as reference the angle in the superconductor, θ S , which is the same for both spin channels, and only consider the anomalous Green's function evaluated at the interface, where x = x = 0.
Next, we define the even and odd parity parts of the anomalous functions as f r µ, x, E, −θ S )]/2, with µ = 0, 3 and the sign change in the angle due to the exchange y ↔ y . Then, to study the frequency dependence of the anomalous function, we use the retarded and advanced functions for positive and negative frequencies, respectively [78], since only there these functions have physical meaning. Consequently, the anomalous function connects both energy branches. Once the anomalous function in Eq. (21) is fully symmetric with respect to spin and parity, it belongs to one of four possible symmetry classes. We label the different classes by their spin state, S for singlet and T for triplet, and their symmetry, even (E) or odd (O), with respect to frequency and space coordinates. For example, even-frequency, singlet, even-parity is ESE and ETO refers to the even-frequency, triplet odd-parity case [63,[65][66][67].
It is now useful to define the short-hand notation a α β = (αE, βθ), with α, β = ± and a standing for either of a S,T , and thus write the energy-and angle-symmetric parts of the amplitude a as Here, the index a E,O (a e,o ) labels the even and odd parts of the amplitude with respect to the energy (angle). From the particle-hole symmetry of the Hamiltonian in Eq. (3) and the micro-reversibility of the scattering matrix, we obtain the condition v ↓,↑ a 1↑,↓ (E, θ) = −v ↑,↓ a * 1↓,↑ (−E, −θ) for |θ S | ≤ θ c = sin −1 √ 1 − X, see the appendix for more details. Once applied to the singlet and triplet amplitudes, this condition reads a S,T (E, θ) = ±a * S,T (−E, −θ), and combined with Eq. (22) reveals that, in the absence of virtual Andreev processes (|θ S | ≤ θ c ), the amplitudes a Ee S , a Oo S , a Eo T and a Oe T are real, while the rest are purely imaginary. Therefore, at the interface we find where we substituted the purely imaginary amplitudes as a = iā, withā being real. It is then straightforward to check that f ESE and f OSO in Eq. (23) are, respectively, even and odd in both frequency and parity [82].
We can now establish a connection between the thermoelectric current and the symmetry of the pairing state. Focusing on the Andreev contribution to Eq. (16), which we have demonstrated to be equal to the total current, we find where we have used that λ σ = vσ/v σ . It is clear from Eq. (24) that a finite polarization (X = 0) couples the spin singlet and triplet states, resulting in a new contribution to the current. To identify the scattering amplitudes with the anomalous Green's function, we first note that the virtual Andreev processes do not contribute to the current and thus the angle average is limited to |θ S | ≤ θ c . Then, using Eqs. (22) and (23) we can write the scattering amplitudes in the absence of virtual processes as Inserting Eq. (25) into Eq. (24), the current becomes In general, the angle-average of the product of functions of opposite parity in the same spin state is zero, i.e., f * ESE f OSO k = f * ETO f OTE k = 0. In F-I-S junctions for arbitrary pairings, see Fig. 5, the thermoelectric effect is then due to the term proportional to X, which mixes even-and odd-frequency terms with different spin states [83]. Therefore, the Andreev thermoelectric effect is only possible if even and odd-frequency pairing coexist [84].

V. BROKEN ROTATIONAL SYMMETRY
In all the previous results, the Andreev thermoelectric effect is due to subgap processes where only Cooper pairs are transmitted into the superconductor, see Eq. (18), for any value of X or Z m . However, the transport properties of unconventional anisotropic pairings depend strongly on the relative orientation between the pairing potential and the junction interface. For example, for d xy - wave pairing, the effective pairings for electron-and holelike quasiparticles are ∆ ± σ = ±ρ∆ 0 sin (2θ − 2α), respectively, where α is the misalignment angle [25,26]. When α is different from one of the highly symmetric orientations, α = 0, π/4, the rotational symmetry for d xy -wave is broken and |∆ + σ | = |∆ − σ |, resulting in a small contribution to the thermoelectric current from quasiparticles.
We show the effect of the misalignment angle in an F-I-S junction with d xy + is-wave pairing in Fig. 6(a). When α = 0, π/4, the transport properties change. For normalized energies E/∆ 0 below ∆ 2 = 1−∆ 1 , the system is fully gapped, the condition σ = 2R N P ↑ A ↑ + P ↓ A ↓ k holds and the antisymmetric part of the conductance is only due to Andreev processes. However, for ∆ 2 < |E/∆ 0 | < 1 quasiparticles can be transferred into the superconductor. Consequently, the subgap condition A σ + B σ = 1 no longer holds and the quasiparticle term I qp becomes finite in Eq. (15). Even though the current now has two contributions, I A and I qp , we are clearly in the situation where I A > I qp , see inset of Fig. 6(a). For the two-dimensional pairings considered here, the quasiparticle current generated by a broken rotational symmetry is not enough to change the sign of the total current around the highly symmetric orientations.
The results from the previous section can be extended to three-dimensional (3d) unconventional superconductors without loss of generality. In higher dimensions, however, the energy dispersion of SABS can become more complicated [85,86]. Recently, the exotic SABS of 3d chiral superconductors have attracted a lot of attention due to the simultaneous presence of line and point nodes [87][88][89][90][91][92][93][94][95][96]. The SABS of 3d chiral superconductors are predicted to originate from the pair potential with k x = k x cos α − k z sin α, k z = k x sin α + k z cos α, α being the misalignment angle from the k x axis, ν a nonzero integer [40,41], and the normalization factors (r ν=0 , r ν=1 , r ν=2 ) = (1, 1/2, 2/ √ 27). For ν = 0, the potential features a line node, while it has two point nodes and a line node for ν = 1, 2. The pairing is in a spintriplet (singlet) state for ν even (odd).
While the case with ν = 0 generalizes the 2d chiral p-wave triplet pairing of the previous section, the cases with finite ν are relevant for heavy fermion superconductors, with ν = 1 and ν = 2 corresponding to the candidate pairing symmetries of URu 2 Si 2 [97][98][99] and UPt 3 [100][101][102][103][104], respectively. The topological origin of the flatband SABS and its fragility against surface misalignment have been clarified in previous works [40], together with the SABS energy dispersion [41]. It was found that a unique characteristic of the ν = 1, 2 pairings is that an external exchange field can turn a zero-bias conductance peak into a dip, due to the complicated dispersion of the SABS [41]. Consequently, an equivalent sign change in the thermoelectric current appears here.
The thermoelectric current for 3d chiral superconductors behaves like the two-dimensional cases for highly symmetric orientations (α = 0, π/2). For different surface misalignments, however, the quasiparticle contribution from broken rotational symmetry can affect the sign of the current, see Fig. 6(b). In the presence of a zeroenergy SABS, however, the quasiparticle contribution is smaller than the Andreev one and the current remains positive. The sign change of the thermoelectric current thus clearly indicates the presence of SABS, even when they have a complicated energy dispersion. Moreover, some orientations, like α = π/2, have a very different zero-energy behavior depending on the value of ν. For ν = 1 and ν = 2, the pairing can be regarded as a chiral p-wave or chiral d-wave, which feature the opposite signs of the current. Consequently, the sign of the thermoelectric current is, again, a good indicator of the presence of SABS, even if the surface states display a complicated dispersion relation.

VI. SUMMARY
Following a microscopic scattering approach, we have shown that the thermoelectric effect in ferromagnetsuperconductor junctions can be entirely dominated by subgap Andreev processes. Consequently, the electric current from a temperature bias changes sign in the presence of SABS from unconventional superconductivity and is thus a sensitive tool to identify these surface states.
When the rotational symmetry of the pairing potential is maintained at the junction interface (that is, |∆ + σ | = |∆ − σ | for each spin channel σ) the transport asymmetry that gives rise to the thermoelectric effect occurs only at subgap energies. The resulting current is dominated by Andreev processes, I(δT ) = 2I A (δT ), with no quasiparticle transfer into the superconductor. A finite polarization induced by a Zeeman exchange field X in an F-I-S junction is a necessary and sufficient condition for this Andreev thermoelectric effect. An additional intermediate ferromagnetic insulator, F-Fi-S junction, enhances the magnitude of the current if the majority spin species is favored (i.e., Z m > 0 for X > 0). We have shown how the exchange field allows for the coupling of even-and odd-frequency Cooper pairs leaking out of the superconductor, which provide the main contribution to the Andreev thermoelectric current. For F-Fi-S junctions, we restricted our analysis to the case of collinear exchange fields in the ferromagnet and ferro-magnetic insulator. For non-collinear fields, a long-range, odd-frequency triplet state emerges [11][12][13]105], which could enhance the thermoelectric current according to Eq. (26).
When the pairing rotational symmetry is broken (|∆ + σ | = |∆ − σ |), quasiparticle transfers into the superconductor are allowed and contribute to the thermoelectric current. Such contributions could change the sign of the current and mask the presence of SABS. However, we have analyzed the thermoelectric effect for two-and three-dimensional pairings with a strong angular dependence and when low-energy SABS are present the current is positive, Andreev-dominated and with a strong magnitude around highly-symmetric orientations. As such, the sign of the thermoelectric current could help identify the pairing state even for cases where the outcome of tunneling spectroscopy can be ambiguous.
In conclusion, exploring the thermoelectric effect in ferromagnet-superconductor hybrid junctions involving unconventional superconductors offers a sensitive method to identify emerging surface Andreev states. Our results can be easily extended to heterostructures, where topological superconductivity is artificially engineered by combining conventional superconductors with topological materials [106].  [22,[72][73][74][75][76][77][78][79][80][81], the retarded Green's function associated with the Hamiltonian in Eq. (3) is constructed from a linear combination of scattering states obtained from Eq. (7) with the appropriate boundary conditions. We now describe the main steps of this method for a generic pairing potential.
First, since the spin channels in Eq. (3) are decoupled, we only need to work with the reduced Hamiltonian where The retarded Green's function thus satisfies the equation of motion withÊ + = (E + i0 + )σ 0 , which can be integrated to obtain the following boundary conditions: A general solution of the BdG equations forĤ σ (x, k) takes the form is the solution for an electron-like (hole-like) quasiparticle with amplitude m α σ (n α σ ), cf. Eq. (7). For the transposed Hamiltonian, we find Notice that both the extra conjugation and the exchange θ + ↔ θ − are a result of the transposition (we take the superconducting phase to be zero for the sake of simplicity). The four possible scattering states for an F-I-S junction are injection of electron (1) or hole (2) from the ferromagnet and injection of an electron-(3) or hole-like (4) quasiparticle from the superconductor. Specifically, we have where the wave vectors in the normal regions are k σ e,h and thus ψ ±k σ e = (1, 0) T and ψ ±k σ h = (0, 1) T , since η ± σ = 0. Finally, the retarded Green's function is a linear combination of asymptotic solutions ofĤ σ (x, k) andĤ T σ (x, k), namely, The coefficients of the linear combination are obtained by imposing the boundary conditions in Eq. (A4). Focusing on the region with x, x < 0, we find from the Nambu-diagonal elements the system of equations resulting in with v σ e,h = k σ e,h /m being the velocity for spin σ electrons and holes, respectively. On the other hand, the Nambu off-diagonal sector provides the detailed balance equations which can be demonstrated using the Wronskian of the scattering processes 1 and 2. For Nambu spinors ψ j (x) = [u j (x), v j (x)] T , the Wronskian is defined as with u = du(x)/dx. For the scattering processes 1 and 2, we obtain The Wronskian is independent of the coordinate x and, thus, equating the two branches of the Wronskian, we obtain one of the detailed balance equations in Eq. (A8). The other one comes from W [φ σ,1 (x), φ σ,2 ]. By plugging the linear coefficients in Eq. (A7) into Eq. (A6), and using the detailed balance relations in Eq. (A8), we obtain the retarded Green's function per spin channel in the normal region aŝ Appendix B: Symmetry analysis of the anomalous Green's function Using as reference the angle in the superconductor, θ S , which does not change with the spin channel, we approximate the wave vectors in the normal region as follows with k σ = k F cos 2 θ S + ρX and v σ = k σ /m. The spatial dependence of the anomalous part is thus approximated as The anomalous part of the spinful Green's function is then given by the functionŝ Next, we compute the advanced anomalous Green's functionĜ a he fromĜ r eh using the propertyǦ a (r, r ) = [Ǧ r (r , r)] † , resulting inĜ At the interface, x = 0, the anomalous Green's functions reduce tô From the detailed balance relations within a spin sector, see Eq. (A10), we find that, if with the minus sign for singlet and plus for triplet. Particle-hole symmetry, for |θ S | ≤ θ c = sin −1 √ 1 − X, imposes that which changes the spin sector. Applying this property to Eq. (B5) relates the retarded and advanced anomalous functions asĜ a he (x, which connects the same scattering process for different spin sectors. Applying these relations to the anomalous Green's function at the interface, we obtain From Eq. (B3) and Eq. (B5), we define the singlet and triplet components of the anomalous Green's functions as The retarded and advanced functions are related by the transformation f a 0,3 (x, x ; E, θ S ) = ±f r 0,3 (x , x; −E, −θ S ), cf. Eq. (B11), which connects the positive and negative energy branches, where these functions have physical meaning. To investigate the frequency dependence of the anomalous function, we define a function equal to the retarded (advanced) anomalous function for positive (negative) frequencies [78,80]. Assuming > 0, we write the energy as E = sgn(E) to find where the sign for negative energies takes into account the different spin-symmetrization of the advanced Green's function. Next, we define the even and odd parity parts of the anomalous function as f µ,± (E) = [f µ (x, x , E, θ S ) ± f µ (x , x, E, −θ S )]/2, with µ = 0, 3. Once the anomalous function in Eq. (B12) is fully symmetric with respect to spin and parity, it belongs to one of four possible symmetry classes. We label the different classes by their spin state, S for singlet and T for triplet, and their symmetry with respect to frequency and space coordinates, even (E) or odd (O). For example, odd-frequency, singlet, odd-parity is OSO and OTE refers to the odd-frequency, triplet even-parity case. Consequently, the spin-singlet, even parity anomalous function reads which is the same for positive and negative energies. One can also check that f 0, We now focus on the local anomalous function with x = x ≡ x 0 . We define the singlet and triplet Andreev scattering amplitudes as a S,T = [a ↑ 1 /v ↑ ∓ a ↓ 1 /v ↓ ]/2 and, using the short-hand notation a α β = (αE, βθ S ), with α, β = ± and a standing for either of a S,T , we can define the energy-and angle-symmetric parts of the amplitude a as Here, the index a E,O (a e,o ) labels the even and odd parts of the amplitude with respect to the energy (angle). We can now apply the condition in Eq. (B9) to the singlet and triplet amplitudes and find that a S (E, Expressing the singlet and triplet amplitudes in terms of the energy and angle-symmetric ones we find Then, using Eq. (B15) and the fact that the symmetric amplitudes are linearly independent, we find that and these amplitudes are thus real. Analogously, we have that and these amplitudes are thus purely imaginary. Therefore, when necessary, we can change the purely imaginary amplitudes as a = iā, withā real. We show an example of the angle dependence of the fully-symmetric scattering amplitudes in Fig. 7. The region between the dashed vertical gray lines corresponds to |θ S | ≤ θ c , where Eqs.
and the thermoelectric current is where we have used that a ↑,↓ 1 = v ↑,↓ (a T ± a S ) from the definition of the singlet and triplet components, λ σ = vσ/v σ , and P ↑(↓) = [1 + (−)X]/2. It is clear in Eqs. (C1) and (C2) that a finite polarization (X = 0) couples the spin singlet and triplet states, resulting in a new contribution to the current.
We can now identify the scattering amplitudes with the anomalous Green's function described above. First, we note that the virtual Andreev processes do not contribute to the current and thus the angle average is limited to |θ S | ≤ θ c . Using Eqs. (B14), (B20) and (B21) we can write the scattering amplitudes in the absence of virtual processes as follows: a S (E) = a Ee S + a Oe S + a Eo S + a Oo S = (a Ee S + iā Oe S ) + (a Oo S + iā Eo a T (E) = a Ee T + a Oe T + a Eo T + a Oo T = (a Oe T + iā Ee T ) + (a Eo where we have used that the amplitudes a Ee S , a Oo S , a Eo T and a Oe T are real, while the rest are purely imaginary and thus can be recast as a = iā, cf. Eqs. (B20) and (B21).
Inserting Eq. (C3) into Eqs. (C1) and (C2), the Andreev conductance reads and the current becomes The angle-average of the product of functions of opposite parity in the same spin state should be zero, i.e., f * ESE f OSO k = f * ETO f OTE k = 0, unless the rotational symmetry is broken by α = 0. That is not the case for the products of different spin states, f * ESE f OTE θ S and f * ETO f OSO θ S , which then result in the thermoelectric current for F-I-S junctions even when α = 0.