Apr 15, 2013 - In this work, we revisit the problem of a diffusive SNS .... H = He e. Ly. Lx. N x=0 y=0. FIG. 1: SNS junction in a transverse magnetic...

2 downloads 44 Views 432KB Size

Dmitri A. Ivanov

arXiv:0907.0632v3 [cond-mat.supr-con] 15 Apr 2013

Institute for Theoretical Physics, ETH Zurich, 8093 Zurich, Switzerland and Institute for Theoretical Physics, University of Zurich, 8057 Zurich, Switzerland We study a diffusive superconductor–normal metal–superconductor (SNS) junction in an external magnetic field. In the limit of a long junction, we find that the form of the dependence of the Josephson current on the field and on the length of the junction depends on the ratio between the junction width and the length associated with the magnetic field. A certain critical ratio between these two length scales separates two different regimes. In narrow junctions, the critical current exhibits a pure decay as a function of the junction length or of the magnetic field. In wide junctions, the critical current exhibits damped oscillations as a function of the same parameters. This damped oscillating behavior differs from the Fraunhofer behavior typical for short or tunnel junctions. In wide and long junctions, superconducting pair correlations and supercurrent are localized along the edges of the junction.

I.

INTRODUCTION

Interplay between magnetism and superconductivity in proximity structures results in rich physics.1,2 One of the remarkable effects arising from magnetism in Josephson junctions is the possibility of a reversal of the critical current (π coupling). It was theoretically predicted3 and then experimentally observed4 in Josephson junctions with a ferromagnetic interlayer. In those junctions, π coupling results from the effective Zeeman field induced by the exchange interaction in the ferromagnet. In the present paper, we revisit another way to reverse the critical current in a Josephson junction by applying an external magnetic field. In this case, the main role of the field is to provide a vector potential coupled to the motion of electrons. This effect was studied in detail in the case of a thin (tunnel) junction. If a transverse magnetic field is applied to such a junction, then the Josephson coupling oscillates along the direction transverse to the field, so that the total current exhibits Fraunhoferlike oscillations as a function of the magnetic flux through the junction.5,6 It has been shown both theoretically7 and experimentally8,9 that the Fraunhofer interference pattern does not appear in narrow and long junctions. In particular, the authors of Ref. 7 have shown numerically how the damped oscillatory behavior (Fraunhofer like) characterizing wide and short junctions is replaced by a monotonic exponential decay in narrow diffusive junctions. They have also identified the length scale over which the transition between the two regimes takes place. In the clean case, studies of similar SNS junctions in a magnetic field were performed in Ref. 10 and, in a slightly different geometry, in Ref. 11. In this work, we revisit the problem of a diffusive SNS junction in a magnetic field considered in Ref. 7 and show that further progress may be achieved in the limit of a

long junction using analytical methods. The external magnetic field introduces a “magnetic” length scale r φ0 , (1) ξH = H where φ0 = h/2e is the (superconducting) flux quantum and H is the external magnetic field. We confirm the findings of Ref. 7 that the structure of the Josephson current depends on the relation between the width of the junction and ξH . For narrow junctions, the critical Josephson current decays exponentially as a function of the length of the junction (or, equivalently, of the magnetic flux through the junction). For wide junctions, the critical Josephson current exhibits damped oscillating behavior (reversing sign), as a function of either the length of the junction or of the total magnetic flux through it. The period of these oscillations is the same as in the Fraunhofer limit but the decay is exponential instead of algebraic. In this wide-junction regime, the superconducting correlations and the supercurrent are localized near the edges of the junction. We find analytically the transition point (the critical width of the junction) where the monotonic decay of the critical current switches over to damped oscillations. The paper is organized as follows. In Section II we describe the model of a diffusive SNS junction in a transverse magnetic field. We compute the superconducting Green function in Section III and the Josephson current in Section IV. In Section V we summarize our results. Several technical details are relegated to the Appendices.

II.

SNS JUNCTION IN A TRANSVERSE MAGNETIC FIELD

We consider a SNS junction in a transverse magnetic field in a quasi-two-dimensional geometry. A coordinate

2 system is introduced in such a way that the x axis is directed along the junction, the field is applied along the z axis, and the SN interfaces are parallel to the yz plane (Fig. 1). The SN interfaces correspond to the coordinates x = 0, Lx . The origin of the y axis is chosen in the middle of one of the SN interfaces and we denote the width of the junction by Ly (so that the edges of the junction correspond to y = ±Ly /2). The system is translationally invariant along the z direction. We assume a uniform magnetic field H directed along the z axis and neglect the screening of the magnetic field by the Josephson currents (which assumes that the dimensions of the junction in the y and z directions are much shorter than the Josephson penetration length).7 For simplicity, we further assume that the London length is short compared to other length scales in the problem and neglect the penetration of the magnetic field in the superconducting electrodes. Then the vector potential of the field may be chosen as A = −yHex , with the phase of the superconducting order parameter (denoted ±χ below) constant along the SN interface in each of the superconducting leads (see Appendix A). This choice of the vector potential was used previously in Ref. 7,12 and is exactly gauge equivalent to the more conventional gauge with A parallel to the SN interfaces considered in Ref. 10. We further assume that the normal layer is strongly disordered and the motion of electrons is diffusive. In this regime, the quasiclassical Green functions (averaged over the fast Fermi oscillations and the momentum directions) are given by the solutions to the Usadel equations:13,14 ˆ gˇ∇ˇ ˆ g − ω [ˆ ~D∇ τ3 , gˇ] = 0 , gˇ2 = ˇ1 . (2) Here D is the diffusion constant, ω = (2n + 1) πT

(3)

is the Matsubara frequency, and τˆ3 is a Pauli matrix in the particle-hole (Nambu) space. The Green function may be parametrized as G F gˇ = (4) F † −G with a real G and complex conjugate F and F † (our definition of F and F † differs from that of Ref. 14 by a ˆ in Eq. (2) includes factor of i). The gradient operator ∇ the vector potential: ˆ = ∇ − ie A[ˆ ∇ τ3 , · ] . ~

(5)

We neglect the Zeeman splitting, which, in the case of the transverse magnetic field, has a typically much smaller effect than the vector-potential term. The Usadel equations (2) are complemented by the boundary conditions at the edges of the junction ∇y gˇ(x, y = ±Ly /2) = 0

(6)

Lx

Ly

S

H = He z x=0 y=0

N

S

ey ex

FIG. 1: SNS junction in a transverse magnetic field.

and by suitable boundary conditions at the SN interfaces (x = 0 and x = Lx ), depending on the transparency of the interface, the levels of disorder in the N and S regions, and the ratio of the Matsubara frequency to the superconducting gap.15 For simplicity, we assume that the superconducting gap is the largest energy scale in the problem and that the SN interfaces are perfectly transparent: then the boundary conditions take the simplest “rigid” form 0 e±iχ gˇ(x = {0, Lx }, y) = . (7) e∓iχ 0 The boundary conditions at the SN interfaces will not be important for our result on the transition between the purely decaying and oscillating-decaying regimes, but will play a role in estimating the pre-exponential factor in Eq. (27) below. There are several energy scales in the problem. We assume that the gap in the superconducting terminals is the largest energy scale. In particular, it is much larger than the energy scale associated with the magnetic length ξH (the “magnetic Thouless energy”) defined as EH =

DeH ~D = 2 ξH π

(8)

This energy scale is, in turn, assumed to be much higher than the temperature T . These assumptions agree with the high-field regime of the model considered in Ref. 7 and of the experiments in Refs. 8,9. The conventional Thouless energy ETh =

~D L2x

(9)

(which is smaller than EH ) plays minor role in our longjunction (or large-field) limit, since the decay of the proximity correlations is mainly determined by the magnetic length ξH . The current density can be calculated from F (x, y) as14 ∞ X 4ieA † † † J = 2πieνDT F ∇F − F ∇F − FF , ~ n=0 (10) where ν is the density of states in the normal metal phase (per one spin projection). The symmetry of translation

3 along the z direction implies that the current remains in the xy plane. The sum is taken over the Matsubara frequencies (3).

2

Im κ

1 III. DECAY OF SUPERCONDUCTING CORRELATIONS IN A LONG JUNCTION

Ly=0 For a sufficiently long junction (a specific condition will be determined below), the anomalous correlations F in the middle of the junction are exponentially weak. Therefore, to the leading order, we may approximate the solution to the SNS-junction problem as a superposition of two solutions to the semi-infinite SN problem16,17 F (x, y) ≈ F∞ (x, y)eiχ + F∞ (Lx −x, −y)e−iχ

(12)

The mode with the slowest decay may be written as x → ∞,

Re κ

2 2.6

−1

−2

where Jtot is the total Josephson current (integrated over the y and z directions). The problem now reduces to finding the mode of F∞ (x, y) with the slowest decay (along the x direction) in the auxiliary problem of a semi-infinite infinite SN junction. At sufficiently large x, the proximity effect is weak, and F∞ (x, y) obeys the linearized version of the Usadel equations (2) [obtained by assuming G ≈ sign(ω) and |F | 1]. For the sake of convenience, in all the equations below we rescale both coordinates x and y in the units of the magnetic length (1). Then the linearized Usadel equation takes the form7 2|ω| 2 (∇x + 2πiy) + ∇2y − F∞ (x, y) = 0 . (13) EH

F∞ (x, y) ∼ e−κx ψ(y) ,

1

(11)

where F∞ (x, y) is the off-diagonal matrix element of the (full nonlinear) solution gˇ∞ for the SN problem with the semi-infinite normal layer. It obeys the same equations (2) with the same boundary conditions at x = 0 and at y = ±Ly /2 and with the boundary at x = Lx replaced by gˇ∞ (x → ∞, y) = sign(ω)ˆ τ3 . This weak-coupling regime, via the expression (10), implies the sinusoidal currentphase relation18 Jtot = Ic sin 2χ ,

Ly=Lc

(14)

where κ is the inverse decay length. By substituting this expression into the linearized Usadel equation, we find that ψ(y) must be a zero mode of the linear operator 2|ω| 2 2 ∇y + (κ − 2πiy) − ψ(y) = 0 (15) EH with the boundary conditions ψ 0 (±Ly /2) = 0. The mode with the slowest decay may be found as the zero mode (15) with the smallest positive real part of κ. One can verify that the real part of κ given by the zeromode condition (15) increases with increasing |ω|. Therefore, under the low-temperature assumption T EH ,

−1 FIG. 2: Effective inverse decay length κ (in the units of ξH ). A purely real κ indicates a monotonic decay. The transition to a damped oscillating regime occurs at Ly = Lc ≈ 0.82 (in the units of ξH ) and corresponds to κ ≈ 2.6.

the main contribution to F∞ (x, y) at large x comes from Matsubara frequencies ω EH . For these frequencies, we may neglect the ω term in Eq. (15). The general solution to the second-order differential equation (15) at ω = 0 can be written in terms of a linear combination of two modified Bessel functions,19 ψ(y) =

p

(iκ + 2πy)2 4π ! (iκ + 2πy)2 + C2 K1/4 . (16) 4π

iκ + 2πy C1 I1/4

The boundary conditions at y = ±Ly /2 fix the ratio C1 /C2 and limit the possible values of κ to a discrete set. As a result, we obtain the values of κ (in the units −1 of ξH ) as a function of Ly (in the units of ξH ). The trajectory of κ in the complex plane, as the width of the junction Ly varies from zero to infinity, is plotted in Fig. 2. In the limit Ly → 0, the 2πiy term in the operator (15) may be neglected, and the spectrum is composed of the non-degenerate eigenvalues κ = nπ/Ly . The “leading” eigenvalue with the smallest real part is thus zero at Ly → 0. At a small finite Ly , the leading eigenvalue κ also becomes finite, but remains purely real. This follows from the combined symmetry of the complex conjugation and the reflection y 7→ −y, which relates the eigenvalues κ and κ∗ . Since the leading eigenvalue is nondegenerate in the limit Ly → 0, by continuity it must remain purely imaginary for sufficiently small Ly .20 At larger Ly , two real eigenvalues may collide and bifurcate to a complex-conjugate pair. This happens at Ly = Lc ≈ 0.82 (see Fig. 2). For Ly > Lc , the contributions of the two modes (corresponding to the complex

4 conjugate pair κ and κ∗ ) decay with the same rate (given by the real part of κ) and must be added together. In the discussion of the wide-junction limit (Section III B), we will show that those modes, in the limit Ly 1, correspond to solutions localized at the two edges of the junction y = ±Ly /2. The critical length Lc separates the regime where the superconducting anomalous Green function F (x, y) decays along the x direction without oscillations (narrow junction, purely real κ) and the regime where the decay of the Green function is damped oscillatory (wide junction, a pair of complex conjugate values of κ).

A.

1

|ψ(y)|/|ψ(Ly/2)|

0.5

Narrow-junction limit

For Ly 1 (in the units of ξH ) we expand the exact solution (16) in powers of Ly and find the wave number κ solving the equation for the boundary conditions at y = Ly /2. This yields the expansion 4π 2 4 π 932π 4 8 L + L + ... . κ = √ Ly 1 + (17) 63 y 218295 y 3 To the lowest order in Ly , the solution to the Usadel equation does not depend on y. In this limit, one can simply average the y 2 term in the Usadel equation (13) and arrive at a pair-breaking term7,8,21 ~D 2 ∇ F (x) = (|ω| + 2Γ) F (x) 2 x

(18)

with Γ=

De2 H 2 L2y ~Dκ2 = . 4 12~

(19)

This result obviously reproduces the first term in Eq. (17).

-Ly/2

y

Ly/2

FIG. 3: Superconducting pair correlations |ψ| normalized to their value at the border of the junction for Ly = 0.25ξH (dash), 0.75ξH (dot), ξH (dash dot) and 2.5ξH (solid).

This zero mode decays quasiclassically towards the opposite edge of the junction, and with an exponential precision we can replace the boundary condition at y = −Ly /2 by the decaying condition at infinity, ψ(y → −∞) = 0. This selects a solution from (16) of the form p (iκ + 2πy)2 ψ ∝ iκ + 2πy K1/4 . (21) 4π Imposing now the boundary condition ψ 0 (Ly /2) = 0 implies an equation on κres : (−iκres )2 = 0. (22) K3/4 4π A numerical solution to this equation gives22

B.

Wide-junction limit

κres ≈ 2.32 − 1.68i . In the limit of a wide junction Ly 1 (as usual, in the units of ξH ), the solution is determined by the two complex conjugate wave vectors κ and κ∗ . We show below −1 that the asymptotic behavior of κ (in the units of ξH ) is κ ≈ iπLy + κres .

(20)

Indeed, in the wide-junction limit, each of the two zero modes (16) is localized near one of the two edges of the junction and decays quasiclassically towards the other edge. The solution localized near y = Ly /2 should therefore have the quasiclassical wave vector in the operator (15) vanishing in that region, which immediately gives the leading asymptotics κ ≈ iπLy (the solution localized at the opposite edge has κ ≈ −iπLy ). To get the subleading term κres , we consider one of those zero modes (say, the one localized near y = Ly /2).

(23)

We illustrate our calculation in Fig. 3, where we plot the zero modes below and above the transition. Below the transition (for Ly < Lc ), the solution is nondegenerate and symmetric, while above the transition (Ly > Lc ) the two zero modes are pushed towards the edges of the junction. The characteristic size of the region near the edge where the proximity correlations are localized is of the order one [from the solution (21)], i.e., ξH in the physical units. IV.

JOSEPHSON CURRENT

In Section III, we have calculated the wave vectors κ of decay of the superconducting correlations in a long SN junction. As we shall see below, the same wave vector

5 κ describes the dependence of the Josephson current on the length Lx of a SNS junction. In particular, the Lx dependence of the Josephson current exhibits two types of behavior: purely decaying (for Ly < Lc ) and decaying with oscillations (for Ly > Lc ). This result is consistent with the previous numerical works,7 which identified these two regimes. To calculate the Josephson current, we substitute the asymptotic behavior (14) of F∞ (x, y) (possibly containing two terms in case of a complex κ), applicable in the middle of the junction, into Eq. (10). It produces the sinusoidal current-phase relation (12). In the “pure decay” phase (Ly < Lc , real κ), the critical current Ic is given by # "Z ∞ Ly /2 X Ic = 8πeνDT Lz (κ − 2πiy)ψ 2 (y)dy e−κLx . n=0

−Ly /2

(24) where Lz is the dimension of the junction along the z direction. Note that this expression is real [since ψ(y) = ψ ∗ (−y) in this regime] and positive (one checks this numerically). In the regime of “decaying oscillations” (Ly > Lc , complex κ), the anomalous Green function F∞ (x, y) contains contributions from two zero modes: ∗

F∞ (x, y) = ψ(y)e−κx + ψ ∗ (−y)e−κ

x

,

(25)

where ψ(y) 6= ψ ∗ (−y) are the two leading zero modes (15). Integrating the critical current along the y and z directions, one finds Ic = 8πeνDT Lz

∞ X n=0

Z

Ly /2

Re

(κ − 2πiy)ψ 2 (y)dy e−κLx .

−Ly /2

(26) Note that in the case of a wide junction, Ly ξH , the localization of the superconducting pair correlations at the edges of the junction results in the localization of the current in the same region. In Eqs. (24) and (26), both κ and ψ(y) depend on the index n of the Matsubara frequency (3). Since Re κ grows with ω, only terms with ω EH contribute in the largeLx limit. Therefore, the leading exponential dependence of Ic on Lx is given by Ic ∼ Re Ic(0) exp[−κLx ] ,

(27)

where κ is now taken at ω = 0 (plotted in Fig. 2). The above formula is equally applicable in the “pure decay” (0) phase (with real κ and Ic ) and in the “decaying oscil(0) lations” phase (with complex κ and Ic ). Re κ and Im κ describe the rates of decay and oscillations of the critical current as a function of Lx , respectively. (0) An accurate calculation of Ic would require knowing the normalization of ψ(y) (which depends on the specific boundary conditions at the NS interface and may require solving the nonlinear Usadel equations near the interface) and the ω dependence of κ. We do not perform such a

0.82 Lx/ξH

Exponential decay Oscillations + exponential decay

Ic

φ/φο =1

φ/φο Effective pair breaking

~1 Fraunhofer

0

Ly/ξH

FIG. 4: The phase diagram of the junction in a magnetic field. In the left region (Ly /ξH < 0.82), the critical current monotonically decays as a function of the field. In the right region (Ly /ξH > 0.82), the current oscillates and decays as a function of the field.

detailed calculation, but only give order-of-magnitude estimates, under the assumption of a perfectly transparent SN interfaces with rigid boundary conditions (7) for the Usadel equations (see Appendix B for details). In the low-temperature regime (below a certain temperature scale T∗ ), many Matsubara frequencies contribute to the total current (24) or (26), while at temperatures above T∗ only the lowest Matsubara frequency is relevant. A simple estimate (see Appendix B) gives ( Ly φ Ly ξH , φ0 ETh = Lx EH , T∗ ∼ √ (28) Lx EH ETh = ξH ETh , Ly ξH , where ETh is the Thouless energy defined in Eq. (9) and φ = HLx Ly is the total flux through the junction. Note that in these estimates and in the estimates below we omit numerical coefficients which may reach an order of magnitude.16,23 For narrow junctions, Ly ξH , we find (see Appendix B) 2 φ I , T T∗ , 0 φ0 (29) Ic(0) ∼ φ I T , T T , 0 ETh

φ0

∗

where I0 is the critical current through the same junction at low temperature (T ETh ) without magnetic field. For wide junctions, Ly ξH , the critical current may be estimated (see Appendix B) as q 3/2 φ Lx I , T T∗ , 0 φ0 Ly Ic(0) ∼ (30) Lx I T , T T . 0 ETh

Ly

∗

We sketch the phase diagram of the junction in Fig. 4 in the coordinates Lx and Ly .

6 In the region Ly < Lc , there is no sign reversal of the Josephson coupling (no π phase). The condition of the applicability of our result (27) is κLx 1, which translates in this region to φ/φ0 1. This is exactly the condition of a “sufficiently long” junction which justifies the weak-coupling approximation (11). At weaker fields (or for shorter junctions), φ/φ0 ≤ 1, the expression (27) is not accurate, but the pair-breaking approximation (18)– (19) remains valid, provided Ly ξH . In the region Ly > Lc , there are sign reversals of the Josephson coupling, and the oscillating exponential decay (27) is applicable as long as Lx ξH (since Re κ is of the order one in this phase). For shorter junctions, Lx ξH , the exponential formula (27) crosses over to the conventional Fraunhofer interference pattern5,6 (see Appendix C for details). In experiments, one usually varies the external field for a junction of fixed dimensions. This corresponds to scanning the phase diagram in Fig. 4 along a diagonal with a fixed ratio Lx /Ly . Then, depending on this ratio, one may observe a crossover from the pure-decay (for Ly Lx ) or Fraunhofer (for Ly Lx ) regime to the exponential oscillating regime as the field increases. For junctions with Ly Lx , the low-field regime (ξH Ly ) corresponds to the “effective spin flip” and “exponential-decay” regions of the phase diagram in Fig. 4. There, the field dependence of the critical current is given [by combining Eqs. (27) and (29)] by α1 φ π φ Ic ∝ exp − √ , (31) φ0 3 φ0 where α1 = 2 for T T∗ and α1 = 1 for T T∗ (note that T∗ itself depends on the magnetic field). The above expression is only applicable in the φ φ0 limit. Note that our T T∗ result agrees with the φ φ0 limit of the more general expression Ic ∝

π φ √ 3 φ0

sinh

π φ √ 3 φ0

(32)

derived in Ref. 24 under the assumption of the sinusoidal current-phase relation and with the result of Ref. 25 in the context of a spin-flip scattering. Upon increasing the magnetic field, at ξH Ly , the same junction crosses over to the “oscillations and decay” region of the phase diagram in Fig. 4. There, the critical current exhibits the decaying-oscillating behavior [found by combining Eqs. (27) and (30)]: α2 φ Lx πφ Lx Ic ∝ exp −2.32 sin − 1.68 + ϕ0 , φ0 ξH φ0 ξH (33) where α2 = 1/2 for T T∗ and α2 = 0 for T T∗ . Here ϕ0 is a phase shift, which we do not compute in this work. The asymptotic expression (33) is applicable at Lx , Ly ξH , which implies, in particular, φ φ0 . Note that while both expressions (31) and (33) decay

exponentially with increasing the field, the expression in the exponent of Eq. (31) is proportional √ to H, while that in the exponent of Eq. (33) only to H. For junctions with Ly Lx , the low-field regime (ξH Lx ) belongs to the Fraunhofer region of the phase diagram in Fig. 4. In this region, a simple analysis of our model produces the usual Fraunhofer pattern of the current-field dependence5,6 with 1/φ decay and oscillations with the period φ0 (see Appendix C). In the same junction, with the increase of the field (ξH Lx ), this Fraunhofer pattern crosses over to the oscillating-decay regime (33). In this crossover, the period of oscillations (as a function of the field) remains the same, only the 1/φ rate of decay gets replaced by an exponential.

V.

SUMMARY AND DISCUSSION

To summarize, we consider a long diffusive SNS junction in an external magnetic field H. We show that depending on the width p of the junction relative to the magnetic length ξH = φ0 /H two different regimes can be observed. For narrow junctions the anomalous Green function F decays monotonically along the junction while for wide junctions exponentially damped oscillations are present. We find that the transition between the two regimes occurs at the width Lc ≈ 0.82ξH . These different behaviors of the proximity correlations translate into two types of the dependence of the Josephson critical current on the magnetic flux through the junction: a monotonic decay for narrow junctions and an oscillating decay for wide junctions. We also show that for wide junctions, both the superconducting pair correlations and the current are concentrated in a small region of size ξH near the edges of the junction. The main finding of the present work, in comparison with previous studies of this problem, is the identification of the damped-oscillating phase for wide and long junctions. This phase may be understood as a “crossover” between the Fraunhofer and the purely damped regimes: the period of oscillations is the same as in the Fraunhofer interference pattern, while the exponentially decaying factor resembles the damped phase. Conceptually, the transition between the two asymptotic regimes for long junctions in our problem is similar to the transition between the two regimes in superconductor–ferromagnet–superconductor junctions with domains studied in Ref. 26. In both systems, the transition between the purely damped and dampedoscillating behavior is related to a bifurcation of the solution to the linearized Usadel equations. Diffusive SNS junctions in a magnetic field have been the subject of recent experimental studies,8,9 and both the Fraunhofer regime and the purely decaying have been observed. The oscillating-exponential regime predicted in our present paper has not been yet identified in experiments. We propose that it may be most conveniently observed in nearly square junctions with Lx & Ly . In such

7 junctions, the purely decaying regime at low fields crosses over to the oscillating-exponential regime at higher (but not very high) fields, see Fig. 4. A robust qualitative feature of this crossover (following from our analysis and observed in the numerical studies of Ref. 7) is the flux of the first reversal of the critical current being larger than the subsequent period of current reversals [which tends to φ0 at large fields, see Eq. (33)]. Such an experiment may, for example, be realized using SNS junctions similar to those studied in Ref. 9, but with intermediate ratios Lx /Ly . If one considers a junction similar to the WAu-Sq sample of Ref. 9, but with Lx = 1.2 µm and Ly = 1.0 µm, then the crossover between the purely decaying and the oscillating-decaying phases occurs around ξH ∼ 1.2 µm (the bifurcation point in Fig. 2), which corresponds to the magnetic field H ∼ 14 G. Using the Thouless energy ETh = 5.3 µeV reported in Ref. 9 for their WAu-Sq sample, we estimate the magnetic energy scale (8) in our example at this magnetic field to be EH ∼ ETh ∼ 60 mK. Therefore, for such a junction, at T ∼ 60 mK (the temperature used in Ref. 9), we expect that our results (assuming T EH ) would hold qualitatively, but not quantitatively. A quantitative agreement would improve at higher fields (deeper in the damped-oscillating region) or at lower temperatures. Acknowledgments

This work was supported by the Swiss National Foundation. We thank P. Ostrovsky, M. Skvortsov, and S. Tollis for helpful discussions.

This approximation has been frequently used in earlier studies of SNS and SN systems in a magnetic field.7,10,12 It amounts to neglecting the y component of the vector −1 potential (A1), which is justified as long as Ay φ0 ξH . This condition may be rewritten as H φ0 /λ2 ∼ Hc1 (the lower critical field). This estimate can also be applied to type I superconductors, where it suffices to assume that the field is much lower than the critical field. Note that the finite value of λ may be taken into account in the renormalization of the x component of the vector potential (A1). In other words, the vector potential (A2) used in the derivations in the main body of the paper involves the effective field Heff = H(Lx + 2λ)/Lx instead of the actual field H. As a result, the flux φ in all our formulas should be understood as the total flux through the junction, including the field penetrating in the superconducting leads. Appendix B: Estimating the magnitude of the critical current (0)

The estimate of the critical current Ic in Eq. (27) depends on whether the temperature T is high or low (either one or many Matsubara frequencies contribute) and on the boundary conditions at the SN interface. In the large-Lx limit assumed in our derivation, the main dependence on the Matsubara frequency comes from the exponent exp[−κLx ] in Eqs. (24) and (26). Therefore the temperature scale T∗ separating the highand low-temperature limits may be estimated from the condition T∗ Re

Appendix A: Gauge choice and neglecting the magnetic penetration length

We fix the gauge of the magnetic field in such a way that the phase of the superconducting order parameter is constant in space in each of the two leads. With this gauge choice, the vector potential in the superconducting leads is proportional to the Meissner screening current. Assuming that the field is uniform in the normal part of the junction and exponentially decaying on the penetration length λ inside the superconductors (Meissner effect), we choose the gauge as follows: x/λ x < 0 , λHe ey , λH 1 − 2x e − yH 1 + 2λ e , y x Lx Lx A= (A1) 0 < x < Lx , −λHe−(x−Lx )/λ e , x > L . y

x

In the main body of our paper, we take the limit λ → 0: then the Meissner screening current is also weak, and the vector potential (A1) reduces to ( −yHex , 0 < x < Lx , A= (A2) 0, x < 0 or x > Lx .

∂κ Lx ∼ 1 . ∂ω

(B1)

In turn, ∂κ/∂ω may be estimated by a variation of Eq. (15): R Ly /2 2 ψ (y)dy ∂κ 1 −Ly /2 = . ∼ R Ly /2 2 ∂ω E hκ − 2πiyi H EH (κ − 2πiy)ψ (y)dy −Ly /2

(B2) The average hκ−2πiyi is taken over the region of y where the mode ψ(y) is located. In the narrow-junction regime Ly ξH , ψ(y) is spread over the whole junction (Fig. 3), and hκ − 2πiyi ∼ κ ∼ Ly [see Eq. (17)]. In the widejunction regime Ly ξH , ψ(y) is localized near the edge of the junction (Fig. 3), and hκ − 2πiyi ∼ κres ∼ 1 [see Eq. (20)]. Combining all these considerations, we arrive at the estimate (28) for T∗ . At T T∗ , only the lowest Matsubara frequency is relevant, and the currents (24) and (26) may be estimated as Z Ly /2 Ic(0) ∼ eνDLz T (κ − 2πiy)ψ 2 (y)dy . (B3) −Ly /2

Assuming the rigid boundary conditions (7), we estimate the overall amplitude of ψ(y) to be of order one, and

8 therefore ( Z Ly /2 (Ly /ξH )2 , 2 (κ − 2πiy)ψ (y)dy ∼ 1, −Ly /2

Ly ξH , Ly ξH . (B4) This produces the high-temperature estimates in Eqs. (0) (29) and (30). For convenience, we express Ic in terms of the critical current I0 ∼ eνDLz ETh

Ly Lx

difference linearly depending on y. The gauge-invariant phase difference across the junction is Z 2e 2e Ax dx = 2χ + yHLx (C1) ϕ(y) = 2χ − ~ ~ (so that 2χ equals the gauge-invariant phase difference in the middle of the junction y = 0). The critical current is then

(B5)

at low temperature (T ETh ) through the same junction at zero magnetic field.16,23 At T T∗ , one needs to sum over many Matsubara frequencies, up to the frequencies of the order T∗ . As a consequence, the same estimate (B3) applies, but with T replaced by T∗ . This immediately produces the lowtemperature estimates in Eqs. (29) and (30).

where j(ϕ) is the Josephson current density without magnetic field at the phase difference ϕ across the junction. At temperatures T ETh , the current-phase relation is sinusoidal,16,23 and the Fraunhofer pattern (C2) takes the simplest form,

Appendix C: Fraunhofer interference pattern in wide and short junctions

For completeness of our discussion, we remark that our model of Section II produces Fraunhofer-like interference pattern5,6 in the limit Lx ξH and Ly ξH . In such wide and short junctions, the gradients along the y direction in the Usadel equation (2) may be neglected, and the resulting system is equivalent to an ensemble of Josephson junctions connected in parallel, with the phase

1 2

3

4

5 6

7

8

9

A. I. Buzdin, Rev. Mod. Phys. 77, 935 (2005). F. S. Bergeret, A. F. Volkov and K. B. Efetov, Rev. Mod. Phys. 77, 1321 (2005). L. N. Bulaevskii, V. V. Kuzii, and A. A. Sobyanin, Pisma Zh. Eksp. Teor. Fiz. 25, 314 (1977) [JETP Lett. 25, 290 (1977)]. V. V. Ryazanov, V. A. Oboznov, A. Yu. Rusanov, A. V. Veretennikov, A. A. Golubov, and J. Aarts, Phys. Rev. Lett. 86, 2427 (2001); V. V. Ryazanov, V. A. Oboznov, A. S. Prokofiev, V. V. Bolginov, A. K. Feofanov, J. Low Temp. Phys. 136, 385 (2004); V. A. Oboznov, V. V. Bolginov, A. K. Feofanov, V. V. Ryazanov, and A. I. Buzdin, Phys. Rev. Lett. 96, 197003 (2006). B. D. Josephson, Rev. Mod. Phys. 36, 216 (1964). A. Barone and G. Paterno, Physics and Applications of the Josephson Effect (Wiley, New York, 1982). J. C. Cuevas and F. S. Bergeret, Phys. Rev. Lett. 99, 217002 (2007); F. S. Bergeret and J. C. Cuevas, J. Low Temp. Phys. 153, 304 (2008). L. Angers, F. Chiodi, G. Montambaux, M. Ferrier, S. Gu´eron, H. Bouchiat, and J. C. Cuevas, Phys. Rev. B 77, 165408 (2008). F. Chiodi, M. Ferrier, S. Gu´eron, J. C. Cuevas, G. Montambaux, F. Fortuna, A. Kasumov, and H. Bouchiat, Phys. Rev. B 86, 064510 (2012).

Ly /2

2e j 2χ + yHLx dy χ ~ −Ly /2 Z πφ/φ0 Ly Lz φ0 j(2χ + φ0 ) dφ0 , (C2) = max χ 2πφ −πφ/φ0 Z

Ic = Lz max

Ic ∝

sin(πφ/φ0 ) . πφ/φ0

(C3)

At lower temperature, the current-phase relation j(ϕ) is non-sinusoidal, and the sine function in Eq. (C3) is replaced by a different, albeit qualitatively similar, periodic dependence with the same period. We remark that this periodicity is the same as in the decaying-oscillating regime (33).

10

11

12

13 14

15

16

17

18

19

V. P. Galaiko, Zh. Eksp. Teor. Fiz. 57, 941 (1969) [Sov. Phys. JETP 30, 514 (1970)]; G. A. Gogadze and I. O. Kulik, Zh. Eksp. Teor. Fiz. 60, 1819 (1971) [Sov. Phys. JETP 33, 984 (1971)]; T. N. Antsygina, E. N. Bratus, and A. B. Svidzinskii, Fiz. Nizk. Temp. 1, 49 (1975) [Sov. J. Low Temp. Phys. 1, 23 (1975)]. G. Mohammadkhani, M. Zareyan, and Y. M. Blanter, Phys. Rev. B 77, 014520 (2008). W. Belzig, C. Bruder, and G. Sch¨ on, Phys. Rev. B 53, 5727 (1996). K. D. Usadel, Phys. Rev. Lett. 25, 507 (1970). N. Kopnin, Theory of nonequilibrium superconductivity (Oxford University Press, 2001). A. Zaitsev, Zh. Eksp. Teor. Fiz. 86, 1742 (1984) [Sov. Phys. JETP 59, 1015 (1984)]; M. Y. Kupriyanov and V. F. Lukichev, Zh. Eksp. Teor. Fiz. 94, 139 (1988) [Sov. Phys. JETP 67, 1163 (1988)]; A. Altland, B. Simons, and D. Taras-Semchuk, Adv. Phys. 49, 321 (2000). A. D. Zaikin and G. F. Zharkov, Fiz. Nizk. Temp. 7, 375 (1981) [Sov. J. Low Temp. Phys. 7, 184 (1981)]. D. A. Ivanov, R. von Roten, and G. Blatter, Phys. Rev. B 66, 052507 (2002). A. A. Golubov, M. Y. Kupriyanov, and E. Il’ichev, Rev. Mod. Phys. 76, 411 (2004). M. Abramovitz and I. A. Stegun, Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical

9

20 21

22

Tables (Dover, New York, 1972). We thank M. Skvortsov for pointing to us this symmetry. B. Crouzy, E. Bascones, and D. A. Ivanov, Phys. Rev. B 72, 092501 (2005). The Macdonald function K3/4 (z) does not have zeros on the principal sheet of the Riemann surface, therefore the solution must have π/2 < | arg(−iκres )| < π. We also impose the condition Re κres > 0, since the mode decays along the x direction. Finally, under these constraints, we select the root with the smallest real part. Technically, the second sheet of the Riemann surface may be accessed with

23

24 25

26

the relation K3/4 (eiπ z) = e−3πi/4 K3/4 (z) − iπI3/4 (z). P. Dubos, H. Courtois, B. Pannetier, F. K. Wilhelm, A. D. Zaikin, and G. Sch¨ on, Phys. Rev. B 63, 064502 (2001). G. Montambaux, arXiv:0707:0411. J. C. Hammer, J. C. Cuevas, F. S. Bergeret, and W. Belzig, Phys. Rev. B 76, 064514 (2007). B. Crouzy, S. Tollis, and D. A. Ivanov, Phys. Rev. B 76, 134502 (2007).