text:
# Superconductive properties of thin dirty SN bilayers

## I Introduction

It is well known that the majority of metallic superconductors is well described by the classical BCS theory of superconductivity.  One of the main qualitative features of the BCS theory is a simple relation between the superconductive transition temperature $`T_c`$ and the low-temperature value of the energy gap for $`s`$-wave superconductors: $`\mathrm{\Delta }(0)=1.76T_c`$. Experimentally, violations of this simple relation are considered as a sign of some unusual pairing symmetry or even of a non-BCS pairing mechanism; many theories of unconventional superconductivity were developed during the last decade, mainly in relation with high-$`T_c`$ materials. Recently, an evident example of such a violation of the BCS theory predictions was found in experiments by Kasumov et al. who studied current-voltage characteristics of a carbon nanotube contact between two metallic bilayers (sandwiches) made of ordinary metals, tantalum and gold. The observed value of the low-temperature Josephson critical current is 40 times larger than the maximum expected (Ambegaokar-Baratoff) value  $`I_c=\pi \mathrm{\Delta }(0)/2eR_{\mathrm{tube}}`$, where the energy gap of the bilayer $`\mathrm{\Delta }(0)`$ is estimated from its transition temperature. The source of such discrepancy is not clear at present. The most recent experiments  demonstrate the existence of intrinsic superconductivity in carbon nanotubes. However, the discrepancy could be due to unusual superconductive properties of the bilayers. The aim of the present paper is to investigate these properties.

An essential feature of the experiment \[\] was that the superconductive layer in the bilayer was very thin ($`d_S/d_N`$ = 5 nm/100 nm = 1/20). In the present paper, we investigate such a bilayer both analytically and numerically, calculating quantities characterizing the superconductivity in this proximity system: the order parameter $`\mathrm{\Delta }`$, the density of the superconducting electrons $`n`$, the critical temperature $`T_c`$, the (mini)gap $`E_g`$ in the single-particle density of states (DOS) and the DOS $`\nu (E)`$ itself, the critical magnetic field $`H_c`$ parallel to the bilayer and the upper critical field $`H_{c2}`$ perpendicular to the bilayer. In our calculations, the parameter controlling the strength of the proximity effect is the (dimensionless) resistance of the SN interface per channel $`\rho _{\mathrm{int}}`$, which is related to the total interface resistance $`R_{\mathrm{int}}`$ as

$$R_{\mathrm{int}}=\frac{R_q\rho _{\mathrm{int}}}{2N_{\mathrm{ch}}},$$
(1)
where $`R_q=h/e^2`$ is the quantum resistance, and $`N_{\mathrm{ch}}=𝒜/(\lambda _F/2)^2`$, with $`\lambda _F`$ being the Fermi wave-length, is the number of channels in the interface of area $`𝒜`$. We choose $`\lambda _F`$ referring to the S layer.

Our results show that the values of the interface resistance can be divided into three ranges: (a) at large resistance, many characteristics of the superconductor \[$`\mathrm{\Delta }`$, $`n`$, $`T_c`$, $`\nu (E)`$, $`H_c`$, $`H_{c2}`$\] are almost unaffected by the presence of the normal layer, i.e., this is the BCS limit (however, we note that $`E_g`$ does not coincide with the order parameter, and even vanishes as $`\rho _{\mathrm{int}}`$ increases); (b) at low resistance, the theory describing the bilayer is of BCS type but with the order parameter $`\mathrm{\Delta }`$ substituted by the minigap $`E_g`$ (for instance, $`E_g=1.76T_c`$, whereas $`E_g\mathrm{\Delta }`$); the original BCS relations are thus severely violated; (c) at intermediate resistance, the behavior of the system interpolates between the above two regimes.

The paper is organized as follows. In Section II, we formulate the standard technique of the Usadel equations for dirty systems,  and introduce a convenient angular parameterization of the quasiclassical Green functions entering these equations. In Section III, we apply the Usadel equations to the thin bilayer that we intend to discuss, and present numerical results for $`\mathrm{\Delta }`$, $`n`$, $`T_c`$, and $`\nu (E)`$. We start an analytical analysis of the Usadel equations for the bilayer by calculating the minigap $`E_g`$ in the density of states in the two limiting cases of high and low interface resistance (Section IV). Then, in Section V, we elucidate the structure of the theory describing the system in the so-called Anderson limit (of relatively low resistance), finding $`\mathrm{\Delta }`$, $`n`$, $`T_c`$, $`E_g`$, and $`\nu (E)`$. The critical magnetic fields $`H_c`$ (parallel to the bilayer) and $`H_{c2}`$ (perpendicular to the bilayer) are calculated in Sections VI and VII, respectively. In Section VIII, we show that the results obtained for the bilayer also apply to more complicated structures such as SNS and NSN trilayers, SNINS and NSISN systems, and SN superlattices. The relationship between our results and the experiment \[\] that has stimulated our research is discussed in Section IX. Finally, we present our conclusions (Section X).

## II Method

### A Usadel equation

Equilibrium properties of dirty systems are described  by the quasiclassical retarded Green function $`\widehat{R}(𝐫,E)`$, which is a $`2\times 2`$ matrix in the Nambu space satisfying the normalization condition $`\widehat{R}^2=\widehat{1}`$. The retarded Green function obeys the Usadel equation

$$D(\widehat{R}\widehat{R})+i[\widehat{H},\widehat{R}]=0.$$
(2)
Here the square brackets denote the commutator, $`D=\mathrm{v}l/3`$ is the diffusion constant with $`\mathrm{v}`$ and $`l`$ being the Fermi velocity and the elastic mean free path, $`\widehat{H}=E\widehat{\sigma }_z+\widehat{\mathrm{\Delta }}(𝐫)`$ with $`E`$ being the energy, whereas $`\widehat{\sigma }_z`$ (the Pauli matrix) and $`\widehat{\mathrm{\Delta }}(𝐫)`$ are given by

$$\widehat{\sigma }_z=\left(\begin{array}{cc}1& 0\\ 0& 1\end{array}\right),\widehat{\mathrm{\Delta }}=\left(\begin{array}{cc}0& \mathrm{\Delta }\\ \mathrm{\Delta }^{}& 0\end{array}\right).$$
(3)
The order parameter $`\mathrm{\Delta }(𝐫)`$ must be determined self-consistently from the equation

$$\widehat{\mathrm{\Delta }}(𝐫)=\frac{\nu _0\lambda }{4}_0^{\omega _D}𝑑E\mathrm{tanh}(\frac{E}{2T})\left[\widehat{R}(𝐫,E)\widehat{R}(𝐫,E)\right]_{\mathrm{o}.\mathrm{d}.},$$
(4)
where the subscript $`\mathrm{o}.\mathrm{d}.`$ denotes the off-diagonal part, $`\nu _0=m^2\mathrm{v}/\pi ^2`$ is the normal-metal density of states at the Fermi level, $`\lambda `$ is the effective constant of electron-electron interaction in the S layer (whereas we assume $`\lambda =0`$ and hence $`\mathrm{\Delta }=0`$ in the N layer), and integration is cut off at the Debye energy $`\omega _D`$ of the S material.

Equation (2) should be supplemented with the appropriate boundary conditions at an interface, which read

$$\sigma _l\left(\widehat{R}_l_𝐧\widehat{R}_l\right)=\sigma _r\left(\widehat{R}_r_𝐧\widehat{R}_r\right)=\frac{g_{\mathrm{int}}}{2}[\widehat{R}_l,\widehat{R}_r],$$
(5)
where the subscripts $`l`$ and $`r`$ designate the left and right electrode, respectively; $`\sigma `$ is the conductivity of a metal in the normal state, and $`g_{\mathrm{int}}=G_{\mathrm{int}}/𝒜`$ (with $`G_{\mathrm{int}}=1/R_{\mathrm{int}}`$) is the conductance of the interface per unit area when both left and right electrodes are in the normal state. $`_𝐧`$ denotes the projection of the gradient upon the unit vector $`𝐧`$ normal to the interface.

The system of units in which the Planck constant and the speed of light equal unity ($`\mathrm{}=c=1`$) is used throughout the paper.

### B Angular parameterization of the Green function

The normalization condition allows the angular parameterization of the retarded Green function:

$$\widehat{R}=\left(\begin{array}{cc}\mathrm{cos}\theta & ie^{i\phi }\mathrm{sin}\theta \\ ie^{i\phi }\mathrm{sin}\theta & \mathrm{cos}\theta \end{array}\right),$$
(6)
where $`\theta =\theta (𝐫,E)`$ is a complex angle which characterizes the pairing, and $`\phi =\phi (𝐫,E)`$ is the real superconducting phase. The off-diagonal elements of the matrix $`\widehat{R}`$ describe  the superconductive correlations, vanishing in the bulk of a normal metal ($`\theta =0`$).

The Usadel equation takes the form

$$\frac{D}{2}^2\theta +\left[iE\frac{D}{2}\left(\phi \right)^2\mathrm{cos}\theta \right]\mathrm{sin}\theta +\left|\mathrm{\Delta }\right|\mathrm{cos}\theta =0,$$
(8)
$$\left(\mathrm{sin}^2\theta \phi \right)=0.$$
(9)
The corresponding boundary conditions are

$$\sigma _l_𝐧\theta _l=g_{\mathrm{int}}\left[\mathrm{cos}(\phi _r\phi _l)\mathrm{cos}\theta _l\mathrm{sin}\theta _r\mathrm{sin}\theta _l\mathrm{cos}\theta _r\right],$$
(11)
$$\sigma _r_𝐧\theta _r=g_{\mathrm{int}}\left[\mathrm{cos}\theta _l\mathrm{sin}\theta _r\mathrm{cos}(\phi _r\phi _l)\mathrm{sin}\theta _l\mathrm{cos}\theta _r\right],$$
(12)
$`\sigma _l\mathrm{sin}^2\theta _l_𝐧\phi _l`$ $`=`$ $`\sigma _r\mathrm{sin}^2\theta _r_𝐧\phi _r`$ (13)
$`=`$ $`g_{\mathrm{int}}\mathrm{sin}(\phi _r\phi _l)\mathrm{sin}\theta _l\mathrm{sin}\theta _r.`$ (14)
The self-consistency equation for the order parameter $`\mathrm{\Delta }(𝐫)`$ takes the form

$$\mathrm{\Delta }=\frac{\nu _0\lambda }{2}_0^{\omega _D}𝑑E\mathrm{tanh}(\frac{E}{2T})Im\left[\mathrm{sin}\theta \right]e^{i\phi }.$$
(15)
The above equations are written in the absence of an external magnetic field. To take account of the magnetic field, it is sufficient to substitute the superconducting phase gradient in the Usadel equations (II B) by its gauge invariant form $`2m𝐯=\phi +2e𝐀`$, where $`𝐀`$ is the vector potential and $`𝐯`$ denotes the supercurrent velocity.

Physical properties of the system can be expressed in terms of the pairing angle $`\theta (𝐫,E)`$. The single-particle density of states $`\nu (𝐫,E)`$ and the density of the superconducting electrons $`n(𝐫)`$ are given by

$`\nu `$ $`=`$ $`\nu _0Re[\mathrm{cos}\theta ],`$ (16)
$`n`$ $`=`$ $`{\displaystyle \frac{2m\sigma }{e^2}}{\displaystyle _0^{\mathrm{}}}𝑑E\mathrm{tanh}({\displaystyle \frac{E}{2T}})Im\left[\mathrm{sin}^2\theta \right],`$ (17)
where $`m`$ and $`e`$ are the electron’s mass and the absolute value of its charge. The total number of single-particle states in a metal is the same in the superconducting and normal states, which is expressed by the constraint

$$_0^{\mathrm{}}𝑑E\left[\nu (𝐫,E)\nu _0\right]=0.$$
(18)
### C Simple example: the BCS case

The simplest illustration for the above technique is the BCS case, when the order parameter $`\mathrm{\Delta }(𝐫)=\mathrm{\Delta }_{BCS}`$ is spatially constant. Its phase can be set equal to zero, $`\phi =0`$. Then the Usadel equations (II B) are trivially solved, and we can write the answer in terms of the sine and the cosine of the pairing angle:

$`\mathrm{sin}\theta _{BCS}(E)`$ $`=`$ $`{\displaystyle \frac{i\mathrm{\Delta }_{BCS}}{\sqrt{E^2\mathrm{\Delta }_{BCS}^2}}},`$ (20)
$`\mathrm{cos}\theta _{BCS}(E)`$ $`=`$ $`{\displaystyle \frac{E}{\sqrt{E^2\mathrm{\Delta }_{BCS}^2}}}.`$ (21)
An infinitesimal term $`i0`$ should be added to the energy $`E`$ to take the retarded nature of the Green function $`\widehat{R}`$ into account, which yields

$$Im\left[\mathrm{sin}^2\theta _{BCS}(E)\right]=\frac{\pi }{2}\mathrm{\Delta }_{BCS}\delta (E\mathrm{\Delta }_{BCS}).$$
(22)
The usual BCS relations are straightforwardly obtained from Eqs. (15), (16), (17) (for simplicity, we consider the case of zero temperature):

$`\mathrm{\Delta }_{BCS}`$ $`=`$ $`2\omega _D\mathrm{exp}\left({\displaystyle \frac{2}{\nu _0\lambda }}\right),`$ (23)
$`\nu _{BCS}(E)`$ $`=`$ $`\{\begin{array}{cc}0,\hfill & \text{if }E<\mathrm{\Delta }_{BCS}\hfill \\ \nu _0{\displaystyle \frac{E}{\sqrt{E^2\mathrm{\Delta }_{BCS}^2}}},\hfill & \text{if }E>\mathrm{\Delta }_{BCS}\hfill \end{array},`$ (26)
$`n_{BCS}`$ $`=`$ $`\pi {\displaystyle \frac{m\sigma }{e^2}}\mathrm{\Delta }_{BCS}.`$ (27)
The critical temperature must be determined from Eq. (15) with vanishing $`\mathrm{\Delta }(T_c)`$; the result is

$$\mathrm{\Delta }_{BCS}(0)=\frac{\pi }{\gamma }T_c^{BCS},$$
(28)
where $`\gamma 1.78`$ is Euler’s constant.

## III Usadel equations for a thin bilayer

Let us consider a SN bilayer consisting of a normal metal ($`d_N<z<0`$) in contact (at $`z=0`$) with a superconductor ($`0<z<d_S`$). We assume that the layers are thin (this assumption will be discussed in Sec. IX) and can be regarded as uniform, which allows us to set the order parameter $`\mathrm{\Delta }`$ equal to a constant in the superconductive layer (we choose its phase $`\phi `$ equal to zero). At the same time, we suppose that electron-electron interaction is absent in the normal layer: $`\lambda =0`$, hence $`\mathrm{\Delta }=0`$, although the superconductive correlations ($`\theta 0`$) exist in the N layer due to the proximity effect. The Usadel equations (II B) take the form

$`{\displaystyle \frac{D_N}{2}}{\displaystyle \frac{^2\theta _N}{z^2}}+iE\mathrm{sin}\theta _N`$ $`=`$ $`0,`$ (30)
$`{\displaystyle \frac{D_S}{2}}{\displaystyle \frac{^2\theta _S}{z^2}}+iE\mathrm{sin}\theta _S+\mathrm{\Delta }\mathrm{cos}\theta _S`$ $`=`$ $`0,`$ (31)
where $`\theta _N`$ and $`\theta _S`$ denote the pairing angle $`\theta `$ at $`z<0`$ and $`z>0`$, respectively.

The boundary conditions (9) reduce to

$$\sigma _N\frac{\theta _N}{z}=\sigma _S\frac{\theta _S}{z}=g_{\mathrm{int}}\mathrm{sin}(\theta _S\theta _N).$$
(32)
Equations (III) can be integrated once, yielding

$`{\displaystyle \frac{D_N}{4}}\left({\displaystyle \frac{\theta _N}{z}}\right)^2iE\mathrm{cos}\theta _N`$ $`=`$ $`f_N,`$ (34)
$`{\displaystyle \frac{D_S}{4}}\left({\displaystyle \frac{\theta _S}{z}}\right)^2iE\mathrm{cos}\theta _S+\mathrm{\Delta }\mathrm{sin}\theta _S`$ $`=`$ $`f_S.`$ (35)
The functions $`f_N(E)`$ and $`f_S(E)`$ are determined from the boundary condition $`\theta /z=0`$ at the nontransparent outer surfaces of the bilayer, which give

$`f_N(E)`$ $`=`$ $`iE\mathrm{cos}\theta _N(d_N,E),`$ (36)
$`f_S(E)`$ $`=`$ $`iE\mathrm{cos}\theta _S(d_S,E)+\mathrm{\Delta }\mathrm{sin}\theta _S(d_S,E).`$ (37)
Let us denote $`\theta _N(E)=\theta _N(d_N,E)`$, $`\theta _S(E)=\theta _S(d_S,E)`$. Because of the uniformity of the layers, the functions $`\theta _N(z,E)`$ and $`\theta _S(z,E)`$ are nearly spatially constant. However, in order to determine them, we should take account of their weak spatial dependence and make use of the boundary conditions at the SN interface. Substituting

$`\theta _N(z,E)`$ $`=`$ $`\theta _N(E)+\delta \theta _N(z,E),`$ (38)
$`\theta _S(z,E)`$ $`=`$ $`\theta _S(E)+\delta \theta _S(z,E)`$ (39)
into Eqs. (III) and linearizing them with respect to $`|\delta \theta _N(z,E)|`$, $`|\delta \theta _S(z,E)|1`$, we find the solution. Finally, boundary conditions at the SN interface lead to

$`i\tau _NE\mathrm{sin}\theta _N(E)`$ $`=`$ $`i\tau _SE\mathrm{sin}\theta _S(E)+\tau _S\mathrm{\Delta }\mathrm{cos}\theta _S(E)`$ (40)
$`=`$ $`\mathrm{sin}\left[\theta _S(E)\theta _N(E)\right],`$ (41)
where we have denoted $`\tau _N=2\sigma _Nd_N/D_Ng_{\mathrm{int}}`$, $`\tau _S=2\sigma _Sd_S/D_Sg_{\mathrm{int}}`$. Using the definition of the interface resistance per channel (1), we can represent these quantities as

$$\tau _N=2\pi \frac{\mathrm{v}_Nd_N}{\mathrm{v}_S^2}\rho _{\mathrm{int}},\tau _S=2\pi \frac{d_S}{\mathrm{v}_S}\rho _{\mathrm{int}},$$
(42)
with $`\mathrm{v}_N`$ and $`\mathrm{v}_S`$ being the Fermi velocities in the N and S layers. The ratio $`\tau _N/\tau _S=\mathrm{v}_Nd_N/\mathrm{v}_Sd_S`$, which is independent of interface properties, can also be interpreted as the ratio of the global densities of states (per energy interval) in the two layers considered,

$$\frac{\tau _N}{\tau _S}=\frac{𝒜d_N\nu _{0N}}{𝒜d_S\nu _{0S}}.$$
(43)
The latter interpretation will prove useful for further analysis.

Having solved the boundary conditions (40), we can determine all equilibrium properties of the system (because knowledge of $`\theta _S`$, $`\theta _N`$ implies knowledge of the retarded Green function $`\widehat{R}`$).

A useful representation of the boundary conditions (40) is obtained as follows. Excluding $`\theta _N(E)`$ from Eq. (40), we arrive at a single equation for the function $`\theta _S(E)`$, which can be written, in terms of $`Z=\mathrm{exp}(i\theta _S)`$, as a polynomial equation

$`iC_6Z^6+C_5Z^5+iC_4Z^4+C_3Z^3+iC_2Z^2+C_1Z+iC_0`$ (45)
$`=0`$
with real coefficients

$`C_6`$ $`=`$ $`\tau _NE\left({\displaystyle \frac{\tau _S}{\tau _N}}\right)^2\left[1+{\displaystyle \frac{\mathrm{\Delta }}{E}}\right]^2,`$ (46)
$`C_5`$ $`=`$ $`\left[1\left(\tau _NE\right)^2\right]\left({\displaystyle \frac{\tau _S}{\tau _N}}\right)^2\left[1+{\displaystyle \frac{\mathrm{\Delta }}{E}}\right]^21,`$ (47)
$`C_4`$ $`=`$ $`\tau _NE\left({\displaystyle \frac{\tau _S}{\tau _N}}\right)^2\left[3\left({\displaystyle \frac{\mathrm{\Delta }}{E}}\right)^2+2{\displaystyle \frac{\mathrm{\Delta }}{E}}1\right],`$ (48)
$`C_3`$ $`=`$ $`22\left[1\left(\tau _NE\right)^2\right]\left({\displaystyle \frac{\tau _S}{\tau _N}}\right)^2\left[1\left({\displaystyle \frac{\mathrm{\Delta }}{E}}\right)^2\right],`$ (49)
$`C_2`$ $`=`$ $`\tau _NE\left({\displaystyle \frac{\tau _S}{\tau _N}}\right)^2\left[3\left({\displaystyle \frac{\mathrm{\Delta }}{E}}\right)^22{\displaystyle \frac{\mathrm{\Delta }}{E}}1\right],`$ (50)
$`C_1`$ $`=`$ $`\left[1\left(\tau _NE\right)^2\right]\left({\displaystyle \frac{\tau _S}{\tau _N}}\right)^2\left[1{\displaystyle \frac{\mathrm{\Delta }}{E}}\right]^21,`$ (51)
$`C_0`$ $`=`$ $`\tau _NE\left({\displaystyle \frac{\tau _S}{\tau _N}}\right)^2\left[1{\displaystyle \frac{\mathrm{\Delta }}{E}}\right]^2.`$ (52)
During further analysis, the choice between the boundary conditions in the forms (40) and (45) will be a matter of convenience.

### A Critical temperature

The critical temperature of the bilayer $`T_c`$ is defined from the condition of vanishing of the order parameter $`\mathrm{\Delta }`$. Near $`T_c`$, the superconducting correlations are very small, $`|\theta |1`$; nevertheless, the self-consistency equation (15) has a nonzero solution $`\mathrm{\Delta }0`$.

Linearizing the boundary conditions (40) with respect to $`\theta _N`$ and $`\theta _S`$, we readily find the solution:

$$\theta _S(E)=i\frac{\mathrm{\Delta }}{E}\left(1\frac{\tau _N}{\tau _S+\tau _Ni\tau _S\tau _NE}\right).$$
(53)
Substituting this into the self-consistency equation (15) and simplifying both its sides by $`\mathrm{\Delta }`$, we obtain an equation determining $`T_c`$, which can be cast into the form

$`\mathrm{ln}{\displaystyle \frac{T_c^{BCS}}{T_c}}={\displaystyle \frac{\tau _N}{\tau _S+\tau _N}}[`$ $`\psi \left({\displaystyle \frac{1}{2}}+{\displaystyle \frac{\tau _S+\tau _N}{2\pi T_c\tau _S\tau _N}}\right)\psi \left({\displaystyle \frac{1}{2}}\right)`$ (55)
$`\mathrm{ln}\sqrt{1+\left({\displaystyle \frac{\tau _S+\tau _N}{\tau _S\tau _N\omega _D}}\right)^2}],`$
where $`\psi (x)`$ denotes the digamma function. A similar formula (except the logarithmic term in the r.h.s.) was obtained by McMillan  (see also Ref. \[\]). The logarithmic term in the r.h.s. takes account of the finiteness of the Debye energy $`\omega _D`$; it becomes important only in the limit of a perfect interface (the Cooper limit), i.e., when $`\tau _S\tau _N\omega _D/(\tau _S+\tau _N)1`$. Equation (55) can be solved numerically over the entire range of $`\rho _{\mathrm{int}}`$ (see Sec. III B); the analytical solution can be found in limiting cases (see Sec. V).

### B Numerical results

The solution of Eq. (45) can be found numerically. To this end, we solve the system of two nonlinear equations for the functions $`ReZ(E)`$ and $`ImZ(E)`$, using the modified Newton method with normalization.

The solution depends on the bilayer’s parameters: the thicknesses of the layers, characteristics of materials constituting the bilayer, and the quality of the SN interface. This dependence enters Eqs. (45), (52) via $`\tau _N`$ and $`\tau _S`$. For numerical calculations, we assume the characteristics of the bilayer to be the same as in the experiment by Kasumov et al.  The superconductive layer is made of tantalum, $`d_S=5`$ nm, and the normal layer is made of gold, $`d_N=100`$ nm. Approximate experimental values of the conductivities are  $`\sigma _S=0.01\mu \mathrm{\Omega }^1\mathrm{cm}^1`$ and $`\sigma _N=1\mu \mathrm{\Omega }^1\mathrm{cm}^1`$. In order to calculate the Fermi characteristics of tantalum and gold, we use the values of the Fermi energy $`E_F(\mathrm{Ta})=11`$ eV, $`E_F(\mathrm{Au})=5.5`$ eV, and the free electrons model.

Once the parameters have been specified, the solution of Eq. (45) depends only on the interface resistance $`\rho _{\mathrm{int}}`$. Having found the function $`Z(E)`$ \[which is equivalent to finding $`\theta _S(E)`$\], we start from the case of zero temperature, $`T=0`$, and study the dependence of the order parameter $`\mathrm{\Delta }`$ and of the superconducting electrons’ density in the S layer $`n_S`$ \[Eqs. (15), (17)\] on $`\rho _{\mathrm{int}}`$. The results are plotted in Fig. 1, where we also show the dependence of the critical temperature $`T_c`$, determined from Eq. (55), on $`\rho _{\mathrm{int}}`$.

The suppression of $`\mathrm{\Delta }`$, $`n_S`$, and $`T_c`$, in comparison to their BCS values in the S layer, is a natural consequence of proximity to the normal metal. At the same time, there is a possibility of BCS-like behavior, which implies the BCS relations between the suppressed quantities and the coincidence of the three curves plotted in Fig. 1. However, the curves split, and the difference between them is largest for relatively small values of $`\rho _{\mathrm{int}}`$. Figure 2 presents the range $`80<\rho _{\mathrm{int}}<150`$ on a larger scale.

Figure 3 shows the temperature dependence of the order parameter $`\mathrm{\Delta }`$ and of the superconducting electrons’ density in the S layer $`n_S`$. Although the smaller $`\rho _{\mathrm{int}}`$ the further it is from the BCS limit (corresponding to $`\rho _{\mathrm{int}}\mathrm{}`$), we observe that at $`\rho _{\mathrm{int}}=90`$ the curves are closer to the BCS behavior than at $`\rho _{\mathrm{int}}=110`$. An explanation of this feature is given in Sec. V.

Finally, the energy dependence of the single-particle density of states in the S and N layers $`\nu _{S,N}(E)`$ is plotted in Fig. 4. The density of states in the bilayer is qualitatively different from the BCS result (26). In particular, there is a minigap $`E_g`$ in the density of states at energies much smaller then $`\mathrm{\Delta }_{BCS}`$, and even much smaller then $`\mathrm{\Delta }`$ in the bilayer. In the next Section, we find this minigap analytically in the two limiting cases of small and large interface resistance $`\rho _{\mathrm{int}}`$. Another feature which can be seen from Fig. 4 is that the order parameter $`\mathrm{\Delta }`$ plays the role of a characteristic energy scale of the system only in the limit of large $`\rho _{\mathrm{int}}`$ (see the case $`\rho _{\mathrm{int}}=2000`$). Otherwise, no peculiarity in the DOS is observed at $`E=\mathrm{\Delta }`$. Some other aspects of the DOS behavior will be discussed in Secs. IV, V.

## IV Minigap in the density of states

In principle, the energy dependence of the single-particle density of states is different in the S and N layers. At the same time, the gap in the DOS is a property of the bilayer as a whole; the gap is spatially independent because there is no localization in the system and each electronic state extends over the whole bilayer.

The presence of the gap thus implies that the density of states vanishes in both layers, when the energy is below the gap:

$$Re[\mathrm{cos}\theta _S]=Re[\mathrm{cos}\theta _N]=0,$$
(56)
leading to $`\theta _S=\frac{\pi }{2}+i\vartheta _S`$, $`\theta _N=\frac{\pi }{2}+i\vartheta _N`$, with real $`\vartheta _S`$ and $`\vartheta _N`$. In this case, Eqs. (40) can be written as

$`\mathrm{tanh}\vartheta _N`$ $`=`$ $`{\displaystyle \frac{\mathrm{sinh}\vartheta _S+\tau _NE}{\mathrm{cosh}\vartheta _S}},`$ (58)
$`\mathrm{cosh}\vartheta _N`$ $`=`$ $`{\displaystyle \frac{\tau _S}{\tau _N}}\mathrm{cosh}\vartheta _S+{\displaystyle \frac{\tau _S\mathrm{\Delta }}{\tau _NE}}\mathrm{sinh}\vartheta _S.`$ (59)
Assuming that $`\mathrm{sinh}\vartheta _S\tau _NE`$ at small energies, from Eqs. (56) we obtain $`\vartheta _N=\vartheta _S`$, and, finally,

$$\mathrm{cos}\theta _S=\mathrm{cos}\theta _N=\frac{E}{\sqrt{E^2E_g^2}},$$
(60)
with $`E_g=\tau _S\mathrm{\Delta }/(\tau _S+\tau _N)`$. This is a BCS-like result \[cf. Eq. (21)\], although the order parameter $`\mathrm{\Delta }_{BCS}`$ is substituted by the minigap $`E_g`$. The assumption is readily checked, and we obtain

$$E_g=\frac{\tau _S}{\tau _N+\tau _S}\mathrm{\Delta },\text{if}\frac{\tau _S\tau _N\mathrm{\Delta }}{\tau _S+\tau _N}1.$$
(61)
Now we proceed to the opposite limit of large interface resistance. Assuming $`\mathrm{sinh}\vartheta _S\tau _NE`$ and $`\mathrm{sinh}\vartheta _S1`$, we solve Eqs. (56) and finally obtain

$`\mathrm{cos}\theta _S`$ $`=`$ $`{\displaystyle \frac{E}{\tau _S\mathrm{\Delta }\sqrt{E^21/\tau _N^2}}}i{\displaystyle \frac{E}{\mathrm{\Delta }}},`$ (63)
$`\mathrm{cos}\theta _N`$ $`=`$ $`{\displaystyle \frac{E}{\sqrt{E^21/\tau _N^2}}}.`$ (64)
The assumption is readily checked, and the result is

$$E_g=\frac{1}{\tau _N},\text{if}\frac{\tau _S\tau _N\mathrm{\Delta }}{\tau _S+\tau _N}1.$$
(65)
Equations (61) and (65) imply that $`E_g`$ is a nonmonotonic function of the interface resistance: with increase of $`\rho _{\mathrm{int}}`$, it first increases at small $`\rho _{\mathrm{int}}`$ \[Eq. (61)\] and then decreases at large $`\rho _{\mathrm{int}}`$ \[Eq. (65)\]. Therefore, $`E_g`$ reaches its maximum at some intermediate value of $`\rho _{\mathrm{int}}`$, corresponding to $`\tau _S\tau _N\mathrm{\Delta }/(\tau _S+\tau _N)1`$, hence $`\rho _{\mathrm{int}}140`$. Numerical results for $`E_g`$ are shown in Fig. 5.

At first sight, vanishing of the minigap $`E_g`$ in the limit of an opaque interface seems to contradict the general tendency to the BCS behavior. However, this contradiction is more apparent than real. Actually, the DOS curve for the S layer does approach the BCS result (26) in this limit, showing the standard peculiarity at $`E\mathrm{\Delta }_{BCS}`$. At the same time, below $`\mathrm{\Delta }_{BCS}`$, the DOS curve sharply drops to very small values (which are still finite in contrast to the BCS case), and turns to zero only at $`E=E_g`$. Simultaneously, the DOS in the N layer approaches the (constant) normal-metal value. The tendency to such behavior is illustrated by Fig. 4, the case $`\rho _{\mathrm{int}}=2000`$.

The results (60), (IV) are valid not only below the minigap but also right above it (in which case the real parts of $`\theta _S`$ and $`\theta _N`$ differ from $`\pi /2`$), providing a comparison between the DOS in the S and N layers. Equation (60) demonstrates the equality of the DOS in the two layers at relatively small $`\rho _{\mathrm{int}}`$ (see Fig. 4, the case $`\rho _{\mathrm{int}}=80`$). Proceeding to the limit of an opaque interface, we should note that the approximation which lead to Eqs. (IV) fails in a narrow vicinity of $`E_g`$ \[this fact does not affect the result for the minigap (65) itself, but the incorrect divergence of the DOS at $`E=E_g`$ disappears\]. Equations (IV) demonstrate that outside this region, at energies right above the minigap, the DOS in the N layer exceeds the DOS in the S layer by the large factor $`\tau _S\mathrm{\Delta }`$ (see Fig. 4, the case $`\rho _{\mathrm{int}}=2000`$).

Finally, we note that the results of the present Section are similar to those obtained by McMillan  (see also Ref. \[\]).

## V Anderson limit

In the limit of relatively low interface resistance (the Anderson limit), the theory describing the bilayer can be developed analytically. The condition defining this limit is $`\tau _S\mathrm{\Delta }`$, $`\tau _N\mathrm{\Delta }1`$.

First of all, we need to determine $`\theta (E)`$ \[or $`Z(E)`$\] solving Eq. (40) \[or Eq. (45)\] over the entire range of energies $`E`$.

In the region $`E>\mathrm{\Delta }`$, the solution of Eq. (45) can be written as $`Z=1+\delta Z`$, with $`|\delta Z|1`$. Keeping terms up to the second order in $`\delta Z`$, we obtain

$$\delta Z=\frac{\mathrm{\Delta }\left(1i\tau _NE\right)}{E\left(\frac{\tau _S+\tau _N}{\tau _S}i\tau _NE\right)}.$$
(66)
This result is general in the sense that it is valid for arbitrary values of $`\rho _{\mathrm{int}}`$.

At $`E<\mathrm{\Delta }`$, the same calculation as for the minigap leads to the result

$$\mathrm{sin}\theta _S=\mathrm{sin}\theta _N=\frac{iE_g}{\sqrt{E^2E_g^2}},$$
(67)
with the minigap $`E_g`$ given by Eq. (61). \[To avoid confusion, we note that under the less strict limitations for $`\tau _S\mathrm{\Delta }`$, $`\tau _N\mathrm{\Delta }`$ used in Eq. (61), the BCS-like results (60) and (67) are valid only up to energies of the order of $`E_g`$.\]

Now $`Im[\mathrm{sin}\theta _S]`$ is readily calculated, and in the case of zero temperature, $`T=0`$, the self-consistency equation (15) can be solved, yielding

$$\frac{\mathrm{\Delta }}{\mathrm{\Delta }_{BCS}}=\left(\frac{\tau _S+\tau _N}{\tau _S}\right)\left[\frac{\mathrm{\Delta }_{BCS}}{2\omega _D}\sqrt{1+\left(\frac{\tau _S\tau _N\omega _D}{\tau _S+\tau _N}\right)^2}\right]^{\tau _N/\tau _S}.$$
(68)
The relation between the order parameter $`\mathrm{\Delta }`$ and the minigap $`E_g`$ is given by Eq. (61), which immediately yields

$$\frac{E_g}{\mathrm{\Delta }_{BCS}}=\left[\frac{\mathrm{\Delta }_{BCS}}{2\omega _D}\sqrt{1+\left(\frac{\tau _S\tau _N\omega _D}{\tau _S+\tau _N}\right)^2}\right]^{\tau _N/\tau _S}.$$
(69)
In the limit of a perfect interface (the Cooper limit), which is defined by the condition $`\tau _S\tau _N\omega _D/(\tau _S+\tau _N)1`$, Eq. (69) reproduces classical Cooper’s result  generalized to the case of different Fermi parameters in the S and N layers:

$$E_g(\rho _{\mathrm{int}}0)=2\omega _D\mathrm{exp}\left(\frac{2}{\nu _0\lambda }\right),$$
(70)
with the effective interaction parameter

$$\nu _0\lambda =\frac{\tau _S}{\tau _S+\tau _N}\nu _{0S}\lambda .$$
(71)
This parameter can be considered as a result of averaging with the weighting factors $`\tau `$, which are proportional to the global normal-metal DOS per unit interval of energy, $`\tau 𝒜d\nu _0`$ (note that the interaction parameter is zero in the N layer).

At the same time, we would like to emphasize that the Anderson limit does not reduce to the Cooper limit with small corrections. On the contrary, due to the relation $`\mathrm{\Delta }\omega _D`$, the Cooper limit’s condition is not satisfied over the most part of the Anderson limit’s validity range; therefore, the minigap $`E_g`$ and the quantities calculated below differ drastically from the Cooper limit expressions.

Now we proceed to calculate the density of the superconducting electrons in the S layer $`n_S`$. On this way, we immediately encounter the problem that the above solution (67) of Eq. (40) at $`E<\mathrm{\Delta }`$ is not accurate enough for our purpose. In fact, as we will see below, the principal contribution to the integral (17) determining $`n_S`$ comes from a narrow region of energies near $`E_g`$. At the same time, Eq. (67) yields $`Im[\mathrm{sin}^2\theta _S]=0`$ and hence no contribution at all from $`E<\mathrm{\Delta }`$. We thus have to calculate a correction to Eq. (67). Assuming this correction to be small, we linearize Eq. (40) and obtain

$$\mathrm{sin}\theta _S=\frac{iE_g}{\sqrt{E^2E_g^2}}+X(E),$$
(72)
with

$$X(E)=\frac{\tau _NE_g}{\left(1+\frac{\tau _S}{\tau _N}\right)\left[1\left(\frac{E_g}{E}\right)^2\right]^2},$$
(73)
which is valid for all $`E<\mathrm{\Delta }`$ except for a narrow vicinity of $`E_g`$.

The accurate consideration of the minigap’s vicinity is possible due to the fact that $`|Z|1`$ in this region. We define a new dimensionless quantity $`ϵ`$ as

$$\frac{E}{E_g}=1+\frac{\left(\tau _S\mathrm{\Delta }\right)^{2/3}}{2\left(1+\frac{\tau _S}{\tau _N}\right)^{4/3}}ϵ,$$
(74)
and consider the region $`|EE_g|E_g`$. The function $`Im[\mathrm{sin}^2\theta _S]`$, which determines $`n_S`$, has a peak at $`ϵ1`$. An analysis of the coefficients (52) shows that only the terms that contain zeroth, first, and third order in $`Z`$ should be retained in Eq. (45). Then, after rescaling

$$Z=\frac{\left(\tau _S\mathrm{\Delta }\right)^{1/3}}{\left(1+\frac{\tau _S}{\tau _N}\right)^{2/3}}Y,$$
(75)
we obtain a cubic equation for the function $`Y(ϵ)`$:

$$4Y^3ϵY+i=0,$$
(76)
which can be solved analytically.

The density of the superconducting electrons $`n_S`$ at zero temperature is now readily calculated:

$`{\displaystyle \frac{n_S}{n_S^{BCS}}}`$ $`=`$ $`{\displaystyle \frac{E_g}{\mathrm{\Delta }_{BCS}}}\{1+{\displaystyle \frac{11\left(\tau _S\mathrm{\Delta }\right)^{5/6}}{2\pi \left(1+\frac{\tau _S}{\tau _N}\right)^{4/3}}}`$ (78)
$`+{\displaystyle \frac{2\tau _S\mathrm{\Delta }}{\pi \left(1+\frac{\tau _S}{\tau _N}\right)^2}}[{\displaystyle \frac{\tau _N}{\tau _S}}+2\mathrm{ln}{\displaystyle \frac{1+\frac{\tau _S}{\tau _N}}{\tau _S\mathrm{\Delta }}}]\},`$
with $`E_g`$ given by Eq. (69). The first term in the curly brackets is the principal one; the two other terms become comparable to unity only near the upper limit of applicability of Eq. (78).

The critical temperature of the bilayer $`T_c`$ is determined by Eq. (55). In the Anderson limit, $`(\tau _S+\tau _N)/\tau _S\tau _NT_c`$, and, with the use of the asymptotic form of the digamma function $`\psi (x)\mathrm{ln}x`$ at $`x1`$, we obtain

$$\frac{T_c}{T_c^{BCS}}=\frac{E_g}{\mathrm{\Delta }_{BCS}},$$
(79)
with $`E_g`$ given by Eq. (69). Interestingly, this result explains a discrepancy in the formulas for $`T_c`$ of a thin bilayer that were found by Cooper  and McMillan.  This discrepancy is discussed in the classical paper by McMillan.  We conclude that both cited results are correct, but their applicability ranges are different, although within the Anderson limit. Cooper’s result corresponds to a perfect interface, $`\tau _S\tau _N\omega _D/(\tau _S+\tau _N)1`$, whereas McMillan’s formula applies in the case $`\tau _S\tau _N\omega _D/(\tau _S+\tau _N)1`$.

Now we can discuss the general structure of the theory describing the bilayer in the Anderson limit. In the limit $`\rho _{\mathrm{int}}0`$, our results for the pairing angle $`\theta `$ (which is constant over the entire bilayer, $`\theta \theta _S=\theta _N`$) yield expressions which can be obtained from the BCS ones \[Eqs. (II C), (22)\] if we substitute the BCS order parameter $`\mathrm{\Delta }_{BCS}`$ by the bilayer’s minigap $`E_g`$. At $`\rho _{\mathrm{int}}>0`$, corrections to this simple result are small while the Anderson limit’s conditions are satisfied. Therefore, we obtain a BCS-type theory with $`E_g`$ substituting $`\mathrm{\Delta }_{BCS}`$ in all formulas.

The results of this Section immediately explain the numerical results in the limit of relatively small $`\rho _{\mathrm{int}}`$, shown in Fig. 2. As we have found, the Anderson limit implies the following relations between the quantities under discussion:

$`E_g`$ $`=`$ $`{\displaystyle \frac{\pi }{\gamma }}T_c,`$ (80)
$`n_S`$ $`=`$ $`\pi {\displaystyle \frac{m\sigma _S}{e^2}}E_g,`$ (81)
which substitute Eqs. (28) and (27). For the Ta/Au bilayer to which the numerical results refer, the Anderson limit is valid at $`\rho _{\mathrm{int}}<80`$ (we see that the values of $`\rho _{\mathrm{int}}`$ can be large although they are relatively small). Therefore, approaching $`\rho _{\mathrm{int}}=80`$, the curves $`n_S/n_S^{BCS}`$ and $`T_c/T_c^{BCS}`$ tend to coincide, and $`\mathrm{\Delta }/\mathrm{\Delta }_{BCS}`$ exceeds them by the large factor $`(1+\tau _N/\tau _S)15`$.

The temperature dependence of $`\mathrm{\Delta }`$ and $`n_S`$, shown in Fig. 3, is quite different at $`\rho _{\mathrm{int}}=90`$ and $`\rho _{\mathrm{int}}=110`$; at $`\rho _{\mathrm{int}}=90`$, the curves are much closer to the BCS behavior. This is also explained by approaching the Anderson limit, where the curves coincide with the BCS ones.

The DOS in the S and N layers coincide in the Anderson limit \[Eqs. (60), (67)\], showing the standard BCS-like peculiarity at $`E=E_g`$. The tendency to such behavior is illustrated by the DOS curve for $`\rho _{\mathrm{int}}=80`$ in Fig. 4.

## VI Parallel critical field

We proceed to calculate the critical magnetic field $`H_c`$ directed along the plane of the bilayer. As it was mentioned in Sec. II B, in the presence of an external magnetic field, the superconducting phase gradient in the Usadel equations (II B) must be substituted by its gauge invariant form, which can be expressed via the supercurrent velocity $`𝐯`$. The spatial distribution of $`𝐯`$ in the bilayer can be found as follows.

Let us direct the $`x`$-axis along the magnetic field $`𝐇`$. The supercurrents $`𝐣=en𝐯`$ are directed along the bilayer and perpendicularly to $`𝐇`$, i.e., $`𝐣=(0,j(z),0)`$ and $`𝐯=(0,v(z),0)`$. Near $`H_c`$, the magnetic field inside the bilayer is uniform, so the vector potential can be chosen as $`𝐀=(0,zH,0)`$. The supercurrent velocity distribution is determined by the equation $`\times 𝐯=e𝐇/m`$. Another essential point is the continuity of $`𝐯`$ at the SN interface, which follows from the continuity of the superconducting phase $`\phi `$ \[see the boundary condition (14)\]. The result is

$$v(z)=v_0\frac{eH}{m}z,$$
(82)
where $`v_0`$ is the supercurrent velocity at the interface, which must be determined from the condition that the total charge transfer across the bilayer’s cross-section is zero:

$$_{d_N}^{d_S}j(z)𝑑z=0,$$
(83)
leading to

$$v_0=\left(\frac{eH}{2m}\right)\frac{n_Sd_S^2n_Nd_N^2}{n_Sd_S+n_Nd_N}.$$
(84)
The density of the superconducting electrons is constant in each layer ($`n_S`$ and $`n_N`$).

Near $`H_c`$, the superconducting correlations are small, $`|\theta |1`$, and the Usadel equation (8) for the paring angle $`\theta (z,E)`$ can be linearized:

$`{\displaystyle \frac{D_N}{2}}{\displaystyle \frac{^2\theta _N}{z^2}}+\left(iE2m^2D_N𝐯^2\right)\theta _N`$ $`=`$ $`0,`$ (86)
$`{\displaystyle \frac{D_S}{2}}{\displaystyle \frac{^2\theta _S}{z^2}}+\left(iE2m^2D_S𝐯^2\right)\theta _S+\left|\mathrm{\Delta }\right|`$ $`=`$ $`0.`$ (87)
At the same time, the second Usadel equation (9) is trivial: its l.h.s. is proportional to

$$\left(\mathrm{sin}^2\theta 𝐯\right)=\mathrm{sin}2\theta \theta 𝐯+\mathrm{sin}^2\theta 𝐯,$$
(88)
where both terms vanish due to the fact that $`\theta `$ is directed along the $`z`$-axis whereas $`𝐯`$ is parallel to the $`y`$-axis.

The pairing angle $`\theta `$ is almost spatially constant in each layer; this allows us to average each of Eqs. (VI) over the thickness of the corresponding layer, obtaining

$`{\displaystyle \frac{\theta _N}{z}}|_{z=0}`$ $`=`$ $`{\displaystyle \frac{2d_N}{D_N}}\left(E_NiE\right)\theta _N,`$ (90)
$`{\displaystyle \frac{\theta _S}{z}}|_{z=0}`$ $`=`$ $`{\displaystyle \frac{2d_S}{D_S}}\left[\left(iEE_S\right)\theta _S+\left|\mathrm{\Delta }\right|\right],`$ (91)
where

$`E_N`$ $`=`$ $`2m^2D_N𝐯^2(z)_N,`$ (92)
$`E_S`$ $`=`$ $`2m^2D_S𝐯^2(z)_S`$ (93)
are $`H`$-dependent energies. Using Eqs. (82), (84), we express them via $`H_c`$ and the densities of the superconducting electrons:

$$E_S=\frac{D_Se^2H_c^2}{6}\left[d_S^2+3d_N^2\frac{n_N^2\left(d_S+d_N\right)^2}{\left(n_Sd_S+n_Nd_N\right)^2}\right],$$
(94)
and $`E_N`$ is obtained by the interchange of all the S and N indices.

Substituting (VI) into the boundary conditions (32) (which should be linearized), we find

$`\theta _N`$ $`=`$ $`\tau _S\left|\mathrm{\Delta }\right|/\{\tau _SE_S+\tau _NE_N+\tau _S\tau _NE_SE_N\tau _S\tau _NE^2`$ (97)
$`iE[\tau _S+\tau _N+\tau _S\tau _N(E_S+E_N)]\},`$
$`\theta _S`$ $`=`$ $`\left(1+\tau _NE_Ni\tau _NE\right)\theta _N.`$ (98)
The order parameter $`\mathrm{\Delta }`$ cancels out from the self-consistency equation (15). However, the resulting equation alone does not suffice for determining $`H_c(T)`$ because it contains $`E_S`$ and $`E_N`$, which are functions of $`n_N/n_S`$. Therefore, to obtain a closed system, we must consider the self-consistency equation together with the equation determining the ratio $`n_N/n_S`$; the latter equation is obtained from Eq. (17). The resulting system of two nonlinear equations for the quantities $`H_c`$ and $`n_N/n_S`$ is

$`\mathrm{ln}{\displaystyle \frac{2\omega _D}{\mathrm{\Delta }_{BCS}}}`$ $`=`$ $`{\displaystyle _0^{\omega _D}}𝑑E\mathrm{tanh}({\displaystyle \frac{E}{2T}}){\displaystyle \frac{Im\theta _S}{\left|\mathrm{\Delta }\right|}},`$ (100)
$`{\displaystyle \frac{n_N}{n_S}}`$ $`=`$ $`{\displaystyle \frac{\sigma _N_0^{\mathrm{}}𝑑E\mathrm{tanh}\left(\frac{E}{2T}\right)Im\theta _N^2}{\sigma _S_0^{\mathrm{}}𝑑E\mathrm{tanh}\left(\frac{E}{2T}\right)Im\theta _S^2}},`$ (101)
with $`\theta _N`$ and $`\theta _S`$ given by Eqs. (VI). The first equation of the system, Eq. (100), can be written via the digamma functions, thus taking exactly the same form as Eq. (135) below (which determines the perpendicular upper critical field) if we denote $`_S=E_S+1/\tau _S`$, $`_N=E_N+1/\tau _N`$.

In the limit $`\rho _{\mathrm{int}}\mathrm{}`$, Eqs. (VI) lead to the BCS result. In this case, the layers uncouple, the density of the superconducting electrons in the N layer vanishes, $`n_N/n_S0`$, and Eq. (100) finally yields

$$\mathrm{ln}\frac{T_c^{BCS}}{T}=\psi \left(\frac{1}{2}+\frac{D_S\left[eH_c^{BCS}d_S\right]^2}{12\pi T}\right)\psi \left(\frac{1}{2}\right),$$
(102)
which determines the parallel critical field $`H_c^{BCS}(T)`$ of a thin superconducting film. Another immediate consequence of Eqs. (VI) is the critical temperature of the bilayer $`T_c`$, which can be found from the condition $`H_c(T_c)=0`$: in this case, Eqs. (VI) reproduce Eq. (55).

The system of equations (VI) can be solved numerically at arbitrary values of the temperature $`T`$ and the interface resistance $`\rho _{\mathrm{int}}`$; the results for $`H_c`$ are presented in Figs. 6, 7.

A remarkable feature of the function $`H_c(\rho _{\mathrm{int}})`$ at zero temperature (Fig. 6) is the steep behavior of $`H_c`$ at $`\rho _{\mathrm{int}}=120`$$`123`$. This feature is due to rearrangement of the supercurrents inside the bilayer, which occurs in the following way. The supercurrent velocity changes across the thickness of the bilayer according to the simple linear law (82). This supercurrent distribution may be characterized by the position of the stationary point $`z_0`$, where the supercurrent velocity is zero: $`v(z_0)=0`$, hence $`z_0=mv_0/eH`$. At large values of the interface resistance $`\rho _{\mathrm{int}}`$, the density of the superconducting electrons in the N layer is very small, $`n_N/n_S1`$, and the supercurrents circulate only in the S part of the system; this case corresponds to

$$z_0=\frac{d_S}{2}.$$
(103)
Then, while decreasing $`\rho _{\mathrm{int}}`$, a shift in $`z_0`$ occurs. Now the supercurrents in the S layer are not compensated (in the sense of the charge transfer); therefore, they must be compensated by the supercurrents in the N layer, which are enhanced due to significant increase in $`n_N`$. This situation corresponds to the beginning of the drop in $`H_c`$. The ratio of the superconducting electrons’ densities grows rapidly, approaching the Anderson limit value $`n_N/n_S=\sigma _N/\sigma _S`$ (see Sec. VI A below); simultaneously, $`z_0`$ tends to

$$z_0=\frac{\sigma _Sd_S^2\sigma _Nd_N^2}{2\left(\sigma _Sd_S+\sigma _Nd_N\right)},$$
(104)
and the steep drop in $`H_c`$ finishes. For the bilayer to which the numerical results refer, $`d_Sd_N`$ and $`\sigma _S\sigma _N`$, so Eq. (104) yields $`z_0d_N/2`$.

This scenario is illustrated by Fig. 8, which has been obtained numerically.

The analytical solution of Eqs. (VI) at zero temperature in the Anderson limit is presented below.

### A $`H_c`$ at zero temperature in the Anderson limit

In the zero-temperature Anderson limit (defined by the conditions $`\tau _SE_S`$, $`\tau _NE_N1`$), the ratio of the superconducting electrons’ densities (101) becomes independent of the magnetic field, $`n_N/n_S=\sigma _N/\sigma _S`$, and the self-consistency equation (100) yields

$$\frac{\tau _SE_S+\tau _NE_N}{\tau _S+\tau _N}=\frac{\mathrm{\Delta }_{BCS}}{2}\left[\frac{\mathrm{\Delta }_{BCS}}{2\omega _D}\sqrt{1+\left(\frac{\tau _S\tau _N\omega _D}{\tau _S+\tau _N}\right)^2}\right]^{\frac{\tau _N}{\tau _S}}$$
(105)
which determines $`H_c`$. This result can be compared to the BCS case, which corresponds to the limit $`\rho _{\mathrm{int}}\mathrm{}`$. In this case, the density of the superconducting electrons in the N layer vanishes, $`n_N/n_S0`$, and the self-consistency equation yields

$$E_S^{BCS}=\frac{\mathrm{\Delta }_{BCS}}{2},$$
(106)
where $`E_S^{BCS}`$ is given by Eq. (94) with $`n_N=0`$. Finally,

$$H_c^{BCS}=\frac{\sqrt{3}\mathrm{\Phi }_0}{\pi \xi _{BCS}d_S},\xi _{BCS}=\sqrt{\frac{D_S}{\mathrm{\Delta }_{BCS}}},$$
(107)
where $`\mathrm{\Phi }_0=\pi /e`$ is the flux quantum, and $`\xi _{BCS}`$ is the correlation length in the dirty limit.

Remarking that the r.h.s. of Eq. (105) is identical to $`E_g/2`$ with the minigap $`E_g`$ given by Eq. (69), we see that equation (105), determining the parallel critical field of the bilayer in the Anderson limit, is obtained from the BCS equation (106) if we substitute the order parameter $`\mathrm{\Delta }_{BCS}`$ by the minigap $`E_g`$ (in accordance with the results of Sec. V) and the $`H`$-dependent energy $`E_S^{BCS}`$ by the corresponding averaged quantity $`(\tau _SE_S+\tau _NE_N)/(\tau _S+\tau _N)`$.

The explicit result for the parallel critical field of the bilayer, obtained from Eq. (105), can be cast into a BCS-like form:

$$H_c=\frac{\sqrt{3}\mathrm{\Phi }_0}{\pi \xi d_{\mathrm{eff}}}.$$
(108)
The bilayer’s correlation length $`\xi `$ is the characteristic space scale on which the order parameter (or the pairing angle $`\theta `$, or the Green function) varies in the absence of the magnetic field. In the Anderson limit (under discussion), the explicit formula for $`\xi `$ is a natural generalization of the BCS expression \[see Eq. (107)\] which implies that $`D_S`$ must be substituted by the averaged diffusion constant $`D`$ and $`\mathrm{\Delta }_{BCS}`$ must be substituted (in accordance with the results of Sec. V) by the bilayer’s characteristic energy scale, the minigap $`E_g`$ \[Eq. (69)\]:

$$\xi =\sqrt{\frac{D}{E_g}},D=\frac{\tau _SD_S+\tau _ND_N}{\tau _S+\tau _N}.$$
(109)
The effective thickness of the bilayer in Eq. (108) is

$`d_{\mathrm{eff}}=[(\sigma _Sd_S+\sigma _Nd_N)(\sigma _Sd_S^3+\sigma _Nd_N^3)`$ (110)
$`+3\sigma _S\sigma _Nd_Sd_N(d_S+d_N)^2]^{1/2}/(\sigma _Sd_S+\sigma _Nd_N).`$ (111)
In the case of equal conductivities, $`\sigma _S=\sigma _N`$, the effective thickness is simply the geometrical one: $`d_{\mathrm{eff}}=d_S+d_N`$. This case corresponds to a uniform density of the superconducting electrons, $`n_S=n_N`$, which implies a continuous distribution of the supercurrents, centered at the middle of the bilayer \[this can be also seen from Eq. (104) which yields $`z_0=(d_S+d_N)/2`$ in the case $`\sigma _S=\sigma _N`$\]. However, in a more subtle situation when the conductivities are different, the density of the supercurrent experiences a jump at the SN interface; this nontrivial supercurrent distribution results in the nonequivalence of $`d_{\mathrm{eff}}`$ to the geometrical thickness of the bilayer.

## VII Perpendicular upper critical field

Now we turn to calculating the upper critical field $`H_{c2}`$ perpendicular to the plane of the bilayer.

As in the case of the parallel critical field, we start with discussing the supercurrent distribution, which is now a function of the sample boundaries in the $`xy`$-plane, perpendicular to the magnetic field $`𝐇`$ (the magnetic field is directed along the $`z`$-axis). The infinite bilayer under consideration can be thought of as a disk of a large radius; let us assume $`x=0`$, $`y=0`$ at the axis of the disk. Then the supercurrent distribution is axially symmetric, and, with the gauge chosen as $`𝐀=\left[\mathrm{𝐇𝐫}\right]/2`$, the superconducting phase must be constant, $`\phi =0`$, which yields a simple result for the supercurrent velocity: $`𝐯=e𝐀/m`$.

Near $`H_{c2}`$, the superconducting correlations are small, $`|\theta |1`$, and the Usadel equations can be linearized:

$`{\displaystyle \frac{D}{2}}\left(i+2e𝐀\right)^2\theta +iE\theta +\mathrm{\Delta }`$ $`=`$ $`0,`$ (113)
$`𝐀\theta `$ $`=`$ $`0.`$ (114)
The second of these equations is trivially satisfied because $`\theta (𝐫)`$ is axially symmetric.

Thus, the Usadel equations reduce to the single Eq. (113) for the pairing angle $`\theta (𝐫,E)`$. Introducing the cylindrical coordinates $`𝐫(z,𝝆)`$ and denoting $`\widehat{𝐏}=i_𝝆+2e𝐀(𝝆)`$, we rewrite this equation as

$`{\displaystyle \frac{D_N}{2}}{\displaystyle \frac{^2\theta _N}{z^2}}{\displaystyle \frac{D_N}{2}}\widehat{𝐏}^2\theta _N+iE\theta _N`$ $`=`$ $`0,`$ (116)
$`{\displaystyle \frac{D_S}{2}}{\displaystyle \frac{^2\theta _S}{z^2}}{\displaystyle \frac{D_S}{2}}\widehat{𝐏}^2\theta _S+iE\theta _S+\mathrm{\Delta }`$ $`=`$ $`0.`$ (117)
We cannot solve these equations straightforwardly because near the upper critical field, the order parameter $`\mathrm{\Delta }(𝝆)`$ is a nontrivial unknown function of the in-plane coordinate $`𝝆`$ (while the $`z`$-dependence is absent due to the small thickness of the bilayer). In this situation, we employ the following approach.

Averaging each of Eqs. (VII) over the thickness of the corresponding layer, we obtain

$`{\displaystyle \frac{\theta _N}{z}}|_{z=0}`$ $`=`$ $`{\displaystyle \frac{2d_N}{D_N}}\left({\displaystyle \frac{D_N}{2}}\widehat{𝐏}^2\theta _NiE\theta _N\right),`$ (119)
$`{\displaystyle \frac{\theta _S}{z}}|_{z=0}`$ $`=`$ $`{\displaystyle \frac{2d_S}{D_S}}\left({\displaystyle \frac{D_S}{2}}\widehat{𝐏}^2\theta _S+iE\theta _S+\mathrm{\Delta }\right).`$ (120)
The averaged pairing angles entering the r.h.s. of Eqs. (VII) are

$`\theta _N(𝝆,E)`$ $`=`$ $`{\displaystyle \frac{1}{d_N}}{\displaystyle _{d_N}^0}𝑑z\theta _N(z,𝝆,E),`$ (122)
$`\theta _S(𝝆,E)`$ $`=`$ $`{\displaystyle \frac{1}{d_S}}{\displaystyle _0^{d_S}}𝑑z\theta _S(z,𝝆,E).`$ (123)
Substituting Eqs. (VII) into the boundary conditions (32) (which should be linearized), we obtain a system of two differential equations for the function $`\theta (𝝆,E)`$:

$`\tau _N\left({\displaystyle \frac{D_N}{2}}\widehat{𝐏}^2\theta _NiE\theta _N\right)`$ $`=`$ $`\tau _S\left({\displaystyle \frac{D_S}{2}}\widehat{𝐏}^2\theta _S+iE\theta _S+\mathrm{\Delta }\right)`$ (124)
$`=`$ $`\theta _S\theta _N.`$ (125)
From the vicinity of the superconductive transition it follows that the pairing angle $`\theta `$ depends on the order parameter $`\mathrm{\Delta }`$ linearly:

$`\theta _N(𝝆,E)`$ $`=`$ $`{\displaystyle \frac{\mathrm{\Delta }\left(𝝆\right)}{\alpha _N\left(E\right)}},`$ (127)
$`\theta _S(𝝆,E)`$ $`=`$ $`{\displaystyle \frac{\mathrm{\Delta }\left(𝝆\right)}{\alpha _S\left(E\right)}},`$ (128)
where the functions $`\alpha _N(E)`$ and $`\alpha _S(E)`$ are spatially independent. Then Eqs. (125) can be rewritten as

$$\frac{D_N}{2}\widehat{𝐏}^2\mathrm{\Delta }(𝝆)=\left[iE+\frac{1}{\tau _N}\left(\frac{\alpha _N\left(E\right)}{\alpha _S\left(E\right)}1\right)\right]\mathrm{\Delta }(𝝆),$$
(130)
$$\frac{D_S}{2}\widehat{𝐏}^2\mathrm{\Delta }(𝝆)=\left[iE+\frac{1}{\tau _S}\left(\frac{\alpha _S\left(E\right)}{\alpha _N\left(E\right)}1\right)+\alpha _S(E)\right]\mathrm{\Delta }(𝝆).$$
(131)
We see that the order parameter must be an eigenfunction of the differential operator $`\widehat{𝐏}^2`$. Moreover, in order to obtain the largest value of $`H_{c2}`$, we should choose the eigenfunction corresponding to the lowest eigenvalue (in complete analogy with Refs. \[\]). The solution of the emerging eigenvalue problem is readily found thanks to its formal equivalence to the problem of determining the Landau levels of a two-dimensional particle with the “mass” $`1/D`$ and the charge $`2e`$ in the uniform magnetic field $`𝐇`$ directed along the third dimension. The lowest Landau level is $`DeH`$; the function $`\alpha _S(E)`$ is straightforwardly determined,

$$\alpha _S(E)=D_SeHiE+\frac{\tau _N\left(D_NeHiE\right)}{\tau _S\left[1+\tau _N\left(D_NeHiE\right)\right]},$$
(132)
and we substitute $`\theta _S(𝝆,E)`$ into the self-consistency equation (15). The order parameter $`\mathrm{\Delta }(𝝆)`$ cancels out, and the resulting equation, which determines $`H_{c2}(T)`$, can be cast into the form

$`\mathrm{ln}{\displaystyle \frac{T_c^{BCS}}{T}}`$ $`=`$ $`{\displaystyle \frac{\tau _N}{\tau _S+\tau _N}}\mathrm{ln}\sqrt{1+\left({\displaystyle \frac{\tau _S+\tau _N}{\tau _S\tau _N\omega _D}}\right)^2}\psi \left({\displaystyle \frac{1}{2}}\right)`$ (135)
$`+{\displaystyle \frac{1}{2}}\left[1+{\displaystyle \frac{_S_N}{\sqrt{\left(_S_N\right)^2+4/\tau _S\tau _N}}}\right]\psi \left({\displaystyle \frac{1}{2}}+{\displaystyle \frac{1}{4\pi T}}\left[_S+_N+\sqrt{\left(_S_N\right)^2+{\displaystyle \frac{4}{\tau _S\tau _N}}}\right]\right)`$
$`+{\displaystyle \frac{1}{2}}\left[1{\displaystyle \frac{_S_N}{\sqrt{\left(_S_N\right)^2+4/\tau _S\tau _N}}}\right]\psi \left({\displaystyle \frac{1}{2}}+{\displaystyle \frac{1}{4\pi T}}\left[_S+_N\sqrt{\left(_S_N\right)^2+{\displaystyle \frac{4}{\tau _S\tau _N}}}\right]\right),`$
where

$`_S`$ $`=`$ $`D_SeH_{c2}+{\displaystyle \frac{1}{\tau _S}},`$ (136)
$`_N`$ $`=`$ $`D_NeH_{c2}+{\displaystyle \frac{1}{\tau _N}}`$ (137)
are $`H`$-dependent energies. The logarithmic term in the r.h.s. of Eq. (135) takes account of the finiteness of the Debye energy $`\omega _D`$; it becomes important only in the limit of a perfect interface (the Cooper limit), i.e., when $`\tau _S\tau _N\omega _D/(\tau _S+\tau _N)1`$.

In the limit $`\rho _{\mathrm{int}}\mathrm{}`$, Eq. (135) yields the classical result of Maki  and de Gennes  for the BCS case,

$$\mathrm{ln}\frac{T_c^{BCS}}{T}=\psi \left(\frac{1}{2}+\frac{D_SeH_{c2}^{BCS}}{2\pi T}\right)\psi \left(\frac{1}{2}\right),$$
(138)
which is valid for bulk superconductors and superconductive layers of arbitrary thickness (when the magnetic field is directed perpendicularly to them). Another immediate consequence of Eq. (135) is the critical temperature of the bilayer $`T_c`$, which can be found from the condition $`H_{c2}(T_c)=0`$: in this case, Eq. (135) reproduces Eq. (55).

Equation (135) can be solved numerically at arbitrary values of the temperature $`T`$ and the interface resistance $`\rho _{\mathrm{int}}`$; the results for $`H_{c2}`$ are presented in Figs. 9, 10.

The analytical solution of Eq. (135) at zero temperature in the Anderson limit is presented below.

### A $`H_{c2}`$ at zero temperature in the Anderson limit

In the zero-temperature Anderson limit (defined by the conditions $`D_SeH_{c2}1/\tau _S`$, $`D_NeH_{c2}1/\tau _N`$), Eq. (135) yields

$$\frac{H_{c2}}{H_{c2}^{BCS}}=\frac{\left(\tau _S+\tau _N\right)D_S}{\tau _SD_S+\tau _ND_N}\left[\frac{\mathrm{\Delta }_{BCS}}{2\omega _D}\sqrt{1+\left(\frac{\tau _S\tau _N\omega _D}{\tau _S+\tau _N}\right)^2}\right]^{\frac{\tau _N}{\tau _S}}$$
(139)
where the zero-temperature BCS value of the upper critical field, as follows from Eq. (138), is

$$H_{c2}^{BCS}=\frac{\mathrm{\Delta }_{BCS}}{2eD_S}=\frac{\mathrm{\Phi }_0}{2\pi \xi _{BCS}^2}.$$
(140)
It is instructive to rewrite the perpendicular upper critical field of the bilayer (139) in the standard BCS-like form

$$H_{c2}=\frac{\mathrm{\Phi }_0}{2\pi \xi ^2},$$
(141)
where $`\xi `$ is the bilayer’s correlation length given by Eq. (109) \[the physical interpretation of this result for $`\xi `$ precedes Eq. (109)\].

## VIII SNS, NSN, SNINS, NSISN, and superlattices

Our results for $`\mathrm{\Delta }`$, $`n_S`$, $`T_c`$, $`E_g`$, $`\nu (E)`$, and $`H_{c2}`$ (i.e., all the results except $`H_c`$) can be directly applied to more complicated structures such as SNS and NSN trilayers, SNINS and NSISN systems, and SN superlattices.

Let us consider, for example, a symmetric SNS trilayer consisting of two identical S layers of thickness $`d_S`$ separated by a N layer of thickness $`2d_N`$. The SN interfaces can have arbitrary (but equal) resistances. As before, the $`z`$-axis is perpendicular to the plane of the structure. This trilayer can be imagined as composed of two identical bilayers perfectly joined together along the N sides. Indeed, the pairing angle $`\theta `$ has zero $`z`$-derivative on the outer surfaces of the bilayers, thus producing the correct (symmetric in the $`z`$-direction) solution for $`\theta `$ in the resulting trilayer. Consequently, the symmetric SNS trilayer has exactly the same physical properties \[$`\mathrm{\Delta }`$, $`n_S`$, $`T_c`$, $`E_g`$, $`\nu (E)`$, $`H_{c2}`$\] as the SN bilayer considered in the present paper. The only point where the above reasoning fails is the calculation of the parallel critical field $`H_c`$. In this case, the combination of the supercurrent distributions in the two bilayers does not yield the correct distribution in the resulting SNS trilayer, which implies that the Usadel equations for the two systems are different.

Evidently, the above reasoning, based on the formal equivalence of the outer-surface boundary condition for the bilayer to the symmetry-caused condition in the middle of the SNS trilayer, also holds for symmetric NSN trilayers (N layers of thickness $`d_N`$, S layer of thickness $`2d_S`$, identical SN interfaces) and SN superlattices (N layers of thickness $`2d_N`$, S layers of thickness $`2d_S`$, identical SN interfaces). Moreover, the same applies to systems composed of two bilayers in nonideal contact with each other: SNINS and NSISN (where I stands for an arbitrary potential barrier), because the presence of a potential barrier does not violate the applicability of the symmetry argument. Thus, all the results obtained for the bilayer (except $`H_c`$) are also valid for these structures.

## IX Discussion

An essential property of the bilayer used throughout the paper is its small thickness. Now we shall argue that the bilayer studied in the experiment by Kasumov et al.  (and to which our numerical results refer) can be considered thin. The Usadel equations (III) imply that the characteristic space scale of the bilayer’s properties variation is $`\sqrt{D_{N,S}/E_0}`$ for the N and S layers, respectively. However, the correct determination of the characteristic energy scale $`E_0`$ is a nontrivial problem. Our results suggest that $`E_0`$ is always smaller than the order parameter $`\mathrm{\Delta }`$: in the BCS limit ($`\rho _{\mathrm{int}}\mathrm{}`$), $`E_0`$ approaches $`\mathrm{\Delta }`$, whereas in the opposite (Anderson) limit, $`E_0`$ is determined by the minigap $`E_g`$ \[see Eq. (61)\]. For the following discussion, it is convenient to write the condition of the small thickness of the bilayer as

$$d_{N,S}\sqrt{\frac{\mathrm{\Delta }}{E_0}}\sqrt{\frac{\mathrm{\Delta }_{BCS}}{\mathrm{\Delta }}}\sqrt{\frac{D_{N,S}}{\mathrm{\Delta }_{BCS}}}.$$
(142)
The individual layers’ thicknesses are $`d_N=`$ 100 nm and $`d_S=`$ 5 nm. The third multiplier (the BCS correlation length) in the r.h.s. of the condition (142) equals 194 nm and 16 nm for the N and S layers, respectively. At the same time, each of the first two multipliers in the r.h.s. of the condition (142) exceeds unity. We can thus conclude that the bilayer can indeed be considered thin.

Now we turn to a possible experimental application of our results. Our results provide a method for determining $`\rho _{\mathrm{int}}`$, a very important parameter of the bilayer which is not directly measurable. By analyzing the experimental  values $`T_c=0.4`$ K and $`H_c=0.1`$ T, we get $`\rho _{\mathrm{int}}111`$ and $`\rho _{\mathrm{int}}121`$, respectively. Within the experimental accuracy of the bilayer’s parameters, the two estimates for $`\rho _{\mathrm{int}}`$ should be considered close. Interestingly, the value $`\rho _{\mathrm{int}}121`$ extracted from the measured value of $`H_c`$ corresponds to the extremely narrow region of the steep drop in $`H_c(\rho _{\mathrm{int}})`$ (see Fig. 6).

Finally, we wish to remark on a peculiarity of real systems which can be relevant when one compares our findings with an experiment. The point is that during the fabrication of a bilayer, the interface between S and N materials cannot be made ideally uniform. In other words, the local interface resistance possesses spatial fluctuations. At the same time, as we have shown, the bilayer’s properties are highly sensitive to the interface quality, which could lead to complicated behavior not reducing to the simple averaging of the interface resistance embodied in $`\rho _{\mathrm{int}}`$. One possibility could be a percolation-like proximity effect. We leave the study of inhomogeneity effects for further investigation.

## X Conclusions

We have studied, both analytically and numerically, the proximity effect in a thin SN bilayer in the dirty limit. The layers were supposed to be thin enough to ensure uniform properties of each layer across its thickness. The strength of the proximity effect is governed by $`\rho _{\mathrm{int}}`$, the resistance of the SN interface per channel.

The quantities calculated were $`\mathrm{\Delta }`$, the order parameter; $`n_S`$, the density of the superconducting electrons in the S layer; $`T_c`$, the critical temperature; $`E_g`$ and $`\nu (E)`$, the minigap in the density of states and the DOS itself; $`H_c`$ and $`H_{c2}`$, the critical magnetic field parallel to the bilayer and the upper critical field perpendicular to the bilayer.

These quantities were calculated numerically over the entire range of $`\rho _{\mathrm{int}}`$. For this purpose, the characteristics of the bilayer were assumed to be the same as in the experiment by Kasumov et al.  that originally stimulated our research (Ta/Au bilayer, $`d_S/d_N=1/20`$). In the limit of an opaque interface, $`\mathrm{\Delta }`$, $`n_S`$, $`T_c`$, $`H_c`$, and $`H_{c2}`$ approach their BCS values. At the same time, $`E_g`$ does not coincide with the order parameter $`\mathrm{\Delta }`$, and $`E_g0`$ when $`\rho _{\mathrm{int}}\mathrm{}`$, although in general, the energy dependence of the DOS in the S and N layers, $`\nu _S(E)`$ and $`\nu _N(E)`$, approaches the BCS and normal-metal results, respectively.

The minigap $`E_g`$ demonstrates nonmonotonic behavior as a function of $`\rho _{\mathrm{int}}`$. Analytical results for the two limiting cases of small and large $`\rho _{\mathrm{int}}`$ show that in the Anderson limit, $`E_g`$ increases with increasing $`\rho _{\mathrm{int}}`$, whereas in the limit of an opaque interface, $`E_g`$ tends to zero. Thus, $`E_g`$ reaches its maximum in the region of intermediate $`\rho _{\mathrm{int}}`$.

The most interesting case of relatively low interface resistance (the Anderson limit) has been considered analytically. The simple BCS relations between $`\mathrm{\Delta }`$, $`n_S`$, $`T_c`$, $`H_c`$, $`H_{c2}`$ are substituted by similar ones with $`E_g`$ standing instead of $`\mathrm{\Delta }`$. The relation between the minigap $`E_g`$ and the order parameter $`\mathrm{\Delta }`$ in this limit is expressed by Eq. (61), implying that in the case where $`\tau _S<\tau _N`$, the BCS relations are strongly violated (by more than the order of magnitude for the above-mentioned Ta/Au bilayer). The DOS in the S and N layers coincide, showing BCS-like behavior with the standard peculiarity at $`E=E_g`$. It should be emphasized that absolute values of $`\rho _{\mathrm{int}}`$ corresponding to the Anderson limit can be large; for the Ta/Au bilayer this limit is already valid at $`\rho _{\mathrm{int}}<80`$.

All the results (except $`H_c`$) obtained for the bilayer also apply to more complicated structures such as SNS and NSN trilayers, SNINS and NSISN systems, and SN superlattices.

## Acknowledgments

This research was supported by the RFBR grant # 98-02-16252.