Square-well model for superconducting pair-potential

We study Andreev reflection in a one-dimensional square-well pair-potential. We discuss the history of the model. The current-phase relation is presented as a sum over Matsubara frequencies. How the current arises from bound and continuum levels is found by analytic continuation. We discuss two limiting cases of the square-well potential, the zero-length well and the infinite well. The model is quantitatively valid in some cases but forms the basis for understanding a wide range of problems in inhomogeneous superconductivity and superfluidity.

Andreev reflection is a fundamental process in superconductivity.The process consists of a particle type excitation being reflected as hole type excitation.In contrast to usual reflection, the momentum change in Andreev reflection is small compared to the particle's momentum, as the hole type excitation moves in the opposite direction than its momentum.Andreev reflection treats the transport of electricity and heat differently.It is the essential process that allows electric current to flow in any inhomogeneous superconductor but it resists the heat current.
The purpose of this article is to study Andreev reflection in a simple model.We model the pair potential by a square-well form, see Fig. 1.This could be compared with the model of square-well potential and its special case the infinite potential well, that are studied in practically all textbooks of quantum mechanics.The study of such models has a twofold justification. 1) Although the model is not precisely realized in Nature, it still forms a reasonable approximation for several applications.2) A model that can be solved precisely allows one to learn physics, which is useful also in cases where the model ceases to be correct quantitatively.1 arXiv:2405.07659v2[cond-mat.supr-con]29 Jun 2024 Fig. 1 The square-well potential for the pair potential ∆(x): it equals ∆ 0 e −iϕ/2 for x < −L/2, 0 for |x| < L/2, and ∆ 0 e iϕ/2 for x > L/2.The figure also shows bound state energies (dashed lines).These are evaluated with parameter values L∆/ℏv F = 4.5 and ϕ = 1 rad.At energies E > ∆ 0 the states form a continuum.
We specify the parameters of the square-well model (Fig. 1).We consider the pair potential ∆(x) that depends on a single coordinate x.The potential is assumed to vanish in the well of width L. Outside of this range, the amplitude is constant, |∆| = ∆ 0 .We allow a phase difference ϕ between the two sides.Because the absolute phase is not important, it is convenient to choose the phase symmetrically, ∆(x) = ∆ 0 e −iϕ/2 for x < −L/2, and ∆ = ∆ 0 e iϕ/2 for x > L/2.It is convenient to call the middle region as normal (N), and the surrounding regions as superconducting banks (S).We can think that x is a coordinate along a quasiparticle trajectory in a bulk metal.Alternatively, we can consider the system to represent a single quantum channel.(For an introduction to quantum channels see Ref. [1].)We neglect all scattering (by interfaces, walls, impurities or lattice vibrations) except the one caused by the pair potential.We count the energies starting from the Fermi level.For brevity, the bank gap amplitude ∆ 0 is denoted by ∆ in the following.
The motivation for Andreev was to understand experimental observations that the heat conductivity in the intermediate state was found much smaller than in the normal or zero-field superconducting state.The intermediate state of type I superconductors, which occurs in the magnetic field, was known to have a lamellar structure consisting of alternating normal and superconducting regions.This naturally leads to thinking that the normal-superconducting interfaces could be responsible for the reduced heat conductivity.The problem, however, is that the force on the electron at such an interface is too weak to cause ordinary scattering of electrons that have momentum near the Fermi momentum.Andreev solved this puzzle by inventing a new type of reflection, where an incident electron is reflected as a hole-type excitation.The latter travels in the opposite direction at a momentum that is only slightly changed from that of the incident electron.In the first publication on the topic [2], he used Gor'kov's equations to solve the reflection problem for a pair-potential step.More generally, he showed that the result for a realistic N-S interface is not essentially different for electron energies slightly above the gap, the relevant ones for heat conduction at a low temperature.This allowed him to estimate the heat conductivity of the electrons at low temperature, in good agreement with experiment.
Andreev studied the problem of the square-well potential (Fig. 1) in Ref. [3].At energies below ∆, an electron is reflected from one step as a hole, which is then reflected back from the other step as an electron.This leads to the formation of discrete energy levels, which we now call Andreev bound states.For low-energy levels, he found equally spaced energy values Here v F is the Fermi velocity, n is an integer and γ is a constant, and it is assumed L ≫ ℏv F /∆.Note that the level spacing is the same as one would obtain in a deep ordinary potential well at energies close to the Fermi energy.Andreev applied this to the intermediate state, where the discreteness of the levels should be visible in thermodynamic properties at temperatures T ≲ ℏv F /L.The next important development was by Kulik [4].He had the idea that the square-well pair-potential could be used as a model of the SNS Josephson junction.For that, he took into account that there could be a phase difference between the superconducting banks.For the bound state energies, he found the equation This formula is valid for all bound state energies and arbitrary L. For low-energy levles, the last term in the parenthesis approaches π.Thus the energy (2) reduces to the Andreev's result (1) but now gives a phase-dependent γ(ϕ) = 1 2 ± ϕ/2π.The levels that increase in energy with increasing ϕ correspond to momentum to the right.This means that there is an electron moving to the right and a hole to the left.For levels that decrease in energy with increasing ϕ, it is the opposite.Thus, depending on the phase difference and the occupations of these levels, there is a current through the system.Thus the SNS system forms a Josephson junction.
Kulik went on to calculate the current.His method was complicated.In particular, the contributions from the bound states and continuum were calculated separately.The result reported in Ref. [4] turned out to be incorrect.In particular, the continuum contribution was found to have a phase dependence of sin ϕ, contrary to later calculations to be discussed below.
The Josephson current was next calculated by Ishii [5].His method was also rather complicated.An essential advantage in his calculation was that he used the imaginary energies ϵ m = πT (2m + 1) (also known as Matsubara frequencies or energies), where m is an integer.This allowed him to write a single compact expression that included both bound state and continuum contributions.His general expression [Eq.(3.5) in [5]] can be written in the form The same result has been obtained in several publications using different approaches [6][7][8][9][10][11][12][13].We comment on some of the calculations.Svidzinskii, Antsygina, and Bratus' solved the problem both by using the Gorkov equations and by using the Bololiubovde Gennes equations, and discussed the problem in the earlier work by Kulik [6,7].
Bezuglyi, Kulik and Mitsai already went on to study more general problems, but commented on the earlier calculations and provided alternative expressions for the current in limiting cases [8].Kulik and Omel'yanchuk used the quasiclassical Green's function satisfying the Eilenberger equation to study the special case L = 0 [14,15].Kupriyanov used the same method in the case of general L and derived the current in the same form as given in (3) [9].He also discussed several limiting cases.A calculation aimed for conducting channels and based on Bololiubov-de Gennes equations is given by Furusaki, Takayanagi and Tsukada [10].A same type of derivation with plenty of explanations was given by Bagwell [11].The derivation using the Eilenberger equation is presented in the review [12] together with numerous studies of the Josephson current in various models.Recently the derivation of (3) has been presented once more in response to false claims about its invalidity [13].In summary, formula (3) was discovered over 50 years ago, and it has formed the starting point for numerous other studies since.Equation ( 3) is a highly convenient expression to calculate the current.The problem in (3) is its rather abstract form as a sum over imaginary energies.In order to gain physical understanding, one should do an analytic continuation of the imaginary energies E = iϵ m to the real E axis.We notice that the denominator of (3) vanishes on the real axis at points This is just another expression for the bound state energies given in (2).We express the sum over ϵ m as an integral and deform the integration path as explained in many textbooks, for example in Ref. [16], and also in Ref. [5].The result is The first term is the contribution from the bound states at |E| < ∆.It is the sum over the energy values satisfying (4).The second term is the contribution from the branch cut of the square root in the continuum |E| > ∆.In both terms a multiplier is the Fermi function f (E) = 1/(e E/T + 1).The derivation given here proves formula (5) for equilibrium occupations of the levels, but the equation is valid also for distributions f (E) differing from 1/(e E/T + 1).A noteworthy feature of the expression ( 5) is that it is a sum and integral over both positive and negative energies.This is a consequence of the Green's function technique, where a hole-type excitation (a particle is removed from the ground state) appears as a negative-energy pole and a particle-type excitation (a particle added to the ground state) as a positive-energy pole.We call the presentation of both negative and positive energies the semiconductor picture, since it makes the superconducting gap look like the gap between the valence and conduction bands in a semiconductor.An alternative to the semiconductor picture is the excitation picture, which shows only the excitation energies, which are positive for both types of excitations.This is used in Fig. 1.In the following, we prefer the excitation picture for two reasons.Firstly, consider a level of a given wave vector that decreases in energy from positive to negative.In the semiconductor picture this goes smoothly, in contrast to the excitation picture, where an excitation of the given momentum disappears and an excitation of the opposite momentum is created.Secondly, the equilibrium current (the supercurrent) can be interpreted to arise mainly from filled negative energy levels since f (E) in ( 5) is small for the positive energy levels (and vanishes at T = 0).The semiconductor picture is a direct generalization of the simple model of the normal state, where the free electron levels are filled up to the Fermi energy.
We study expression (5) for the current by evaluating the different terms numerically at zero temperature, T = 0.The results are presented in Figures 2-4.The difference between the three figures is the length of the normal region.The left panels in the figures show the energies of the bound states as functions of the phase difference ϕ.These are conveniently obtained by plotting ϕ as a function of E (4) and then exchanging the axes.Because of the symmetry of the energy levels with respect to E = 0, we discuss only the levels at negative energies, which are filled in equilibrium at T = 0.With increasing ϕ the levels degenerate at ϕ = 0 shift up or down in energy depending on their direction of momentum, as illustrated in Ref. [13].In the case of Fig. 2, there are two bound states for small ϕ.At a critical ϕ c = 2 the negativemomentum level merges with the continuum, so that only the positive-momentum bound state survives in the range ϕ c < ϕ < π.
The contribution of a bound state to the current ( 5) is given by the derivative of its energy with respect to the phase.These are evaluated in the middle panels of Figs.2-4.In the case of Fig. 2, we see that the currents of the two filled bound states are in opposite directions.Their currents cancel at ϕ = 0.With increasing ϕ the negative contribution decreases leaving a net positive current (in units of 2e∆/ℏ, where e is the electron charge).At ϕ > ϕ c the negative-momentum level has moved to the continuum so that the bound-state current comes from the single positive-momentum level.
The right panels of Figs.2-4 show the continuum current (5) by orange lines.We see that it vanishes only at special points ϕ = nπ with integer n.In the case of Fig. 2, it has a maximal negative value at the cusp at ϕ = ϕ c , where the negative-momentum bound state disappears.The total current is given by the blue line.This equals the sum of the bound-state currents (middle panel) and the continuum contribution (orange line).We see that the cusp in the continuum current just cancels the vanishing of bound-state current of the negative-momentum level at ϕ c .The total current does not exhibit any special behavior at ϕ c .Similar conclusions are obtained in all Figs.2-4.
The insensitivity of the current on the bound states near the gap edge was related to scattering theory in Refs.[6,7].It can also be argued that the current is determined by levels near the Fermi level and thus should be independent of the details at the relatively high energy of the gap edge [13].Further plots of the different contributions to the current can be found in Refs.[11,17].We see from the plots above that the zero-temperature current is discontinuous at ϕ = π.The current in the neighborhood of this point is This current comes solely from the two levels that cross the Fermi level at ϕ = π, as the contributions from all other bound states cancel each other and the continuum contribution vanishes.The quantity L * can be interpreted as the effective length of a bound state, where the latter term in L * describes the penetration of the bound state in the superconducting banks [18].The square-well pair-potential model has a few limiting cases.The case of the zerolength potential well, L = 0, was first studied by Kulik and Omel'yanchuk [14,15].It is easy to calculate from equations ( 4) and ( 5) the bound state energies E = ±∆ cos(ϕ/2) and the current The current comes from the two bound states only as the continuum contribution vanishes.
The second limiting case we consider is the infinite potential well.For an ordinary potential, this model is well presented in books of quantum mechanics because it has a fully analytic solution.In the case of the pair potential we take the limit ∆ → ∞.It is easy to calculate from equations ( 4) and (3) that the bound states have energies The infinitely deep potential well is a good approximation for low-energy levels, E n ≪ ∆ in a long junction, L ≫ ℏv F /∆ because in this limit the bound state energies (2) become independent of ∆.This model is equivalent to the one used by Bardeen and Johnson [18] (except that they use L * instead of L).
We study the current (8).At T = 0 the sum goes to an integral.Evaluating this gives Infinitely-deep pair-potential well (8) evaluated at temperatures 0, 1,. . ., 5 in units of ℏv F /2πL.In both panels, an increasing temperature corresponds to lower current, but the scales of the current and the temperature differ essentially.
and repeated periodically with period 2π.This result was found by Ishii (the linear dependence on ϕ in Ref. [5] and the coefficient in Ref. [19]).We see that the linear dependence is rather well satisfied in the selected examples of the general model, Figs.2-4.
An interesting feature of the current ( 9) is that it corresponds to motion at a velocity v s = (ℏ/2m)ϕ/L, which can be interpreted as the superfluid velocity.Thus the normal region seems not to differ from a superconductor, as long as the critical velocity v c = πℏ/2mL is not exceeded.Expressed in other words, the system behaves like it would be Galilean invariant, although it is not because the locations of the N-S interfaces are fixed [18].
The temperature dependence of the current can be calculated from (3) in the general case.Examples in two limiting cases are shown in Figure 5.
Let us consider the application of the square-well pair-potential model to experimentally realizable systems.The main deficiency of the model in several applications is that the pair-potential is not self-consistent.For example, consider a normal section in a uniform-thickness superconducting wire.Because ordinary scattering was neglected, the normal region induces a proximity effect on the superconducting banks.As a result, the pair potential is reduced in the banks near the normal region.This affects especially the levels near the gap edge, which can penetrate deeper into the banks.There are, however, some cases where the square-well model can be justified quantitatively.
• The limit of a long junction, L ≫ ℏv F /∆.Here the current is determined by lowenergy levels, E n ≪ ∆.The energies (2) of these levels are not dependent on ∆.Thus we expect this to be the case also for ∆(x) modified by the proximity effect.• The case of a weak link separating two wider superconducting leads.In this case the order parameter in the banks is determined by quasiparticle trajectories staying in the leads with a minimal contribution from the transmitting trajectories [14].• The case that the Fermi velocity in the banks is much larger than in the normal region.This means that only a few of the trajectories incident on the junction will be transmitted, and thus the proximity effect on the banks is small [9].
Another possible point of worry in the square-well pair-potential model could be that current conservation is strictly valid only in the normal region.There is no current deep in the superconducting regions of constant phase.This problem can be corrected by adding a spatially varying phase in the superconducting regions.In the three cases where the proximity effect in the banks can be neglected, the current in the normal region is always small compared to the critical current in the superconducting regions.Thus the induced phase gradient in the superconducting regions is small and does not affect the calculation of the bound states and the current.This means that current conservation is not a problem that is separate from the self-consistency problem.Thus attempts to restore current conservation without solving the self-consistency of the pair potential are misguided.
The square-well pair-potential model can be used as a qualitative model in cases where its assumptions are not strictly valid.For example, consider a flux line in a type II superconductor and a quasiparticle trajectory passing the vortex line at distance b (impact parameter).There is phase change along the trajectory.We can model this using the square-well model.The phase along the trajectory changes by 2π as a function of b changing from −∞ to +∞.There has to be a bound state whose energy changes from positive to negative energy, and an opposite momentum state that shifts from negative to positive energy.The zero crossing takes place as the trajectory goes through the vortex line, b = 0.This gives a qualitative semiclassical description of the Caroli-de Gennes-Matricon vortex-core states.
The model can also be used to demonstrate the difference between an ideal normal metal and a superconductor.Consider an occupation of the levels that the positive momentum levels are filled up to an energy E > 0 and correspondingly, there is a deficit in filling the negative momentum levels.Such a state would not decay if we strictly assumed the pure limit.The current in such a state can be considerably higher than that obtained with equilibrium occupation.This is not supercurrent since, for a filled positive momentum level, there are empty negative momentum levels at lower energy.Inelastic scattering can reduce the distribution and the current to their equilibrium values.
There are a great number of generalizations of the square-well pair-potential model.The most obvious is to allow for a self-consistent pair potential, see Refs.[20,21] for some early calculations.A second important generalization is to consider the normal scattering of electrons.This couples the levels with opposite momentum and results in that the crossings of the opposite momentum levels at ϕ = nπ (integer n) are replaced by avoided crossings [8,11,[22][23][24].This is illustrated in Fig. 6.There are other generalizations to non-equilibrium phenomena, different geometries, unconventional superconductors/superfluids etc.Some of the topics are considered in Refs.[12,25,26].In fact, all studies of inhomogeneous superconductors can be seen as generalizations of the square-well pair-potential model.
I will finish this article with a personal view of how I learned and used the squarewell pair-potential model.First, an impurity in the vortex core acts similarly as in a SNS junction at phase difference ϕ = π: it increases the superconductor's condensation energy.This leads to a pinning force on the vortex to the impurity [27,28].This is a single-impurity version of the effect known earlier that impurity scattering makes the  [24].The parameters are the transmission probability 0.9 at both barriers, L∆/ℏv F = 5.7, zero temperature, and the phase acquired in one round trip between the two barriers equal to π/2.[Note the misprint of an extra factor sin χ in formula ( 16) of Ref. [24].] superconducting coherence length shorter.Second, an impurity or a wall in a p-wave superfluid acts like a SNS junction as the phase on the quasiparticle trajectory jumps in scattering [29][30][31].Third, we have studied both equilibrium and non-equilibrium currents in pinholes, the analog of point contacts in superconductors, of superfluid 3 He [32][33][34].This allows one to explain π states observed experimentally in an array of pinholes.In summary, similarly to the role of the square-well potential model to learn quantum mechanics, the square-well pair-potential model is an analytically solvable model that allows us to learn inhomogeneous superconductivity and fermion superfluidity.

Fig. 2 Fig. 3
Fig.2Numerical results in the square-well pair-potential model with a short N region.(left) The bound state energies as functions of the phase difference ϕ from (4).(middle) The contributions of the filled bound states (E < 0) to the current are given by the derivatives of the bound state energies(5).(right) The total current (blue) (3) and the contribution of the continuum (orange) [integral term in(5)].The parameters are L∆/ℏv F = 1 and T = 0.

)Fig. 4
Fig. 4 Same as Figs. 2 and 3 but L∆/ℏv F = 5.7.This is an example of a bound state emerging from the continuum with increasing ϕ.

Fig. 6
Fig.6Same as Fig.4but including barriers at both NS interfaces.The noteworthy difference to Fig.4is the opening of gaps in the energy spectrum at ϕ = 0 and π.The figure is based on the calculation by Galaktionov and Zaikin[24].The parameters are the transmission probability 0.9 at both barriers, L∆/ℏv F = 5.7, zero temperature, and the phase acquired in one round trip between the two barriers equal to π/2.[Note the misprint of an extra factor sin χ in formula(16) of Ref.[24].]