Yajiang Chen

A A Shanenko
email: arkady.shanenko@ua.ac.be

F M Peeters

Quantum cascades in nano-engineered superconductors: geometrical, thermal and paramagnetic effects

The effect of a parallel magnetic field on the orbital motion of electrons in high-quality superconducting nanowires resulting in a superconductor-to-normal transition which occurs through a cascade of jumps in the order parameter as a function of the magnetic field. Such cascades originate from the transverse size quantization that splits the conduction band into a series of subbands. Here, based on a numerical solution of the Bogoliubov-de Gennes equations for a hollow nanocylinder, we investigate how the quantum-size cascades depend on the confining geometry, i.e., by changing the cylinder radius R and its thickness d we cover the range from the nanowire-like to the nanofilm-like regime. The cascades are shown to become much less pronounced when increasing R/d, i.e., when the nanofilm-like regime is approached. When the temperature is non-zero they are thermally smoothed. This includes the spin-magnetic-field interaction which reduces the critical (depairing) parallel magnetic field H c, but does not have any qualitative effect on the quantum cascades. From our calculations it is seen that the paramagnetic limiting field H par significantly exceeds H c, even in extremely narrow nanocylinders, i.e., when R, d are down to a few nanometers, and H c, is only about 10% larger when switching-off the spin-magnetic-field interaction in this case. Both characteristic fields, H c, and H par , exhibit pronounced quantum-size oscillations. We demonstrate that the quantum cascades and the quantum-size oscillations survive in the presence of surface roughness.

Introduction

Three basic stereotypes are associated with downsizing of superconducting devices. First, it is assumed that significant structural imperfections are inevitable in the fabrication of ultrasmall metallic specimens. Second, a lack of size control is usually thought to accompany miniaturization. And, last but not least, it is well known that superconducting correlations suffer from fluctuations in lower dimensions and, hence, the common expectation is that superconductivity will be suppressed in extremely small devices. Furthermore, the role of fluctuations is enhanced in the presence of significant disorder and so the previous investigations of superconductivity in ultrathin films and narrow wires have been mainly focused on dissipative phenomena caused by the presence of pronounced superconducting fluctuations. In particular, phase-slip fluctuations have been investigated at length in experiments with amorphous superconducting wires [1], and the results were found to be in good agreement with existing theoretical models [2][3][4][5] supporting the quantum-phase-slip scenario. Another important aspect was the study of fluctuations next to the point of the superconductor-to-normal quantum-phase transition induced by various pair-breaking perturbations (including a parallel magnetic field in nanowires) in quasi-1D and quasi-2D dirty superconducting systems [6].

However, recent breakthroughs in nanofabrication introduced an essentially new motive in the study of superconductivity in nanosized systems. In particular, high-quality nanoscale superconducting systems are now experimentally attainable, such as Sn and Al nanowires (both made of strongly coupled grains [7][8][9] and single-crystalline [10][11][12]), with diameters down to 8-10 nm, and single-crystalline superconducting Pb and In nanofilms [13][14][15][16][17][18][19][20][21][22], with thickness down to a few monolayers. Structural imperfections were minimized so that such nanowires and nanofilms did not show signatures of suppression of superconductivity due to disorder (see, e.g., [15]). For example, polycrystalline Al and single-crystalline Sn nanowires exhibit a clear increase of the critical temperature when decreasing the nanowire width [7][8][9][10][11][12], as opposed to the typical results for amorphous MoGe wires [1]. Use of the so-called quantum-growth mode [23,24] (invoking the impact of quantum-size effects that overcome the strain contribution for soft metals like Pb) for metallic nanofilms allowed experimentalists to control the film thickness on an atomic scale. A striking finding is that superconductivity is very robust against fluctuations in high-quality nanosized superconductors. For instance, atomically uniform flat islands made of Pb and In on a silicon substrate were found to be superconducting even for a thickness of only a single atomic layer [22], i.e., the worlds thinnest superconductor. Moreover, the reported temperature dependence of the excitation gap shows no significant deviations from the mean-field result [21,22]. In high-quality nanowires the role of superconducting fluctuations is, of course, more significant. However, a substantial corruption of the superconducting state due to phase fluctuations occurs only for diameters less than ≈8 nm (see [8]). We also note that typical variations of the diameter along the nanowire could be reduced to 10-15%.

The question arises as to what is the main physical mechanism underlying the pronounced difference between those high-quality nanowires and nanofilms and the previously investigated superconducting films and wires. The mean-free path in high-quality nanofilms and nanowires is approximately equal to or larger than the thickness D (see, for instance, [8,9,15]), i.e., D. So, the perpendicular motion of electrons cannot be considered in a diffusive manner, contrary to strongly disordered films and wires with D. As the perpendicular motion is governed by quantum confinement, the conduction band splits up into a series of subbands (see, e.g., [16,21]). This leads naturally to multiband superconductivity similar to, e.g., the two-band superconductivity as found in novel superconductors such as MgB 2 (see, e.g., [26]) and the iron-pnictides (see, e.g., [27]). Note that scattering of electrons on imperfections results in a broadening of the perpendicular single-electron levels. When increasing disorder, such a broadening will eventually wash out any signatures of the multiple-subband structure, and the electronic motion will become diffusive.

The formation of multiple subbands has a significant effect on the superconducting properties when the intersubband energy spacing δ sub ≈ h2 2m e π 2 a 2 (with a the nanowire/nanofilm thickness) becomes of the order of the bulk energy gap bulk . For conventional superconductors, say, aluminum or tin, bulk = 0.3-0.7 meV (see, e.g., textbooks [START_REF] Fetter | Quantum Theory of Many-Particle Systems[END_REF][START_REF] De Gennes | Superconductivity of Metals and Alloys[END_REF]) and, so, δ sub exceeds bulk when a 30-40 nm. Before being eventually suppressed in lower dimensions, superconductivity in high-quality nanofilms/nanowires exhibits an almost unexplored regime governed by the formation of multiple subbands.

The presence of multiple subbands causes many effects that may not be observed in bulk superconductors, e.g., quantum-size superconducting resonances and related quantum-size oscillations of the basic superconducting properties [START_REF] Blatt | [END_REF][31][32][33][34][35][36], the formation of Andreev states induced by quantum confinement [37,38], giant size-dependent drops of the Bardeen-Cooper-Schrieffer (BCS) coherence length [25] etc. In the present paper we focus on another interesting effect due to multiple subbands in nanoscale superconductors, i.e., quantum cascades induced by a magnetic field. It was recently shown [39] that the superconducting-to-normal transition in quantum nanowires driven by a parallel magnetic field occurs as a cascade of jumps in the order parameter as a function of the magnetic field. The magnetic field does not 'kill' superconductivity simultaneously in all relevant single-electron subbands and the field-induced decay of superconducting correlations occurs step by step, i.e., first in one subband, then in another and so on. This results in a cascade of transitions. At zero temperature the cascade structure was predicted [39] to be pronounced in nanowires with diameters 15 nm, which is an experimentally reachable domain. Here we extend [39] and consider the impact of the confining geometry and of thermal effects on these quantum cascades. Additionally, we take into account the spin-magnetic-field interaction and study quantum-size oscillations of the paramagnetic limiting field [41]. Our work is based on a numerical solution of the Bogoliubov-de Gennes (BdG) equations for a hollow metallic nanocylinder. Varying its inner radius R and its wall thickness d, we can switch from a nanowire-like (R/d 1) to a nanofilm-like (R/d 1) scenario and back [40]. We also investigate how the cascades are influenced by surface roughness, which is the main issue of minor disorder in high-quality nanoscale superconductors.

Our paper is organized as follows. In section 2 we outline the formalism of the BdG equations for a hollow metallic nanocylinder in the presence of a parallel magnetic field (with only the interaction of the orbital-magnetic moments with the external magnetic field). Further in this section we discuss how the scenario of the superconductor-to-normal transition, i.e., quantum-size cascades of jumps in the order parameter as a function of the magnetic field, depends on the confining geometry and on the multiple-subband structure induced by quantum confinement. Then we investigate thermally activated smoothing of the cascades. In section 3 the spin-magnetic-field interaction is included and we briefly discuss the spin-generalized BdG equations for a hollow nanocylinder. Based on this generalization, we investigate how the inclusion of the spin-magnetic-field interaction influences our results. By switching-off the potential energy associated with the orbital-magnetic moments, we further consider the paramagnetic breakdown of superconductivity in the presence of quantization of the perpendicular electron motion. Effects of quantum confinement on the paramagnetic limiting field, including its quantum-size oscillations, are discussed here. In section 4 we investigate smoothing of the cascades by surface roughness and discuss the effects of fluctuations. Conclusions are given in section 5.

Orbital-magnetic moments and quantum-size cascades

Formalism

Previously, we showed that the effect of a parallel magnetic field on the electron orbital motion in cylindrical superconducting nanowires results in an intriguing scenario of a superconductor-to-normal transition at zero temperature, i.e., quantum cascades [39]. Here we generalize our previous study to the case of a metallic hollow nanocylinder in a parallel magnetic field, and we investigate how the confining geometry (by changing the cylinder radius R and its wall thickness d) influences these quantum-size cascades. Effects of quantum confinement are pronounced at low temperatures and, therefore, we expect thermal smoothing of the cascades, which will also be studied in the present section. For the moment we neglect spin-magnetic-field interaction and related paramagnetic effects. Their contribution will be investigated in section 3.

Bogoliubov-de Gennes equations.

As the translational invariance in the radial direction of the nanocylinder is broken, the superconducting order parameter will be position dependent, i.e., = (r). A powerful way to implement a mean-field treatment in such a case is by invoking the BdG equations, which, in the most general form, can be written as [START_REF] De Gennes | Superconductivity of Metals and Alloys[END_REF]

H e * -H * e u ν v ν = E ν u ν v ν , (1)
with E ν the quasiparticle energy, u ν (r) and v ν (r) the particle-and hole-like wavefunctions. In the clean limit, the single-electron Hamiltonian reads H e = 1 2m e (-i h∇ -e c A) 2 µ, with m e the electron-band mass (which we put equal to the free-electron mass) and µ the chemical potential. The BdG equations are solved in a self-consistent manner, together with the self-consistency relation

(r) = g ν u ν (r)v * ν (r) 1 -2f ν , (2)
where f ν is the quasiparticle Fermi distribution, g > 0 is the coupling constant. The summation in equation ( 2) is over the states with |ξ ν | < hω D , where ω D is the Debye frequency and ξ ν is the single-electron energy given by

ξ ν = d 3 r u * ν (r)H e | A=0 u ν (r) + v * ν (r)H e | A=0 v ν (r) . (3)
Note that use of H e instead of H e | A=0 in equation ( 3) (and so in the selection procedure |ξ ν | < hω D ) will produce practically the same numerical results. However, the choice of H e | A=0 is more convenient for analytical calculations (see equation (14) together with the related discussion). We also note that the usual account for only positive quasiparticle energies in equation ( 2) is not relevant in the presence of the depairing induced by a magnetic field (or a supercurrent).

Quasiparticles with negative energies are relevant in this case because they survive at zero temperature and, thus, cause a reconstruction of the ground state. Such a reconstruction is a result of a depletion of the superconducting condensate due to the depairing. To take account of the depairing reconstruction, one should select the states for which the quasiparticle energy becomes positive only in the limit A → 0, i.e., this implies that the physical solution satisfies (see, for instance, [39,[START_REF] Swidzinsky | Spatially Inhomogeneous Problems in the Theory of Superconductivity[END_REF][START_REF] Zagoskin | Quantum Theory of Many-Body Systems: Techniques and Applications[END_REF][START_REF] Peeters | Nanoscale Superconductivity[END_REF]53])

E ν | A→0 > 0. (4)
A further insight concerning the choice of the physical states can be obtained from our consideration of Anderson's recipe below (see equation ( 12)).

Due to the grand-canonical formulation of the BdG equations, it is convenient to measure all the relevant energies in the problem with respect to µ. Given the mean-electron density n e , the chemical potential µ is calculated from

n e = 2 V ν d 3 r f ν |u ν (r)| 2 + (1 -f ν )|v ν (r)| 2 , (5)
with V the system volume. For elemental superconductors, like, e.g., aluminum, the pairing energy is typically much smaller than the Fermi energy E F . Hence, changes in the chemical potential due to the superconducting order/temperature are not important when solving equation (5) and, therefore, we set temperature to zero and neglect the superconducting order in equation (5). This gives µ = E F in a numerical solution of the BdG equations.

Anderson's recipe.

Following the seminal paper by Anderson [46] it was known that the pairing of the time-reversed states is a good approximation in the absence of a magnetic field. For nanoscale systems, e.g., superconducting nanowires and nanofilms, the use of Anderson's recipe is well justified, i.e., corrections were found to be even less than a few per cent (see, e.g., [40]). Even for a broken time-reversal symmetry it was shown [39] that Anderson's approximation still provides very good results for a cylindrical superconducting nanowire in a parallel magnetic field, because the single-electron wavefunctions remain unaltered in the linear (leading) order in the vector potential A. The same holds for a hollow nanocylinder and, so, below we use Anderson's approach. When working with the BdG equations, Anderson's approximation means that we look for a solution of the BdG equations when the particle-like u n (r) and hole-like v n (r) amplitudes are proportional to the corresponding single-electron wavefunction, i.e., (ρ, ϕ, z are the cylindrical coordinates)

u ν (r) = U jmk ψ jm (ρ) e imϕ √ 2π e ikz √ L , v ν (r) = V jmk ψ jm (ρ) e imϕ √ 2π e ikz √ L , (6)
with ν = {j, m, k}, with j = 0, 1, 2, . . . the radial quantum number, m = 0, ±1, ±2, . . . the azimuthal quantum number, and k the wavevector of the quasi-free electron motion in the z-direction (periodic boundary conditions are applied in this direction with a unit cell of length L = 50 µm). Diamagnetic contributions to A can be entirely ignored in the case of interest, which makes it possible to consider the system in a uniform external magnetic field. In our case of a parallel magnetic field H = (0, 0, H ) it is convenient to adopt the Coulomb gauge, i.e., A = H ρ 2 e ϕ , with e ϕ the unit vector in the azimuthal direction. This allows us to choose u ν (r) and v ν (r) in the factorized form of equation (6).

The problem can be further simplified due to a negligible role of the term ∝ A 2 in the single-electron Hamiltonian, which reads 1 8 m e 2 ρ 2 , with = |e|H /m e c. In particular, for a nanowire with a diameter of about 10 nm the critical (depairing) magnetic field H c, is about 1-3 T, see [39]. For this case 1 8 m e 2 ρ 2 10 -4 meV bulk = 0.25 meV. We also note that the corresponding magnetic length is a H = √ h/m e ≈ 25 nm, which is significantly larger than the radius of 5 nm. Then, when keeping only terms linear in A, ψ jm (ρ) remains the same as in the absence of a magnetic field. This is why considering only pairing of time-reversed states is still a very good approximation. Hence, ψ jm is obtained from the solution of the second-order differential equation

ξ jmk ψ jm (ρ) = h2 2m e - 1 ρ ∂ ∂ρ ρ ∂ ∂ρ + m 2 ρ 2 + k 2 -µ ψ jm (ρ). (7)
The radial part of the single-electron wavefunction ψ jm (ρ) can be expressed (see, e.g., [40]) in terms of the Bessel functions of the first and second kind with two coefficients for each set of quantum numbers (j, m). These two coefficients together with the relevant single-electron energy ξ jmk can be obtained by inserting the expression of ψ jm (ρ) into the following quantum confinement boundary conditions

ψ jm (R) = ψ jm (R + d) = 0. ( 8
)
Using these energies and based on equations ( 6) and (7), from the BdG equations we can obtain (U jmk , V jmk and jm are real)

E jmk U jmk = ξ + jmk U jmk + jm V jmk , E jmk V jmk = jm U jmk -ξ - jmk V jmk , (9)
where (µ B is the Bohr magneton)

ξ ± jmk = ξ jmk ± mµ B H , (10)
and the subband-dependent pairing potential jm = j,-m is of the form

jm = R+d R dρ ρ (ρ)ψ 2 jm (ρ), (11)
with (ρ) = * (ρ). The nontrivial physical solution of the above equations is given by

E jmk = ξ 2 jmk + 2 jm + mµ B H , (12)
where the sign '+' in front of the square root is selected by virtue of equation ( 4) and corresponds to the physical solution. In the limit A → 0 we recover the usual choice of the quasiparticle energies in the BCS model. Use of equations ( 9) and ( 12) together with the usual constraint [START_REF] De Gennes | Superconductivity of Metals and Alloys[END_REF][START_REF] Peeters | Nanoscale Superconductivity[END_REF] U 2 jmk + V 2 jmk = 1, makes it possible to express U jmk and V jmk in terms of ξ jmk and jm and, then, to rewrite equation (2) in the form of the BCS-like self-consistent equation

j m = - jmk J j m ,jm jm 2 ξ 2 jmk + 2 jm 1 -2f jmk , (13)
with the pair-interaction matrix element given by

J j m ,jm = - g 2πL R+d R dρ ρ ψ 2 j m (ρ)ψ 2 jm (ρ). (14)
Due to equation (3), we should make the summation in equation ( 13) over the states with |ξ jmk | < hω D . This is similar to the well-known regularization procedure for ultraviolet divergence in the bulk gap equation (see, for instance, [START_REF] Fetter | Quantum Theory of Many-Particle Systems[END_REF]).

Based on a numerical solution of equation (13), one can simply calculate the position-dependent order parameter from

(ρ) = g 4πL jmk ψ 2 jm (ρ) jm ξ 2 jmk + 2 jm 1 -2f jmk , (15)
which follows from equations ( 2) and (6).

Parameters of calculations.

In our calculation, we consider aluminum as the material of the hollow nanocylinder: hω D /k B = 375 K with k B the Boltzmann constant, gN(0) = 0.18 with N(0) the bulk DOS at the effective Fermi level chosen as E F = 0.9 eV. It is necessary to note that to obtain the correct period for the quantum-size superconducting oscillations in nanoscale systems with a simplified parabolic-band approximation (based on the band mass m e ), one should use an effective Fermi level rather than the real one (see a detailed discussion in the second paper of [33]). In general, the effective Fermi level depends on the complicated interplay between the crystalline and confinement directions, and its possible values for elemental superconductors (when ignoring the shift-up in E F for extremely small d) can be found in [START_REF] Shanenko | [END_REF]. For aluminum E F = 0.9 eV is justified from a good agreement with experimental data on aluminum superconducting nanowires [35].

Effect of the confining geometry

We start our discussion by considering the impact of the confining geometry on the scenario of the superconductor-tonormal transition at zero temperature. The critical field H c, is set in our numerical study as the magnetic field at which the spatially averaged order parameter, i.e.

¯ = 2 d(2R + d) R+d R dρ ρ (ρ), (16)
drops below 0.01 bulk , where bulk = 0.25 meV for the parameters in use. Our previous investigation [39] of the critical magnetic field in cylindrical nanowires demonstrated that H c, is significantly enhanced as compared to the bulk value H c,bulk = 0.01 T (see the textbook [START_REF] De Gennes | Superconductivity of Metals and Alloys[END_REF]) and exhibits pronounced quantum-size oscillations. Another finding of [39] is that the scenario of the superconductorto-normal transition in sufficiently narrow nanowires is very different from both the first-order transition in bulk for type-I superconductors (we ignore the formation of vortices in nanowires) and the smooth crossover in superconducting wires with diameters [10,11,[START_REF] Han | [END_REF] 20-200 nm (such a crossover is typical of mesoscopic samples, see [START_REF] De Gennes | Superconductivity of Metals and Alloys[END_REF]49]). As follows from the results of [39], the superconductor-to-normal transition in high-quality metallic nanowires with diameters less than 10-15 nm occurs as a cascade of jumps in the order parameter as a function of the parallel magnetic field.

2.2.1.

Quantum-size oscillations of the critical magnetic field. Based on a numerical solution of the BdG equations for a hollow nanocylinder in the presence of a parallel magnetic field, we find a rather complex scenario of the quantum-size oscillations of the critical magnetic field H c, . In figure 1 To show details about quantum-size cascades we plot the magnetic field dependence of the spatially averaged order parameter ¯ for different values of R and d in panels (b)-(f) of figures 1-3. Note that the choice of R = 0.1 nm for the nanowire-like regime can hardly be realized in practice. This choice is not crucial and is only an illustration for the nanowire-like regime: our results practically do not change with R for R < 0.1-0.3 nm and, so, we are at the limit R → 0.

From figure 1(a), we can see that the critical magnetic field can reach values as high as 20-25 T for resonant points located around d = 1.5 nm. This is by three orders of magnitude larger than H c,bulk = 0.01 T and an order of magnitude higher than the estimate H (GL) c, = 4λH c,bulk /(R+d) based on the result of the Ginzburg-Landau theory for a cylinder with diameter 2(R + d) (see [50,51], with λ ≈ 50 nm for aluminum in the clean limit). Beyond the superconducting resonances, H c, falls down to 0.5-2 T, which is close to H (GL) c, (the same was found in our previous paper for cylindrical superconducting nanowires [39]).

The reason for pronounced oscillations of H c, is as follows. Due to a systematic shift-down of the single-electron energies with increasing R and d, the bottoms of the single-electron subbands cross one by one the Fermi surface, which results in quantum-size superconducting oscillations with pronounced resonant enhancements each time when the bottom of a new subband approaches the Fermi surface. Since the seminal paper by Blatt and Thompson [START_REF] Blatt | [END_REF] (see, also, [31][32][33][34][35]) it is known that when the bottom of a single-electron subband forming due to quantization of the perpendicular electron motion is situated in the vicinity of the Fermi level, the superconducting quantities are enhanced due to an increase in the relevant DOS. Hence, the pair-condensation energy increases at resonant points and so does the critical magnetic field. In particular, as seen from figure 1(b), ¯ is about 2.05 meV (at H = 0) at the resonant point d = 1.44 nm, where H c, ≈ 24.5 T. Beyond this resonance, say, at d = 1.4 nm, the spatially averaged order parameter falls down to ≈ bulk , i.e., by an order of magnitude, and the corresponding critical field is also down by a factor of 30, i.e., to 0.7 T. We stress that there is no simple relation like H c, /H c,bulk ∝ ¯ | H =0 / bulk , as can be seen from figures 1(b)-(f). The reason is that H c, is not simply determined by the value of ¯ at H = 0. It also depends on the orbital-magnetic moments of those states that make major contributions to the superconducting characteristics. The smaller the absolute value of the relevant azimuthal quantum numbers, the larger the critical magnetic field. In figure 2(a) (R = 1 nm) the maximal critical magnetic field is about 9 T at d = 1.75 nm, while for R = 8 nm (see figure 3(a)) the maximum of H c, is about 0.75 T at d = 1.75 nm. As seen, H c, exhibits an overall decrease with increasing R and d. The reason for this is two-fold. First, the resonant enhancements of ¯ are significantly reduced for larger R and d (the superconducting resonances are eventually washed out for large enough R and d, see, for more details, [35,[START_REF] Shanenko | [END_REF]). Second, a contribution of subbands with large azimuthal quantum numbers (|m| 1) increases with R and d. Such a contribution is less stable with increasing H due to the term mµ B H next to the square root in the quasiparticle energy.

Comparing figures 1(a), 2(a) and 3(a), one finds a qualitative change in the character of the quantum-size oscillations in H c, . From our study of a hollow nanocylinder in the absence of a magnetic field, we learned that there are two distinctive regimes of the quantum-size oscillations in ¯ as a function of d: (i) an almost chaotic sequence of peaks for R/d 1, i.e., the nanowire-like behavior, and (ii) regular oscillations with an overall decrease (when increasing d) of the magnitude of resonant enhancement for R/d 1, i.e., the nanofilm-like regime. We observe a similar picture for the critical magnetic field: irregular oscillations in figure 1(a) that change into an almost regular sequence of enhancements of H c, with period d = λ F /2 ≈ 0.7 nm (λ F is the Fermi wavelength) when R/d 1. Switching between these two regimes is governed by the centrifugal energy h2 2m e m 2 ρ 2 , see equation (7). The role of the centrifugal energy is more important for small R. In this case the difference in energy between the bottoms of two neighboring single-electron subbands strongly oscillates with changing relevant quantum numbers. This leads to an irregular introduction of superconducting resonances and, as a result, to irregular oscillations of H c, (and ¯ ) with d, as seen in figure 1(a). Unlike R → 0, the centrifugal energy becomes infinitesimally small for R → ∞. So, for sufficiently large R (and d) the single-electron spectrum is almost reduced to

h2 2m e π 2 (j+1) d 2
(i.e., it is nearly degenerate with respect to m) and one can recover the regular nanofilm-like appearance of the superconducting resonances with the above-mentioned period d = λ F /2. It is worth noting that the nanowire-like and nanofilm-like regimes coexist at R = 1 nm, as seen from figure 2(a). In particular, for d = 1.5-2.0 nm (R/d ∼ 1) the overall behavior of H c, in figure 2(a) resembles that of figure 3(a). This is also seen from the fact that the maximum of H c, in both figures is at d = 1.75 nm. However, the pattern of the oscillations in figure 2(a) is still more complex as compared to figure 3(a) due to many smaller resonant enhancements that are nearly washed out at R = 8 nm. It is interesting that a similarity between figures 2(a) and 3(a) is still visible even for d = 3.5-4.0 nm (i.e., at R/d = 0.2-0.3 at R = 1.0 nm).

2.2.2.

Quantum cascades in the order parameter. Now let us discuss the results given in panels (b)-(g) of figure 1, where ¯ is shown versus H for R = 0.1 nm and different values of d. Similarly to the case of a cylindrical nanowire, we observe clear signatures of the formation of a cascade of jumps in ¯ as a function of H for the hollow nanocylinder. In particular, as seen from figure 1(b), for R = 0.1 nm and d = 1.44 nm we obtain a cascade of four jumps in the dependence of ¯ on H . The first, second and third jumps are rather small (they appear like kinks in the figure). However, the fourth jump is pronounced: ¯ drops from ≈1.7 meV down to zero.

What is the physics underlying these cascades? As already mentioned above, due to the size quantization the conduction band splits up into a series of subbands. The formation of multiple subbands results in the appearance of different quasiparticle branches. Any quasiparticle branch corresponds to the single-electron states with the same orbital-magnetic moment. Orbital-magnetic moments acquire a potential energy in the presence of a magnetic field, which results in a magnetic-field-dependent shift of the quasiparticle energy. A jump in the order parameter as a function of the magnetic field can arise when one of the relevant quasiparticle branches touches zero. From this point on, such a branch supplies the system with the states with a negative quasiparticle energy. For these quasiparticles f jmk = 1 at T = 0 and, as seen from equation (2), they survive even at zero temperature. It means that we face a reconstruction of the ground state which is nothing more than a signature of depairing in the corresponding subband. A jump is pronounced when such a subband makes a significant contribution to the order parameter. Due to the multiple-subband structure, depairing induced by a magnetic field does not occur at the same time in all single-electron subbands. It first takes place in subbands that are most sensitive to a magnetic field, i.e., the subbands with large orbital-magnetic moment (large modulus of the azimuthal quantum number). This results in a cascade of jumps (and kinks).

The first jump in figure 1(b) is related to the onset of depairing in the subband (j, m) = (0, -4). It is small because the contribution of this subband to ¯ (at H = 0) is only about 4.5%. We note that before depairing starts in subband (j, m) = (0, -4), ¯ does not exhibit any change with H : all quasiparticle energies are positive, f jmk = 0 for all relevant quantum numbers and, as a result, the order parameter does not depend on H . The second and third small jumps in panel (b) mark the onsets of depairing in subbands (j, m) = (0, -3) (contribution of about 6% at H = 0) and (j, m) = (0, -2) (contribution ≈ 5% at H = 0), respectively. The eventual drop to zero measures H c, ≈ 24 T and occurs when the bottom of the quasiparticle branch associated with subband (j, m) = (1, -2) touches zero. This is one of the resonant subbands (together with (j, m) = (1, 2)) and its contribution to ¯ is about 30% at H = 0 and approaches 40% in the vicinity of H c, (the growth is due to the depairing in subbands (j, m) = (0, -2) and (j, m) = (0, -3)). The question can arise why the third and fourth jumps in ¯ at d = 1.44 nm (see figure 1(a)) are situated at different magnetic fields in spite of the fact that there are the same azimuthal quantum numbers for the associated quasiparticle branches (j, m) = (0, -2) and (1, -2). The point is that the subband-dependent order parameters jm s are generally different due to the Andreev mechanism in nanoscale superconductors (for more detail, see [37,38,40]) and in the case of interest we have 0,-2 < 1,-2 . Our previous study of magnetically induced quantum cascades in superconducting cylindrical nanowires [39] and our present investigation of a hollow nanocylinder for R → 0 (R = 0.1 nm) show that it is almost a typical situation when the onset of depairing in a resonant subband 'kills' the superconducting solution and, hence, marks the critical field H c, (see panels (b), (d), and (f) of figure 1 for the resonant points d = 1.44, 2.12 and 3.06 nm, respectively). Note that the cascade pattern is changed in the presence of many resonant subbands, as can be seen from figure 1(e), d = 2.78 nm. For this resonant thickness there are four resonant subbands with (j, |m|) = (1, 5) (a contribution of about 23% in ¯ at H = 0) and (0, 5) (16%). Two of them, with (j, m) = (1, -5) and (0, -5), are responsible for two jumps of about 0.2-0.25 meV in ¯ .

We also remark that unlike the resonant points, the scenario of the magnetically induced breakdown of superconductivity outside superconducting resonances exhibits, as a rule, a pronounced jump down with a subsequent tail of ¯ slowly approaching zero with increasing H , see, e.g., figure 1(c). Why is there a pronounced jump in ¯ as a function of H that is still visible for non-resonant conditions in spite of the fact that there are no resonant subbands? Let us consider what happens during the decay of a superconducting resonance with increasing d. The bottom of a resonant subband moves down in energy and leaves the Debye window |ξ jmk | < hω D . As a result, the DOS at the Fermi level is reduced and a resonant enhancement disappears. However, the contribution of the former resonant subband(s) into ¯ remains larger as compared to other single-electron subbands. For instance, a large jump in ¯ in figure 1(c) measures the onset of depairing in the single-electron subband (j, m) = (2, -1), whose contribution to ¯ is about 13% at H = 0, which appears to be sufficient in order to lead to a significant jump in ¯ but it is not sufficient to 'kill' the superconducting solution without a remaining after-jump tail. Such a tail is not an artifact of our numerical calculations. In particular, the after-jump tail in figure 1(c) was self-consistently calculated to within an accuracy of about 10 -6 meV.

When switching from the nanowire-like regime to the nanofilm-like behavior by increasing the ratio R/d, see figures 2(b)-(f) and 3(b)-(f), we find significantly smaller jumps in ¯ as a function of H . In particular, figures 2(b)-(f) show patterns of the superconductor-to-normal transition for R = 1.0 nm for the five resonant thicknesses d = 1.8, 2.12, 2.8, 3.07 and 3.77 nm, respectively. It is instructive to compare these patterns with those of figures 1(b), (d)-(f) (results for the non-resonant thickness d = 1.9 nm are given in figure 1(c)). As seen, the quantum-size cascades are still well pronounced but jumps become smaller for R = 1.0 nm (and the numbers of the jumps increases). In addition, there are clear signatures of after-jump tails in ¯ in panels (e) and (f) of figure 2. This is similar to the non-resonant picture shown in figure 1(c). The reason for all the features mentioned above is a rise in the number of relevant subbands that make a contribution to the superconducting order parameter (i.e., entering the Debye window). The larger R, the smaller the centrifugal energy in equation (7), which results in a decrease of the interband energy spacing. So, more and more single-electron subbands come in play and the role of any particular subband diminishes. For instance, the typical eventual drop to zero in ¯ at H = H c, , present in most cases in figures 1 and 2, breaks into a chain of almost insignificant jumps in figures 3(b)-(f), notwithstanding whether d is close to a superconducting resonance (as d = 3.1 nm in figure 3(e)) or not (as d = 1.9, 2.04, 2.64 and 3.9 nm shown in figures 3(b)-(f), respectively). Further increases in R/d lead to a smoothing of the jumps together with the accompanying cascade structure. Here it is worth noting that the microscopic basis of the resonant enhancement for large R/d is different to those for small R/d. Indeed, for a superconducting resonance at a small ratio R/d we have, as a rule, one single-electron subband whose bottom is situated in the Debye window |ξ jmk | < hω D . For large R/d, the bottoms of many subbands are condensed in the vicinity of the level

h2 2m e π 2 (j+1) d 2
due to an extremely small centrifugal energy. When such a level occurs in the Debye window, we obtain a nanofilm-like superconducting resonance. Due to the contribution of many different azimuthal numbers, the response of the system to the external parallel magnetic field is no longer as dramatic as for small R/d.

Off-equilibrium superconductor-to-normal transition.

In this subsection we make some additional remarks on the superconductor-to-normal transition in the presence of the quantum-size cascades. In figure 4 we show the Gibbs free energy (measured from G 0 , with G 0 the Gibbs free energy of the normal state) versus H , panel (a), together with the corresponding dependence of ¯ on H , panel (b). As seen, jumps in ¯ are accompanied by discontinuities in G/G 0 (see also clear hysteretic patterns in the inserts). It is well known that when the thermodynamic potential acquires a discontinuity, the system loses its thermodynamic stability at this point. Thus, any jump in ¯ is an off-equilibrium transition associated with depairing in one of single-electron subbands. Even in the presence of smoothing (e.g., induced by an increase in R or d), the system can remain thermodynamically unstable in the vicinity of the point of a smoothed drop. This can significantly extend the domain where the quantum-size cascades of jumps in the order parameter can be observed experimentally. It is similar to depairing in the presence of a uniform supercurrent in bulk. Here the supercurrent first increases with the superfluid velocity but, after reaching its maximal value, starts to decrease quickly to zero. The regime of a continuous decrease of the supercurrent with increasing superfluid velocity cannot be realized experimentally because of a thermodynamic instability [52]. In practice, after reaching the maximal value of the supercurrent, the system will pass from the superconducting state to the normal one with a jump similar to our jumps in ¯ as a function of H (see also our recent paper about current-induced cascades [53]).

Thermal smoothing

In section 2.2, we investigated how the quantum-size cascades in the order parameter as a function of a parallel magnetic field change with the confining geometry. In particular, it was shown that cascades becomes much less pronounced with increasing ratio R/d. The reason is a significant increase in the number of relevant single-electron subbands and, so, a decrease of the relative contribution of each particular subband. As a result, the number of jumps in a cascade increases but the jumps become smaller. Then, they are converted into kinks in the dependence of ¯ on H and, finally, the kinks also disappear with a further increase in R/d.

In the present section we investigate the smoothing of the quantum-size cascades due to thermal effects. Switching on temperature means that in addition to quasiparticles with negative energy, there appears a contribution to ¯ from quasiparticles with positive energy E jmk ∼ k B T. One can expect that this will result in a smoothing of the jumps but, unlike the case of increasing R/d, the number of jumps (or smoothed jumps) in a cascade will not change with temperature. In particular, figure 5 illustrates how the cascade structure of the superconductor-to-normal transition is smoothed by temperature for the nanowire-like regime R/d 1 (R = 0.1 and d = 1.44 nm): (a) the contour plot of ¯ as a function of T/T c and H /H c, , with T c the critical temperature at zero magnetic field; (b) ¯ as a function of H /H c, for T/T c = 0, 0.3, 0.6, and 0.95; (c) ¯ as a function of T/T c for H /H c, = 0, 0.6, and 0.73. Here T c /T c,bulk = 11.5, with T c,bulk = 1.65 K for our parameters, and H c, = 24.5 T. Our numerical analysis shows that the jump of ¯ to zero (that measures the critical magnetic field as calculated at T = 0) survives up to T * /T c = 0.58 (see the solid curve in figure 5(a)). Physically, we can assume that T * is mainly determined by the magnitude of the jump in ¯ at zero temperature, i.e., ¯ | jump ≈ 1.74 meV. A naive expectation suggests that the ratio ¯ | jump /k B T * can be close to the BCS parameter 1.76, i.e., the ratio of the energy gap to the critical temperature. Interestingly, this expectation is supported by our numerical results. In figure 6 5 but for R = 8 nm and d = 2.12 nm. The dashed curve in panel (a) represents H c, (T) = 1.43H c, (1 -T/T c ) 1/2 ; the solid line in (a) shows the position of an abrupt decay of ¯ to zero by a cascade of jumps of about 0.02-0.03 meV. Panel (b) shows our results for T/T c = 0, 0.65, and 0.97 (T c /T c,bulk = 1.9 K), panel (c) displays the data for H /H c, = 0, 0.65, and 0.9 (with H c, = 0.234 T).

smoothed into continuous drops even below T * . Hence, T * measures the temperature above which the cascade structure is almost washed out due to thermal effects. We find that T * strongly oscillates with changing R and d and exhibits an overall decrease with increasing R and d, which is similar to other superconducting characteristics. This can be seen from figure 7, where our numerical results for R = 8 nm and d = 2.12 nm are given: (a), the contour plot of ¯ ; (b), the data for ¯ as a function of H /H c, for T/T c = 0, 0.65, and 0.97 (with H c, = 0.234 T and T c /T c,bulk = 1.88); (c), ¯ versus T/T c for H /H c, = 0, 0.65, and 0.9. As follows from our discussion in section 2.2.2, only small jumps survive in ¯ at such R and d. For instance, jumps of about 0.02-0.03 meV are still visible at T = 0 in the vicinity of H = H c, in figure 7(b). They are smoothed out when reaching T * = 0.08T c . So, the ratio T * /T c decreases from figures 5-7 and, in addition, T c also reduces. It is of interest to discuss the behavior of the temperature-dependent critical magnetic field H c, (T) for T → 0 (i.e., the temperature dependence of the position of the final drop to zero in ¯ as a function of H ). By numerically fitting, we find (H c, (0) = H c, )

H c, (T)/H c, (0) = 1 -αT/T c (T → 0), ( 17
)
with α ≈ 1.53 for R = 0.1 nm and d = 1.44 nm. As seen, the temperature-dependent contribution to H c, (T) is linear in T for T → 0, which has nothing to do with the temperature dependence of the thermodynamic critical field in a bulk superconductor with the leading-order temperature correction ∝ (T/T c ) 2 (see [START_REF] Fetter | Quantum Theory of Many-Particle Systems[END_REF]). In turn, numerical fits for R = 1.0 nm, d = 3.07 nm and R = 8.0 nm, d = 2.12 nm give α ≈ 1.0 and 0.73, respectively. The factor α in equation ( 17) decreases to zero with increasing R and d (our numerical results for other values of R and d shows that such a decrease is accompanied by size-dependent oscillations). The difference in the low-temperature behavior of H c, (T) as compared to the bulk thermodynamic field is directly related to the temperature dependence of ¯ for H > H * . Here, due to the presence of the large jump in ¯ as a function of H , we obtain a similar jump in ¯ as a function of T, as seen from the results for H /H c, = 0.73 in figure 5(c). Only for H < H * the ordinary BCS-like dependence of the order parameter on the temperature is recovered, see figure 5(c), the data for H /H c, = 0, 0.6, and figures 6(c) and 7(c).

For T → T c , H c, (T) goes to zero, which can be fitted to the expression

H c, (T)/H c, (0) = β(1 -T/T c ) 1/2 (T → T c ), (18)
with β ≈ 1.17, 1.02 and 1.44 for the parametric sets (R, d) = (0.1, 1.44) nm, (1.0, 3.07) nm, and (8.0, 2.12) nm, respectively. Despite the regular BCS-like temperature dependence of the order parameter ¯ ∝ (1 -T/T c ) 1/2 in the domain H /H c, 1, the critical field H c, (T) again does not exhibit any similarity with the bulk thermodynamic field (the latter is proportional to 1 -T/T c (T → T c ), see the textbook [START_REF] Fetter | Quantum Theory of Many-Particle Systems[END_REF]). The reason is that the magnetic field easily penetrates into a nanoscale superconductor, which is completely different from a bulk type-I superconductor. In the case under investigation, H c, (T) is always close to H at which the resonant subband (or the former resonant subband for the non-resonant case) begins to supply the system with negative-energy quasiparticles (thermal effects do not alter this trend but result in a smoothing of jumps in the cascade). So, based on equation (12), we can expect that H c, (T) ≈ jm /(µ B m), where j and m < 0 are the quantum numbers of a subband making a major contribution to the superconducting quantities. As jm ∼ ¯ ∝ (1 -T/T c ) 1/2 , we recover the result of equation (18). Similar to the factor α in equation ( 17), the factor β in equation ( 18) varies with R and d, as well. However, our numerical analysis does not show any signature of a systematic reduction of β. Unlike β, the factor α in equation ( 17) exhibits a significant overall decrease to zero with increasing R and d. It means that the low-temperature trend of equation ( 17) disappears together with the quantum-size cascades for large enough R and d. Equation ( 18) is still applicable until d approaches the magnetic penetration depth. From this point on we can no longer consider a homogeneous magnetic field and neglect the contribution of the diamagnetic currents.

Paramagnetic effects

Formalism

The spin-generalized BdG equations read [START_REF] Ketterson | Superconductivity[END_REF]

     H e,↑ 0 0 * -H * e,↓ 0 0 0 0 H e,↓ 0 0 * -H * e,↑           u N,↑ v N,↓ u N,↓ -v N,↑      = E N      u N,↑ v N,↓ u N,↓ -v N,↑      , ( 19
)
where N = (ν, σ ), with ν = {j, m, k} (the same as in section 2) and σ = ±1 the spin-like discrete index, and the spindependent Hamiltonian H e,α (α = {↑, ↓}) is given by

H e,α = 1 2m e -i h∇ - e c A 2 + α -µ, (20)
with ↑ = -↓ = µ B H . The self-consistency relation is now written in the form

(r) = g N [u N,↑ (r)v * N,↓ (r)(1 -f N ) (21)
+ u N,↓ (r)v * N,↑ (r)f N ]. ( 22
)
Making use of Anderson's recipe (see section 2), we put

u N,α (r) = U N,α ψ jm (ρ) e imϕ √ 2π e ikz √ L , v N,α (r) = V N,α ψ jm (ρ) e imϕ √ 2π e ikz √ L , (23)
where the factors U N,α and V N,α can be chosen real. Here it is necessary to note that the constraint U 2 jmk + V 2 jmk = 1, as used for the BdG equations in the absence of the spin-magnetic energy in section 2, is now replaced by

α U 2 N,α + V 2 N,α = 1. ( 24
)
When inserting equations ( 23) into (19), equation ( 19) is reduced to four linear algebraic equations, which form two systems of two equations with two variables, i.e.,

E N U N,↑ = ξ + jmk + ↑ U N,↑ + jm V N,↓ , E N V N,↓ = jm U N,↑ -(ξ - jmk + ↓ )V N,↓ (25)
and where ξ ± jmk and jm are defined by equations ( 10) and (11). Equations ( 25) have a nontrivial solution when E N = ξ 2 jmk + 2 jm + (m + 1)µ B H (here a sign '+' in front of the square root is taken for the physical solution, as usual). However, this quasiparticle energy is not the same as that for the nontrivial solution of equations ( 26)

E N U N,↓ = ξ + jmk + ↓ U N,↓ -jm V N,↑ , E N V N,↑ = -jm U N,↓ -ξ - jmk + ↑ V N,↑ , (26)
E N = ξ 2 jmk + 2 jm + (m -1)µ B H . It means that we obtain two different excitation branches, i.e., E jmk,σ = ξ 2 jmk + 2 jm + (m + σ )µ B H , (27)
where for σ = 1 we have U N,↑ , V N,↓ = 0 while U N,↓ = V N,↑ = 0; and for σ = -1 we obtain U N,↑ = V N,↓ = 0 and U N,↓ , V N,↑ = 0. Now, solving equations ( 25) and ( 26) together with the constraint of equation ( 24) for both quasiparticle branches, we arrive at the following spin-generalized BCS-like equation:

j m = - jmk J j m ,jm jm (1 -σ f jmk,σ ) 2 ξ 2 jmk + 2 jm , (28)
with J j m ,jm given by equation (14). Following the cut-off rule of equation ( 3), the summation in equation ( 28) runs over the states with |ξ jmk | < hω D . The position-dependent superconducting order parameter is now of the form

(ρ) = g 4πL jmk ψ 2 jm (ρ) jm (1 -σ f jmk,σ ) ξ 2 jmk + 2 jm . (29)

Spin-magnetic-field interaction

Similar to equation ( 13), the magnetic field H influences the solution of equation ( 28) only through the quasiparticle occupation number f jmk,σ . It means that jm is not sensitive to H for f jmk,σ = 0, i.e., for positive quasiparticle energies at zero temperature. So, we expect that the spin-magneticfield interaction will not qualitatively change the cascade scenario of the superconductor-to-normal transition. However, including spin-magnetic moments, we obtain an extra contribution into the quasiparticle energy that splits the excitation spectrum into the two different branches. This can transform some of the jumps in the order parameter as a function of H into smaller jumps. Another quantitative effect of the spin-magnetic-field interaction concerns the magnitude of the critical magnetic field. When the relevant azimuthal quantum numbers are large enough (i.e., for large R and d), the spin-magnetic-field potential energy cannot significantly alter H c, . However, for small R and d, when a major contribution to the order parameter comes from states with |m| < 10-20, one expects a persistent effect of spin-magnetic moments on H c, .

Let us first examine the dependence of the zerotemperature critical magnetic field on d as shown in figures 8(a) (R = 0.1 nm) and 9(a) (R = 1.0 nm). As was expected, qualitative patterns of the quantum-size oscillations in H c, do not change when including the spin-magnetic-field interaction. In particular, we again observe the same irregular sequence of resonant enhancements of H c, for the nanowire-like regime at R = 0.1 nm, see figure 8(a). The only difference is that H c, is generally smaller by 10%-20% in the presence of the spin-magnetic-field interaction (this is clearly seen at the resonant points). However, there are a few exceptions, for example, H c, increases when the spin-magnetic moments are included at, e.g., d = 1.85 nm. At first sight this looks strange. However, the underlying reason is rather simple. The spin-magnetic-field interaction splits any resonant subband into two smaller subbands. This not only transforms the final jump to zero in ¯ into smaller drops but, sometimes, also leads to the appearance of an after-jump tail in ¯ (see figure 1(c) and the related discussion). If such a tail is significant, we obtain an increase in H c, rather than a decrease, due to the definition of H c, used in the present work. Recall that H c, is set in our numerical study as the magnetic field at which ¯ (H ) drops below 0.01 bulk . The decrease in H c, , as induced by the spin-magnetic-field interaction, is less pronounced for R = 1.0 nm, see figure 9(a). Here H c, almost does not exhibit any change as compared to figure 2(a), except the domain d < 2 nm, where a general decrease of about 5-10% is observed. This is quite expected because of a significant increase in the absolute value of the relevant azimuthal quantum numbers when passing from R = 0.1 nm to 1.0 nm. Note that due to an even more significant increase in |m| for R = 8 nm, we did not find changes that are more than 1% in the corresponding H c, as compared to our results without the spin-magnetic-field interaction.

In figures 8(b)-(f) we show how the spatially averaged order parameter ¯ depends on H when including the spin-magnetic-field potential energy in the nanowire-like regime for d = 1.44, 1.7, 2.12, 2.78, and 3.06 nm (only resonant points are given). As discussed in the first paragraph of this subsection, including the spin-magneticfield interaction not only reduces, as a rule, H c, but also breaks jumps in ¯ into smaller jumps. This is clearly the case, as seen from the comparison of our results in figure 8 with those in figure 1. For instance, in figure 8(b) the order parameter decays through two large jumps while in figure 1(b) only one significant jump appears. In figure 1(e) there are two pronounced jumps in ¯ as a function of H (here two single-electron subbands with (j, m) = (0, -5) and (1, -5) are responsible for the decay of ¯ ). So, it is natural to expect that these two jumps are replaced by four smaller jumps in the presence of the spin-magnetic-field interaction, see figure 8(e). It is interesting that such a trend of doubling the number of jumps is not fully supported by the results given in figures 9(b)-(f), where ¯ as a function of H is shown for the resonant thicknesses d = 1.8, 2.12, 2.8, 3.07, and 3.77 nm. Note that there are also examples (see panels (c) and (d)) when including the spin-magnetic moments, contrary to expectation based on figures 8(b)-(f), results in the appearance of a single pronounced jump substituting for several smaller ones. This kind of behavior can be observed only for large enough R and d in the presence of many resonant subbands. Here the spin-Zeeman splitting of different but close to each other quasiparticle branches does not result in the formation of clearly distinguished subbranches any more. Instead, we can obtain bunches of almost degenerate subbranches: we note that the term µ B H (m + σ ) in equation ( 27) has the same value, e.g., for m = -3, σ = 1 and for m = -1, σ = -1. This complex interplay can lead to a significant rearrangement of jumps in ¯ and, sometimes, even to larger jumps that arise due the spin-Zeeman splitting. Another deviation from, say, the doubling rule mentioned above is for resonant enhancements that develop owing to the contribution of states with |m| = 1. In this case we will obtain a jump in ¯ driven by a downward shift of the quasiparticle branch associated with m = -1. The spin-Zeeman splitting replaces m by m + σ . For σ = 1, E jmk,σ does not include the Zeeman term any more. It appears only in E jmk,σ for σ = -1. As a result, one jump in the absence of the spin-magnetic-field interaction will be replaced by a smaller jump in its presence.

Quantum-size oscillations of the paramagnetic limiting field

It is well known that the upper limit for the critical magnetic field is controlled by the paramagnetic breakdown of the Cooper pairs [41]. To investigate how such a breakdown is sensitive to the quantization of the electron motion, here we continue our numerical study of the spin-generalized BdG equations but we switch off the orbital-magnetic interaction, i.e., the term µ B mH is removed from equation (27). As the states with negative azimuthal quantum numbers are responsible for the cascades of jumps in ¯ as a function of H , one can easily see that removing the aforementioned term makes the system less sensitive to H : the resulting paramagnetic limiting field H par will be generally larger than H discussed in section 3.2. We plot H par as calculated from the spin-generalized BdG equations versus d for R = 0.1 and 1.0 nm, see figures 10(a) and (b). Comparing the results in figure 10(a) with those of figure 8(a), we observe that H par is larger than H c, by a factor of about 2-3. The difference between H par and H c, is enhanced for larger R. In particular, H par in figure 10(b) differs from H c, by a factor of about 3-4 for the resonant points and by an order of magnitude for non-resonant thicknesses. Such a significant difference agrees with the fact that the spin-magnetic-field interaction brings only small corrections (∼10%) in addition to the effect of the orbital-magnetic moments, as found in previous subsection. Due to the size quantization of the electron motion, H par exhibits pronounced size-dependent oscillations with significant resonant enhancements over the Clogston-Chandrasekhar field H CC (the dashed line in figures 10(a) and (b)), i.e., the paramagnetic limiting field in bulk H CC = bulk /( √ 2µ B ) ≈ 3.6 T. However, when a superconducting resonance decays due to a change in d, H par can drop even below H CC .

It is instructive to get a simple analytic estimate for the paramagnetic limiting field for a superconducting nanowire. Based on equation (27) and the above discussion of our numerical results, we can find H par ≈ jm (H = 0)/µ B , where j and m are the quantum numbers associated with a resonant subband (or a former resonant subband as for the non-resonant case). In most cases jm ≈ ¯ and, so we obtain

H par ≈ ¯ (H = 0)/µ B . (30)
There are two points that makes equation ( 30) different from the Clogston-Chandrasekhar result. First, the translational invariance in the radial direction is broken and, so, the spatially averaged order parameter appears in equation [START_REF] Blatt | [END_REF] rather than the uniform pairing energy as in the expression for the Clogston-Chandrasekhar field. Second, the factor 1/ √ 2 is not present in equation [START_REF] Blatt | [END_REF]. At first sight it can lead to some inconsistency because it seems that H par of equation [START_REF] Blatt | [END_REF] does not approach H CC when, e.g., d → ∞. However, this is not true because the result of equation ( 30) is derived for R, d λ, with λ the magnetic penetration depth, i.e., when a magnetic field easily penetrates the system.

Despite the rather simplified character of equation ( 30), it provides a reasonable estimate of H par . This can be seen from figures 10(c)-(e), where the dependence of ¯ on H is shown for R = 0.1 nm and d = 2.12 nm (c), R = 0.1 nm and d = 1.7 nm (d), R = 1.0 nm and d = 2.8 nm (e), and R = 1.0 nm and d = 3.07 nm (f). For panel (c), equation [START_REF] Blatt | [END_REF] gives 19 T, which is only 10% higher than H par = 17.3 T shown in the figure. It is worth noting that the Clogston-Chandrasekhar-like estimate produces ¯ (H = 0)/( √ 2µ B ) = 13 T, which is 30% lower than the result from a numerical solution of the spin-generalized BdG equations. For panel (d), use of equation [START_REF] Blatt | [END_REF] results in the estimate H par ≈ 26.6 T, which is about 10% lower than our numerical result H par = 28.9 T. In figures 10(e) and (f) we have H par = 14.3 and 11.5 T, respectively. Here, equation [START_REF] Blatt | [END_REF] gives 15.7 and 12.7 T. As seen, corrections to the estimate given by equation [START_REF] Blatt | [END_REF] are generally within 10%. Concluding this subsection, it is worth noting that the quantum-size cascades survive in the dependence of ¯ on H even in the absence of the orbital-magnetic interaction. They appear because the quasiparticle branches with σ = -1/2 are not degenerate, due to the difference in jm . As already mentioned in section 2, this difference is a manifestation of the formation of Andreev states induced by quantum confinement (see [40]).

Surface roughness and fluctuations

Smoothing due to the surface roughness

We must mention that up to now our analysis was limited to the ideal case of a hollow nanocylinder without surface roughness and in the absence of any diameter variation along the wire. Typical experimental variations of the diameter along the wire are relatively slow and can be taken into account via an adiabatic treatment, resulting in a variation of the superconducting order parameter along the longitudinal coordinate. As a result, along the wire there will be a set of domains with an enhanced superconducting condensate (they correspond to resonant values of the diameter) interconnected by weak superconducting links (governed by non-resonant diameters). The evolution of each domain in a parallel magnetic field will exhibit a cascade decay of the local superconducting order parameter similar to that in an ideal nanocylinder.

Understanding the impact of surface roughness requires a more involved analysis. In high-quality superconducting nanowires one should distinguish the perpendicular and longitudinal motion. As the electron mean-free path is close to or even larger than the nanowire width, the perpendicular motion cannot be considered as diffusive. In this case an adequate manner to take into account the surface roughness is to invoke available models for accounting such imperfections in semiconductor nanowires, see, e.g., [54,55]. These models introduce corrections to the single-electron energies due to the scattering of electrons on impurities. Below we use the simplified approach of [START_REF] Han | [END_REF], where the single-electron energies in the presence of surface imperfections, i.e., ξ (im) jmk , are approximated as (for simplicity, let us consider R = 0, i.e., a cylinder with diameter D = 2d)

ξ (im) jmk = ξ jmk + s jmk dξ jmk dD δD, (31)
with ξ jmk given by equation ( 7) and s jmk a random variable uniformly distributed in [-0.5, 0.5]. Equation ( 31) accounts for the uncertainty in the single-particle kinetic energy due to fluctuations in the confining dimension.

We performed a numerical study of the BdG equations with equation (31). For illustration, in figure 11 the spatially averaged order parameter is shown versus H for R = 0 and d = 3.17 nm (i.e., for a cylindrical nanowire with the mean diameter D = 6.34 nm) for the ideal case of δD = 0 (triangles) and in the presence of the surface roughness with δD = 4 Å. As could be expected, surface roughness smooths the size-dependent oscillations of the DOS at the Fermi level and, as a result, makes resonant enhancements less pronounced. In particular, for δD = 4 Å the critical magnetic field H c, in figure 11 is decreased by almost a factor of 3 and ¯ (taken at H = 0) becomes smaller by a factor of 2. However, even in the presence of such significant surface imperfections the cascade decay of ¯ survives. Our choice of D = 6.34 nm is not crucial for the basic conclusions, and similar results with clear signatures of remaining quantum cascades can be obtained for other diameters, see figure 12.

Superconducting fluctuations

Another important point to address concerns the superconducting fluctuations in metallic nanowires. As mentioned in section 1, our previous numerical investigation [39] shows that the quantum cascades in ¯ as a function of H can be significant in superconducting nanowires with diameters 15 nm. However, it is known that quantum-phase slips proliferate in superconducting nanowires with D 8 nm, see, e.g., [8,9,56], which corrupts the superconducting state, leading to a long tail in the resistance below the critical temperature. This may question the relevance of our study for such narrow nanowires. However, in fact, our results are useful even for extremely narrow nanowires with the superconducting state influenced by quantum fluctuations. The point is that the superconducting coherence length ξ (governed by the mean-field order parameter) is an important quantity controlling the rate of quantum-phase slips and, so, the residual resistance significantly below T c [2][3][4]56]. In particular, when assuming that the normal state resistance of the nanowire R N is close to the resistance quantum given by R Q = π h2 2e 2 (it is the case for narrow enough nanowires), one obtains for the residual resistance R res (see, e.g., [3,8,56])

R res R N ∼ L ξ e -L/ξ ,
where L is the nanowire length. So, if the mean-field order parameter exhibits a cascade decay in the parallel magnetic field, a similar cascade structure will appear in ξ as function of H and, in turn, in the dependence of R res on H . It is also of importance to discuss thermal fluctuations. They are not extremely pronounced in high-quality superconducting nanowires. In particular, from [8] we can learn that the resistance drops rather fast below T c even for aluminum nanowires with thickness ∼ 10 nm. Such a drop exceeds an order of magnitude when T decreases from T c to 0.9 T c , see figure 2 in [8]. It means that we still have a well-pronounced transition with a thermal broadening of about δT = 0.1 T c . This broadening can provide us with an estimate of the Ginzburg-Levanyuk parameter (Gi ∼ δT T c ≈ 0.1), which measures the impact of thermal fluctuations around T c . In the situation when the conventional Ginzburg-Landau formalism is not well justified for superconducting nanowires (the nanowire diameter can be much smaller than the BCS coherence length), this gives a reasonable guess about the role of thermal fluctuations. Such an estimate of Gi shows that the mean-field treatment is still applicable for low enough temperatures even in high-quality nanowires with diameters down to 5-10 nm. It is worth noting that a naive calculation of Gi in superconducting nanowires with width ∼ 10 nm based on conventional Ginzburg-Landau theory in the clean limit, see [25], produces a similar result, i.e., Gi ∼ 0.1 for diameters ≈ 10 nm and Gi ∼ 0.2 for diameters ≈ 5 nm. We also note that high-quality Al and Sn superconducting nanowires discussed in [8][9][10][11] are far from the dirty limit of strongly disordered superconducting nanowires, e.g., made of MoGe, see [1]. In disordered nanowires the electron mean-free path is significantly smaller than the nanowire diameter, and the electron motion is diffusive. As a result, the superconducting temperature is systematically reduced with decreasing nanowire diameter. In addition, the superconducting fluctuations play a more significant role. For the transition temperature of about 3.5 K in MoGe superconducting nanowires with width ∼10 nm the thermal broadening of the resistance transition is ≈1 K, see [1]. This makes it possible to expect that Gi ≈ 0.3, which is by a factor of three larger than the above estimate for aluminum nanowires. It is in agreement with a common expectation that the role of fluctuations in disordered systems is more significant.

Conclusions

In conclusion, we investigated how quantum-size cascades appearing in the scenario of a superconducting-to-normal transition in metallic nanowires with increasing parallel magnetic field can be influenced by changing the confining geometry, by thermal smoothing and by the spin-magneticfield interaction. Our study was based on a numerical solution of the BdG equations for a hollow metallic nanocylinder. Varying its inner radius R and its wall thickness d, we observed a qualitative change in the cascades when passing from the nanowire-like regime, i.e., R/d 1, to the nanofilm-like regime, i.e., d R. In particular, while being pronounced for R/d 1, jumps in the order parameter as function of an external parallel magnetic field are significantly smoothed for R d. The reason for such a smoothing is that an increase in R results in a significant increase of the number of single-electron subbands that contribute to the order parameter. Hence, the relative contribution of each subband is strongly reduced and, as a result, the appearance of negative quasiparticle energies in one of subbands will no longer have an appreciable effect on the solution.

We also found that a smoothing of the jumps in the cascades occurs when the temperature is non-zero. However, here the reason for smoothing is different. The point is that at finite temperature there is a non-zero population of quasiparticles already at zero magnetic field, i.e., there is thermally activated depairing. Therefore, the effect of any additional depairing mechanism, e.g., due to the negative-energy excitations appearing in the magnetic field, is weakened because of such thermally activated depairing.

We found that taking into account the spin-magneticfield interaction increases the number of small jumps in a given cascade, which is due to the spin-Zeeman splitting of the quasiparticle energies. In most cases the critical magnetic field is smaller in the presence of the spin-magnetic-field interaction. However, interestingly, due to the complex interplay of the orbital and spin-magnetic-field effects, the critical magnetic field can sometimes be enhanced. The reason is that pronounced jumps in quantum-size cascades occur when a quasiparticle branch associated with the resonant subband touches zero at some magnetic field. In this case the order parameter, as a rule, immediately drops down to zero, which measures H c, . The spin-magnetic-field interaction splits each quasiparticle branch into two parts, which, as a rule, divides the initial pronounced jump in two smaller ones. Sometimes, the second of them can be situated at a higher magnetic field as compared to the initial jump and, thus, we obtain a higher critical field. In addition, the spin-Zeeman splitting reduces the effect of each subband on the superconducting quantities, which can lead to smoothing of some of the jumps in the cascades. This can also contribute to rare enhancements of H c, when taking into account the spin-magnetic moments. Due to significant orbital motion in a hollow nanocylinder, we found that the paramagnetic limiting field H par is essentially larger than H c, calculated when including orbital and spin-magnetic moments. This can be different for less symmetric confining geometries, e.g., for a rectangular-cross-section nanowire, and deserves further investigations. The paramagnetic limiting field H par was shown to oscillate significantly with R and d due to the quantum-size effects. For the resonant points, H par is significantly enhanced as compared to the bulk Clogston-Chandrasekhar field H CC due to resonant enhancements of the pairing energy. We note that the spin-orbit coupling, as is known [START_REF] Ketterson | Superconductivity[END_REF], can reduce the effects of the paramagnetic pair breaking. However, due to a minor contribution of the paramagnetic mechanism, the role of spin-orbit coupling will not be essential in our system and we leave this subject for future investigations.

Figure 1 .

 1 Figure 1. Hollow aluminum nanocylinder in a parallel magnetic field for R = 0.1 nm at zero temperature (neglecting the spin-magnetic-field interaction): (a) the critical magnetic field H c, as a function of the nanocylinder wall thickness d; (b)-(f) quantum-size cascades in the spatially averaged order parameter ¯ as a function of H for d = 1.7, 2.12, 2.78, 3.06 and 3.81 nm, respectively.

  (a) the critical (depairing) magnetic field H c, is shown as a function of d for R = 0.1 nm and d = 1.4-4.0 nm (this choice represents the nanowire-like regime with R/d 1). The results of figure 1(a) can be compared with those of figure 2(b) (H c, as a function of d for R = 1.0 nm and d = 1.4-4.0 nm) and figure 3(a) (H c, versus d for R = 8.0 nm and d = 1.4-4.0 nm, the nanofilm-like regime with R/d 1).

Figure 2 .

 2 Figure 2. The same as figure 1 but for R = 1 nm. Panels (b)-(f) show results for d = 1.8, 2.12, 2.8, 3.07 and 3.77 nm, respectively.

Figure 3 .

 3 Figure 3. The same as in figures 1 and 2 but for R = 8 nm. Panels (b)-(f) represent d = 1.9, 2.04, 2.64, 3.1 and 3.9 nm, respectively.

Figure 4 .

 4 Figure 4. The Gibbs free energy of the superconducting state (at T = 0) of the hollow nanocylinder for R = 0.5 nm and d = 1.36 nm as a function of the parallel magnetic field.

Figure 5 .

 5 Figure 5. Effect of temperature on the quantum-size cascades for R = 0.1 nm and d = 1.44 nm (no spin-magnetic-field interaction): (a) the contour plot of the spatially averaged order parameter ¯ (in millielectronvolt) as a function of H /H c, (with H c, = 24.5 T) and T/T c (with T c /T c,bulk = 11.4); (b) ¯ versus H /H c, for T/T c = 0, 0.3, 0.6 and 0.95; (c) ¯ as dependent on T/T c for H /H c, = 0, 0.6 and 0.73. The solid curve in (a) shows the position of the discontinuous jump in ¯ down to zero (ending at (T * , H * )). The dashed curve in (a) represents the temperature-dependent critical field H c, (T) = 1.17H c, (1 -T/T c ) 1/2 .

  we have ¯ | jump /k B T * ≈ 1.85. A similar result is obtained for R = 1.0 nm and d = 3.07 nm, see figure 7, where panel (a) shows the contour plot of ¯ as dependent on T/T c and H /H c, (here T c /T c,bulk = 2.8 and H c, = 1.99 T); panel (b) represents the data for ¯ as dependent on H /H c, for T/T c = 0, 0.3, 0.6, and 0.95; (c) gives ¯ as a function of T/T c for H /H c, = 0, 0.6, and 0.73. Here ¯ | jump ≈ 0.086 meV and T * = 0.12T c , which results in ¯ | jump /k B T * ≈ 1.81. Smaller jumps in the cascade are

Figure 6 .

 6 Figure 6. The same as in figure 5 but for R = 1 nm and d = 3.07 nm. Panel (b) shows data for T/T c = 0, 0.74, and 0.95 (with T c /T c,bulk = 2.8 K), panel (c) shows results for H /H c, = 0, 0.45, and 0.75 (with H c, = 1.99 T). Here the dashed curve corresponds to H c, (T) = 1.02H c, (1 -T/T c ) 1/2 .

Figure 7 .

 7 Figure 7. The same as in figure5but for R = 8 nm and d = 2.12 nm. The dashed curve in panel (a) represents H c, (T) = 1.43H c, (1 -T/T c ) 1/2 ; the solid line in (a) shows the position of an abrupt decay of ¯ to zero by a cascade of jumps of about 0.02-0.03 meV. Panel (b) shows our results for T/T c = 0, 0.65, and 0.97 (T c /T c,bulk = 1.9 K), panel (c) displays the data for H /H c, = 0, 0.65, and 0.9 (with H c, = 0.234 T).

Figure 8 .

 8 Figure 8. Including the spin-magnetic-field potential energy: (a) H c, as a function of d for R = 0.1 nm; (b)-(f) quantum-size cascades in ¯ as a function of H for d = 1.7, 2.12, 3.06 and 3.81 nm, respectively.

Figure 9 .

 9 Figure 9. The same as figure 8 but for R = 1 nm. Panels (b)-(f) shows results for d = 1.8, 2.12, 1.8, 3.07 and 3.77 nm.

Figure 10 .

 10 Figure 10. Paramagnetic breakdown of superconductivity: (a) and (b) quantum-size oscillations of the paramagnetic limiting field H par as a function of d (the dashed line represents the Clogston-Chandrasekhar bulk field H CC ); (b)-(f) the spatially averaged order parameter ¯ versus H (in the presence of only spin-magnetic-field interaction) for (R, d) = (0.1, 2.12) nm, (0.1, 1.7) nm, (1.0, 2.8) nm, and (1.0, 3.07) nm, respectively.

Figure 11 .

 11 Figure 11. ¯ as a function of H in the presence of surface roughness for R = 0 and D = 2 d = 6.34 nm: triangles are for the ideal case with δD = 0; squares represent the results obtained from equation (31) with δD = 4 Å.

Figure 12 .

 12 Figure 12. Effect of surface roughness on ¯ as a function of H for R = 0 and various mean diameters D = 2 d = 8.4 nm (a), 8.8 nm (b) and 9.7 nm (c). In each panel triangles represent the ideal cylindrical nanowire and squares show results calculated from equation (31) with δD = 2 Å.

Acknowledgments

This work was supported by the Flemish Science Foundation (FWO-VI) and the ESF-AQDJJ network. MDC acknowledges the support of the EU Marie Curie IEF Action (Grant Agreement No. PIEF-GA-2009-235486-ScQSR).


Mihail D Croitoru
email: mikhail.croitoru@uni-bayreuth.de

Manuel Houzet

Alexandre I Buzdin

Field-Direction Dependence of the Upper Critical Field in Organic Superconductors

Keywords: FFLO state, Anisotropy, Critical field

niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

and singlet pairing is mainly limited by the Zeeman effect. Nevertheless, this limitation can be surmounted by allowing a spin-triplet state or a spatially modulated spin-singletstate (limit is upshifted) [15,16]. Whether a spin part of the superconducting order parameter singlet or triplet is currently under discussion because the NMR experiment with (TMTSF) 2 PF 6 salt at T c and under pressure showed the absence of the Knight shift, supporting the triplet scenario of pairing [17], while the 77 Se NMR Knight shift in recent experiment with (TMTSF) 2 ClO 4 reveals a decrease in spin susceptibility χ s consistent with singlet spin pairing [18]. 13 C NMR measurements with κ-(BEDT-TTF) 2 Cu(NCS) 4 evidenced for a Zeeman-driven transition within the superconducting state and stabilization of inhomogeneous, spatially modulated phase [19].

Recently the field-amplitude and field-angle dependence of the superconducting transition temperature T c (H ) of the organic superconductor (TMTSF) 2 ClO 4 in magnetic field applied along the conduction planes have been reported [14]. The authors observed an upturn of the curve of the upper critical field at low temperatures. This behavior has often been discussed in connection with the possibility of the FFLO state formation [20][21][START_REF] Shimahara | The Physics of Organic Superconductors and Conductors[END_REF][START_REF] Buzdin | [END_REF]. Moreover, an unusual inplane anisotropy of H c2 in the high-field regime was observed, which authors interpret as an evidence of the FFLO state formation. Motivated by the experimental findings, in this work we investigate the influence of the modulating superconducting phase on the in-plane anisotropy of H c2 in layered conductors.

We consider here a quasi-two-dimensional superconductor in magnetic field applied parallel to the plane of the most conductive layers (xy-plane). The system single-electron spectrum is approximated by

ξ p = E (2) p + ε(p z ), where E (2) p = p x 2m x + p y
2m y , and the tight-binding approximation is used to describe the electron motion along the z-direction: [24]. The coupling between layers, t , satisfies the condition

1 ε(p z ) = 2t cos(p z d)
T 2 c0 E F t T c0 .
Here T c0 is the critical temperature of the system at H = 0 with neglected interlayer coupling and E F is the Fermi energy. The Eilenberger equation describing such systems acquires the form [START_REF] Kopnin | Theory of Nonequilibrium Superconductivity[END_REF] L 2 f ω (n p , r, p z ) = Δg ω (n p , r, p z ), (1) where g ω and f ω are the quasiclassical Green functions, operator

L 2 = Ω n + 1 2 v F ∇ x,y + t sin(p z d)[e iQr -e -iQr ]
with Ω n ≡ ω nih and Q = πd φ 0 (-H y , H x ) with φ 0 = π c/e. Here h ≡ μ B H is the Zeeman energy, v F is the 2D Fermi velocity and d is the interlayer distance. In the following, we assume that the system is near the transition, therefore g ω (n p , r, p z ) sign(ω n ), and that the mean free path inside the layers are much larger than the corresponding coherence length.

Restricting ourselves first to the case when paramagnetic limit (FFLO state) is absent, performing averaging over the Fermi surface, and using the iterative procedure with respect to f ω (n p , r, p z ) we obtain the extended Lowerence-Doniach equation (in the sequel, MLD equation)

Δ(r) ln T T c0 = πT n L MLD (ω n )Δ(r) (2) with L MLD (ω n ) ≡ 1/(4ω 3 n ){ v 2 F x ∂ 2 x + v 2 Fy ∂ 2 y } - t 2 [1 -cos(2Qr)]/ω 3 n + [ v 2 F x Q 2 x + v 2 Fy Q 2 y ] × t 2 [1 - 7 cos(2Qr)]/4ω 5 n .
In order to obtain the iterative scheme convergent, we have required Δ(r) v F ∇ x,y Δ(r)/2πT c0 , which implies that the characteristic scale of the order parameter variations should be much smaller than ξ 0 = v F /2πT c0 . The MLD equation contains the term proportional to (v F Q) 2 t 2 , which is absent in the Lowerence-Doniach equation. In the case of v F Q T c0 , one may neglect the last term and obtain the standard Lowerence-Doniach equation.

The angle resolved highest magnetic field at which superconductivity can nucleate in a sample is obtained as

H c2 (α, T )| lim = H c2 ( π 2 )| lim sin 2 (α) + m x m y cos 2 (α) , ( 3
)
where α is the angle between the direction of the external field and x-axis.

Here H c2 ( π 2 )| lim=κ I = m x v F T c0 dt φ 0 2π (1 - T T c0 ), when H t πdv F φ 0 (κ I ). The MLD equation in the regime of magnetic fields, tφ 0 /πdv F H √ tT c0 φ 0 / πdv F (κ II ) gives rise to H c2 ( π 2 )| κ II = 7ζ(3)t 2 8π 3 dT c0 m x T -T ct φ 0 . In the regime κ III of the magnetic fields, φ 0 πdv F √ tT c0 H φ 0 πdv F T c0
, the MLD equation describes the beginning of the reentrant superconductivity regime [START_REF] Lebed | [END_REF]27] with the upper critical field defined as

H c2 π 2 κ III = 28 31
ζ(3) ζ(5) T ct √ T -T ct dt √ 2m x φ 0 . ( 4
)
Here, one can see that T > T ct , where T ct is introduced via the relation ln(

T c0 /T ct ) = 7ζ(3) 4π 2 t 2 T 2 c0 .
The previous treatment can be extended to extremely large magnetic fields, provided that one can omit the assumption that the characteristic scale of the order parameter variations is much smaller than ξ 0 . For this purpose, we need to reconsider the solution of the Eilenberger equation. According to the 2D Eilenberger equation, the magnetic field induced potential has the form V (r) = 2it sin(p z d) sin(Qr), i.e., it is periodic in real space. Since the system under consideration is in an arbitrary magnetic field, the solution of Eq. ( 1) can be written without any loss of generality as a Bloch function

f ω (n p , r, p z ) = e iqr m e imQr f m (ω n , n p , p z ), (5)
or f ω (n p , r, p z ) = e iqr f ω (n p , r, p z ), where f ω (n p , r, p z ) has the periodicity of 2π/Q. In Eq. ( 5), we take into account the possibility for the formation of the pairing state (k+ q 2 ↑; -k + q 2 ↓) with the modulation vector q ∼ 2μ B H/ v F . Because of the form for f ω (n p , r, p z ) of Eq. ( 5), one can write Δ(r) as

Δ(r) = e iqr m e i2mQr Δ 2m . ( 6
)
From symmetry considerations, it follows that Δ -2m = Δ 2m . Substituting Eqs. ( 5) and ( 6) into Eq. ( 1) and retaining in the sum only the zero and the first term, since we want to investigate the regime H T c0 πdv F φ 0 (T c0 v F Q), when q = 0 one gets after averaging over the Fermi surface ln

T c0 T = πT n t 2 ω 2 n ω 2 n + v F Q 2 2 -1/2 . ( 7
)
For extremely large values of the external magnetic field, we can use ω n ∼ T c0 v F Q and obtain Eq. ( 3) with

H c2 ( π 2 )| κ I V = t 2 2dv F φ 0
(T c0 -T ) . Therefore, summarizing the above discussion, we can conclude that in the non-modulated phase the in-plane anisotropy of the upper critical magnetic field can be defined by Eq. ( 3) and H c2 ( π2 )| lim should be taken accordingly, as described in Table 1.

As it is known [START_REF] Buzdin | [END_REF][START_REF] Saint-James | Type II Superconductivity[END_REF], the FFLO state only appears at

T < T * 0.56T c0 or H > H * 1.06T c0 /μ B . Therefore, v F Q T c0 or √ tT c0 (v F Q).
This allows to retain only the zero and first term in Eq. ( 5). Expanding up to the second order in t/T c0 and making use of the average scheme for the Fermi surface as, . . .

≡ π/d -π/d d dp z 2π
2π 0 dϕ(. . .), one obtains [START_REF] Kopnin | Theory of Nonequilibrium Superconductivity[END_REF] ln

T c0 T = πT n 1 ω n - 1 L n (q) + t 2 2 1 L n (q + Q)L 2 n (q) + 1 L n (q -Q)L 2 n (q) , (8)
Table 1 The in-plane upper critical magnetic field H c2 (π/2) Regime Limits of validity

H c2 ( π 2 )| lim I H φ 0 πdv F t mx v F T c0 dt φ 0 2π (1 -T T c0 ) II φ 0 πdv F t H φ 0 πdv F √ tT c0 7ζ(3)t 2 8π 3 dT c0 mx T -Tct φ 0 III φ 0 πdv F √ tT c0 H φ 0 πdv F T c0 28 31 ζ(3) ζ(5) Tct √ T -Tct dt √ 2mx φ 0 IV H φ 0 πdv F T c0 t 2 2dv F φ 0 (T c0 -T )
where the following notation is introduced:

L n (q) = ω n + ih sign(ω n ) + i v F q
2 . We assume that the optimum direction of the FFLO modulation vector is fixed by the symmetry of the Fermi surface and is solely determined in the paramagnetic limit [START_REF] Buzdin | [END_REF]30]. The magnitude of the FFLO modulation vector is defined according to [31]. Having in our disposal the absolute value and the direction of the FFLO modulation vector, we can calculate the contribution of the orbital effect. Introducing the temperature of the Pauli limited case by ln( T c0

T cP ) = πT cP n { 1 ω n -1 L n (q) }, we can write the equation for the onset of the superconductivity as ln

T cP T = AπT n T n (ω, q, +Q) + T n (ω, q, -Q), (9)
where

T n (ω, q, ±Q) = t 2 2 2π 0 dϕ 2π L -1 n (q ± Q)L -2 n (q)
, and A is a function taken at T cP (H ) and it will be displayed elsewhere [32]. This equation determines H -T phase diagram for layered superconductors in the presence of the paramagnetic and the orbital effects in the limits: tT c0 (v F Q) 2 and t T c0 .

The summation over the Matsubara frequencies of the orbital terms were performed numerically. In the sum over the Matsubara frequencies, N = 10 4 terms were used. This number suffices for the convergence. The results presented in this section are in the dimensionless units. Let us start with a far-fetched situation and investigate the anisotropy of the upper critical field in the absence of the FFLO modulation vector for H > H * . The correction to the pure Pauli limited (PL) H -T phase diagram due to the orbital effects is then described by ln

T cP T = AπT t 2 n ω n (a 2 + ω 2 n -3h 2 ) (ω 2 n + h 2 )z , with z ≡ a 4 + 2a 2 (ω 2 n -h 2 ) + (ω 2 n + h 2 ) 2 and a = v F Q/2.
At low temperature, the sum in r.h.s can become negative. However, coefficient A is also negative, meaning that T cP > T . This is seen in the polar plot in the upper panel in Fig. 1, which illustrates the magnetic field angular dependence of the normalized onset temperature in the non-modulated superconducting phase, T c /T cP , taken at H/H c0 = 0.96 and H/H c0 = 0.9, where H c0 is the critical magnetic field at T = 0 in Pauli limited 2Dsuperconductors. Calculations were performed for the Fermi Fig. 1 The in-plane field-direction dependence of T c /T cP at H/H c0 = 0.96 (red line) and H/H c0 = 0.9 (green line) for m x = 5m y , q = 0 (upper panel), and q = 0 (lower panel), respectively. The Fermi velocity parameter η = 3.45, and the q-vector is along the x-axis (Color figure online) velocity v F = 1.25 × 10 7 cm/s, or introducing dimensionless velocity parameter η = v F πd/φ 0 μ B , for η = 3.45, and for the in-plane mass-anisotropic situation, when m x = 5m y . Hereinafter, t/T c0 = 0.2. In the polar plot, the direction of each point seen from the origin corresponds to the field direction and the distance from the origin illustrates to the normalized critical temperature. In the case of uniform Fig. 2 The in-plane field-direction dependence of T c /T cP at H/H c0 = 0.96 (red line) and H/H c0 = 0.90 (green line) for m x = 5m y and q = 0. The Fermi velocity parameter η is 0.31 (upper panel) and 0.63 (lower panel) (Color figure online) (non-modulated) pairing, a small dip is observed at angles ±90 • from the x-axis. The developments of these dips is due to the ellipsoidal Fermi line, corresponding to an ellipsoidal electronic dispersion relation. Concluding this part, one can infer that when magnetic field is rotated in the xy-plane, the upper critical field is smallest at angles ±90 • , i.e., when it is parallel to the direction of the lightest electron mass.

Hitherto the possibility of the formation of the FFLO state was not considered in our numerical calculations. When abating this constraint, the in-plane anisotropy of H c2 changes essentially. The lower panel in Fig. 1 shows the polar plot of the magnetic field direction dependence of the normalized superconducting transition temperature, taken at H/H c0 = 0.96 and T /T c0 = 0.9, for Fermi velocity parameter η = 3.45 and mass anisotropy m x = 5m y . Hereinafter the q vector is along the x-axis. Comparing this result with that in the upper panel, we see that letting a pairing with the finite momentum in the layered system makes ΔT cP strongly anisotropic. Namely, the position of the maximum of the onset temperature is rotated from α = 0 • , and 180 • typical for q = 0, to α = ±70 • , ±110 • for q = 0 and H/H c0 = 0.9. The minimum in T c (α) is rotated from α = ±90 • for q = 0 to α = ±50 • , ±130 • for q = 0 and H/H c0 = 0.9.

Let us now investigate the influence of the absolute value of the Fermi velocity on the in-plane anisotropy of the upper critical field in layered superconductors. The upper and lower panels of Fig. 2 display T c (α) for the Fermi velocity parameter η = 0.31 and η = 0.63, respectively. One can immediately infer that absolute value of the Fermi velocity parameter essentially influences the field-direction dependence of the superconducting onset temperature. It is seen that the position of the maximum value of T c (α) in the case of the smallest considered Fermi velocity parameter is shifted back towards α = ±90 • . The same tendency is observed for η > 4 (not shown here).

In conclusion, in this work we have analyzed the fieldamplitude and the field-direction dependence of the onset temperature of the superconductivity in layered conductors. In such superconductors, the interplay of the system dimensionality, orbital pair-breaking, Pauli paramagnetism, and a possible formation of the FFLO state can lead to an anomalous in-plane anisotropy of the upper critical field.

Acknowledgements

We acknowledge the support by the European Community under a Marie Curie IEF Action (Grant Agreement No. PIEF-GA-2009-235486-ScQSR), the French Project "SINUS" ANR-09-BLAN-0146, and thank A. S. Mel'nikov for fruitful discussions.