Strong Landau-quantization effects in high-magnetic-field superconductivity of a two-dimensional multiple-band metal near the Lifshitz transition

TL;DR

This paper explores how Landau quantization causes giant oscillations and reentrant superconductivity in a two-dimensional multiband superconductor near a Lifshitz transition under high magnetic fields.

Contribution

It demonstrates the significant impact of Landau quantization on superconductivity near a Lifshitz transition, revealing reentrant behavior and oscillations in the transition temperature.

Findings

Giant oscillations in transition temperature T_C2(H) due to Landau levels.

Reentrant superconducting regions emerge at specific magnetic fields.

Reentrance is suppressed by Zeeman spin splitting but can be restored under resonance conditions.

Abstract

We investigate the onset of superconductivity in magnetic field for a clean two-dimensional multiple-band superconductor in the vicinity of the Lifshitz transition when one of the bands is very shallow. Due to small number of carriers in this band, the quasiclassical Werthamer-Helfand approximation breaks down and Landau quantization has to be taken into account. We found that the transition temperature TC2(H) has giant oscillations and is resonantly enhanced at the magnetic fields corresponding to full occupancy of the Landau levels in the shallow band. This enhancement is especially pronounced for the lowest Landau level. As a consequence, the reentrant superconducting regions in the temperature-field phase diagram emerge at low temperatures near the magnetic fields at which the chemical potential matches the Landau levels. The specific behavior depends on the relative strength of the intraband and interband pairing interactions and the reentrance is most pronounced in the purely interband coupling scenario. The reentrant behavior is suppressed by the Zeeman spin splitting in the shallow band, the separated regions disappear already for very small spin-splitting factors. On the other hand, the reentrance is restored in the resonance cases when the spin-splitting energy exactly matches the separation between the Landau levels. The predicted behavior may realize in the gate-tuned FeSe monolayer.

Equations

$$\mathcal{H}=$$

$$-$$

$$\mathcal{H}\approx {\mathcal{H}}_{\mathrm{HF}}$$

$$\displaystyle-\sum_{\alpha\beta}\bar{\Delta}^{\alpha}(\bm{r})U^{-1}_{\alpha\beta}\Delta^{\beta}(\bm{r})\Big{\}},$$

$${\hat{\Sigma }}^{\alpha }(\mathbf{r})=\left[\begin{array}{cc}0& {\Delta }^{\alpha }(\mathbf{r})\\ {\bar{\Delta }}^{\alpha }(\mathbf{r})& 0\end{array}\right]$$

$${\mathcal{G}}^{\alpha }(\mathbf{r},{\mathbf{r}}^{\prime };\tau )=-⟨T_{\tau }[{\psi }_{\alpha }(\mathbf{r},\tau ){\psi }_{\alpha ,}^{†}(\mathbf{r},0)]⟩,$$

$$[\mathrm{i}{\omega }_n+{\hat{\xi }}_{\mathbf{k}}^{\alpha }\widehat{{\tau }^z}-{\mu }_{z}H-{\hat{\Sigma }}^{\alpha }(\mathbf{r})]{\mathcal{G}}_{{\omega }_n}^{\alpha }(\mathbf{r},{\mathbf{r}}^{\prime })={\hat{\tau }}^0\delta \left(\mathbf{r}-{\mathbf{r}}^{\prime }\right),$$

$${\Delta }^{\alpha }(\mathbf{r})=-T\sum _{{\omega }_n=-\infty }^{\infty }\sum _{\beta }{U}_{\alpha \beta }{F}_{{\omega }_n}^{\beta }(\mathbf{r}),$$

$$F_{{\omega }_n}^{\alpha }(\mathbf{r})\approx \int d{\mathbf{r}}^{\prime }{K}_{\alpha }(\mathbf{r},{\mathbf{r}}^{\prime };{\omega }_n){\Delta }^{\alpha }({\mathbf{r}}^{\prime })$$

$$K_{\alpha }(\mathbf{r},{\mathbf{r}}^{\prime };{\omega }_n)=-G_{0,{\omega }_n,+}^{\alpha }(\mathbf{r},{\mathbf{r}}^{\prime })G_{0,-{\omega }_n,-}^{\alpha }(\mathbf{r},{\mathbf{r}}^{\prime }),$$

$${\hat{\Lambda }}^{-1}\left[\begin{array}{c}{\Delta }_0^h\\ {\Delta }_0^e\end{array}\right]\approx \sum _{{\omega }_n>0}^{\Omega }\frac{2\pi T}{{\omega }_n}\left[\begin{array}{c}\left[\frac{1}{2}+{\eta }_h({\omega }_n)\right]{\Delta }_0^h\\ {\Delta }_0^e\end{array}\right],$$

$${\Lambda }_{0,e}^{-1}$$

$${\Lambda }_{0,h}^{-1}$$

$${{\rm Y}}_C\equiv \sum _{{\omega }_n>0}\frac{2\pi T_C}{{\omega }_n}{\eta }_h({\omega }_n)=\frac{2}{\pi }\sum _{n=0}^{\infty }\frac{{\tan }^{-1}\left(\frac{{\mu }_h/T_C}{\pi (2n+1)}\right)}{2n+1},$$

$$\left[{\hat{\Lambda }}^{-1}-{\hat{\Lambda }}_0^{-1}\right]\left[\begin{array}{c}{\Delta }_0^h\\ {\Delta }_0^e\end{array}\right]=0,$$

$$\det (\hat{W})=W_{11}{W}_{22}-W_{12}{W}_{21}=0.$$

$${\Lambda }_{0,e}^{-1}=$$

$$+$$

$$G_{0,{\omega }_n,\pm }^{\alpha }(\mathbf{r},{\mathbf{r}}^{\prime })={\mathrm{e}}^{2\mathrm{i}{\varphi }_{\mathbf{A}}(\mathbf{r},{\mathbf{r}}^{\prime })}{g}_{0,\pm }^{\alpha }(|\mathbf{r}-{\mathbf{r}}^{\prime }|,{\omega }_n)$$

$$K_{\alpha }(\mathbf{r},{\mathbf{r}}^{\prime };{\omega }_n)=-\exp \left(-\mathrm{i}\frac{{\left[\mathbf{r}\times {\mathbf{r}}^{\prime }\right]}_z}{l^2}\right)g_{0,+}^{\alpha }(\rho ,{\omega }_n)g_{0,-}^{\alpha }(\rho ,-{\omega }_n)$$

$${\Delta }^{\alpha }(\mathbf{r})={\Delta }_0^{\alpha }\exp \left(-\frac{r^2}{2l^2}\right).$$

$$F_{{\omega }_n}^{\alpha }(\mathbf{r})=-{\Delta }^{\alpha }(\mathbf{r})$$

$$\times \int {\Phi }_{\mathbf{r}}(\rho )\rho \mathrm{d}\rho g_{0,+}^{\alpha }(\rho ,{\omega }_n)g_{0,-}^{\alpha }(\rho ,-{\omega }_n){\mathrm{e}}^{-\frac{{\rho }^2}{2l^2}}$$

$$\displaystyle\Phi_{\bm{r}}(\rho)=\int^{2\pi}_{0}\mathrm{d}\theta\exp\Big{[}-\frac{r\rho}{l^{2}}\mathrm{e}^{\mathrm{i}\theta}\Big{]}=2\pi.$$

$$F_{{\omega }_n}^{\alpha }(\mathbf{r})=-\pi N_{\alpha }{\lambda }_{{\omega }_n}^{\alpha }{\Delta }^{\alpha }(\mathbf{r}),$$

$${\lambda }_{{\omega }_n}^{\alpha }={\int }_0^{\infty }\frac{2\rho \mathrm{d}\rho }{N_{\alpha }}{g}_{0,+}^{\alpha }(\rho ,{\omega }_n)g_{0,-}^{\alpha }(\rho ,-{\omega }_n){\mathrm{e}}^{-\frac{{\rho }^2}{2l^2}}.$$

$$g_0^e(\rho ,{\omega }_n)=\int \frac{d\mathbf{k}}{(2\pi {)}^2}\frac{\exp (\mathrm{i}\mathbf{k}\mathbf{\rho })}{\mathrm{i}{\omega }_n+{\xi }_{\mathbf{k}}^e}.$$

$$g_0^e(\rho ,{\omega }_n)g_0^e(\rho ,-{\omega }_n)=-\int \frac{d\mathbf{k}}{(2\pi {)}^2}\frac{d{\mathbf{k}}^{\prime }}{(2\pi {)}^2}\frac{\exp [\mathrm{i}\left(\mathbf{k}-{\mathbf{k}}^{\prime }\right)\mathbf{\rho }]}{\left(\mathrm{i}{\omega }_n-{\xi }_{\mathbf{k}}^e\right)\left(\mathrm{i}{\omega }_n+{\xi }_{{\mathbf{k}}^{\prime }}^e\right)}.$$

$$g_0^e(\rho ,{\omega }_n)g_0^e(\rho ,-{\omega }_n)=\pi N_e\int \frac{d\mathbf{q}}{(2\pi {)}^2}{⟨\frac{\exp (\mathrm{i}\mathbf{q\rho })}{|{\omega }_n|+\mathrm{i}{\mathbf{v}}_e\mathbf{q}/2}⟩}_e$$

$${\lambda }_{{\omega }_n}^e=4\pi \int \frac{d\mathbf{q}}{(2\pi {)}^2}{\int }_0^{\infty }\rho \mathrm{d}\rho {⟨\frac{\exp (\mathrm{i}\mathbf{q\rho }-\frac{{\rho }^2}{2l^2})}{2|{\omega }_n|+\mathrm{i}{\mathbf{v}}_e\mathbf{q}}⟩}_e.$$

$${\lambda }_{{\omega }_n}^e=4\pi {\int }_0^{\infty }\mathrm{d}s{\int }_0^{\infty }\rho \mathrm{d}\rho$$

$$\times \int \frac{d\mathbf{q}}{(2\pi {)}^2}{⟨\exp \left[-s\left(2|{\omega }_n|+\mathrm{i}{\mathbf{v}}_e\mathbf{q}\right)+\mathrm{i}\mathbf{q\rho }-\frac{{\rho }^2}{2l^2}\right]⟩}_e$$

$$=4\pi \underset{0}{\overset{\infty }{\int }}\mathrm{d}s\underset{0}{\overset{\infty }{\int }}\rho \mathrm{d}\rho \delta \left(\mathbf{\rho }-s{\mathbf{v}}_e\right){⟨\exp \left[-2s|{\omega }_n|-\frac{{\rho }^2}{2l^2}\right]⟩}_e.$$

$${\lambda }_{{\omega }_n}^e=2{\int }_0^{\infty }\mathrm{d}s{\mathrm{e}}^{-2s|{\omega }_n|}{⟨\exp \left(-\frac{v_e^2{s}^2}{2l^2}\right)⟩}_e.$$

Full text

I Introduction

The rich field of multiple-band superconductors has been reincarnated by the unexpected discovery of superconductivity in the magnesium diboride at 40 K Nagamatsu et al. (2001); *Liu:PhysRevLett2001; *Xi:RPP.2008; *Budko:PhysC.2015 and received a further powerful boost from the discovery of several families of iron-based superconductors (FeSCs), see, e. g., reviews Paglione and Greene (2010); *Hirschfeld:RPP.2011; *Chubukov:AnnCMP3.2012; *Hosono:PhysC514.2015. The Fermi surfaces of these materials are composed of several disconnected nonequivalents parts. These parts not only have different electronic properties but, also, in the superconducting state they may have different gaps causing many peculiar properties of these materials.

In contrast to the magnesium diboride, FeSCs are semimetals: their band sizes are rather small with the typical Fermi energies $\epsilon_{F}\lesssim$ 50 meV and Fermi velocities $v_{F}\lesssim$ 107 cm/s. In addition, the band Fermi energies can be shifted by doping or pressure. As a result, several FeSC compounds can be driven through Lifshitz transitions at which top or bottom of one band crosses the Fermi level and the corresponding Fermi pocket vanishes. Examples include Ba1-xKxFe2As2 near $x\approx 0.8$ Xu et al. (2013) and LiFe1-xCoxAs for $x\lesssim 0.1$ Miao et al. (2015). In the first case the electron band at the $M$ point is shifted above the Fermi level and in the second case one of three hole bands at the $\Gamma$-point sinks below the Fermi level.

A special case is realized in the simplest compound FeSe. Discovery of superconductivity in the FeSe single layer grown on SrTiO3 substrate with record high transition temperature for FeSCs, $T_{C}\!\gtrsim\!55$ K, has been a major breakthrough in the field Wang et al. (2012); Liu et al. (2012). 111The in-situ four-terminal transport measurement give the transition temperature above $100$ K Ge et al. (2015). This claim, however, was not yet confirmed by any other technique. The bulk material has the tetragonal-to-orthorhombic transition at 87 K which is followed by the superconducting transition at 9 K. Its Fermi surface is composed of one hole pocket and two electron pockets which have very small sizes with $\epsilon_{F}\!\sim\!10-20$ meV Mukherjee et al. (2015). The electronic and superconducting properties of the tetragonal FeSe monolayer on SrTiO3 are very different from the bulk material. The optimally-doped state has only electron bands and the hole band is sinked $\sim 80$ meV below the Fermi level Liu et al. (2012) meaning that the single layer is strongly electron-doped with respect to the bulk crystal. This doping is probably caused by oxygen-vacancies diffusion in SrTiO3 during annealing. Such difference implies that at the intermediate electron doping level FeSe goes through the Lifshitz transition at which the hole band at the $\Gamma$-point is depleted. This transition has been indeed observed in K-dosed FeSe thin films Wen et al. (2016). Such electronic structure is also realized in the intercalated compound $({\mathrm{Li}}{\mathrm{Fe}})\mathrm{OH}\mathrm{Fe}\mathrm{Se}$ with $T_{C}\!=\!40$ K Niu et al. (2015). Transport measurements have been done on the monolayer protected by the FeTe capping layers Sun et al. (2014); Zhang et al. (2014); Zhou et al. (2016), which reduces the transition temperature down to $\sim 23-25$K. The upper critical field of such system has been found to be around 50 Tesla. In a controlled way, the FeSe monolayer can be doped using $K$ coating Shi et al. (2016). It was shown that such coating causes the second Lifshitz transition at which the electron band emerges at the $\Gamma$-point which promotes strong enhancement of $T_{C}$. Also, it was found that the gating of small-size FeSe crystals induces the surface superconductivity with $T_{C}\!=\!48$K Lei et al. (2016). It is likely that in this case the surface region acquires the band structure similar to the FeSe monolayer.

The ubiquity of shallow bands and Lifshitz transitions in FeSCs motivated several recent theoretical studies devoted to the influence of such bands on superconducting pairing Bang (2014); Chen et al. (2015); Chubukov et al. (2016), see also related general considerations Innocenti et al. (2010); Valentinis et al. (2016). One can distinguish two basic scenarios Chen et al. (2015): (i) The shallow band is essential for superconductivity. In this case the superconducting state vanishes when this band is depleted. (ii) The Cooper pairing is dominated by deep bands and superconducting gap is induced into the shallow band via pair-hopping interactions. In this case the superconducting temperature changes only weakly at the Lifshitz transition. It was also demonstrated in Ref. Koshelev and Matveev (2014) that in the case of the second scenario the superconductivity actually smears the Lifshitz transition in the thermodynamics sense but, nevertheless, the density of states changes qualitatively when the shallow band is depleted.

In this paper we investigate the influence of shallow bands on the onset of superconductivity in the magnetic field. The upper critical field, $H_{C2}$, is one of the key characteristics of type-II superconductors. In most materials superconductivity is destroyed by the orbital effect of the magnetic field. In the case of weak impurity scattering the orbital upper critical field, $H_{C2}\equiv H_{C2}^{\mathrm{O}}$, scales inversely proportional to Fermi velocity squared, $H_{C2}^{\mathrm{O}}\propto v_{F}^{-2}$  Helfand and Werthamer (1966); *Werthamer:PRev147.1966; Kogan and Prozorov (2012), meaning that the orbital effect diminishes with decreasing the band size. The temperature dependences of $H_{C2}$ and its anisotropy may be strongly influenced by multiple-band structure Dahm and Schopohl (2003); *Gurevich:PhysRevB.2003; *Golubov:PRB.2003; Kogan and Prozorov (2012).

The superconductivity is also destroyed by the Zeeman spin splitting induced by the magnetic field. Without the orbital effect in a single-band material the superconducting state is destroyed when the spin-split energy $\mu_{z}H$ exceeds $\Delta/\sqrt{2}$, where $\mu_{z}$ is magnetic moment and $\Delta$ is the energy gap. This gives the paramagnetic limit, $H_{P}=\Delta/\sqrt{2}\mu_{z}$. In most materials the spin-splitting effects are weak in comparison with the orbital ones, $H_{C2}^{\mathrm{O}}\ll H_{P}$. In the case when both orbital and spin effects are present, the relative contribution of the spin splitting is usually characterized by the Maki parameter $\alpha_{M}=\sqrt{2}H_{C2}^{\mathrm{O}}/H_{P}$, which in clean single-band materials can be evaluated as $\alpha_{M}=\pi^{2}\Delta/4\epsilon_{F}$, where $\epsilon_{F}$ is the Fermi energy. This means that the role of spin effects enhances in small Fermi surfaces. The spin splitting also dominates in two-dimensional and layered materials when the magnetic field is applied along conducting planes.

The standard theory of $H_{C2}$ is based on the quasiclassical approximation which neglects the Landau quantization Helfand and Werthamer (1966); *Werthamer:PRev147.1966; Kogan and Prozorov (2012). This theory works with very high accuracy for overwhelming majority of superconductors because at $H\!\sim\!H_{C2}$ the cyclotron frequency $\omega_{c}$ is typically much smaller than $\epsilon_{F}$. Nevertheless, the effects of Landau quantization on the behavior of $H_{C2}$ and related superconducting properties in single-band materials were first studied in the seminal papers Rajagopal and Vasudevan (1966); Gunther and Gruenberg (1966); *Gruenberg:PRev176.1968 and later have been worked out in great detail Tešanović et al. (1989); *Tesanovic:PRB43.1991; Rieck et al. (1990); Maniv et al. (1992); MacDonald et al. (1992); Norman et al. (1992a, b); Vavilov and Mineev (1997); Mineev (2000); *Champel:PhilMagB81.2001; Zhuravlev and Maniv (2012), see also reviews Rasolt and Tešanović (1992); Maniv et al. (2001). It was predicted that in clean materials the quantization may dramatically influence the low-temperature behavior of the upper critical field. The density of states is sharply enhanced when the chemical potential crosses the Landau levels, $\mu\!=\!\omega_{c}(\ell\!+\!1/2)$, at the magnetic fields $H=H_{\ell}$. This enhancement is beneficial for superconductivity. It was actually demonstrated that in an ideally clean single-band superconductor without Zeeman spin splitting the transition temperature is always finite at $H=H_{\ell}$ Gruenberg and Gunther (1968). Such resonant enhancements of the transition temperature are especially pronounced in two-dimensional case MacDonald et al. (1992); *Norman:PhysicaC196.1992; Norman et al. (1992b); Maniv et al. (1992). This would mean that in conventional clean materials superconductivity should persist up to fields much higher then the quasiclassical orbital upper critical field. Moreover, in the extreme quantum limit the local maximums of transition temperature were predicted to increase with the magnetic field Gunther and Gruenberg (1966); *Gruenberg:PRev176.1968; Tešanović et al. (1989); *Tesanovic:PRB43.1991; Rieck et al. (1990). In most superconducting materials, however, this limit requires magnetic fields above 100T, which is beyond practical accessibility. In addition, this ultra-high-field reentrant superconductivity is easily destroyed by impurity scattering and Zeeman spin splitting Gunther and Gruenberg (1966); *Gruenberg:PRev176.1968; Rieck et al. (1990), unless the spin-splitting energy exactly matches the Landau-level spacing Maniv et al. (1992); Norman et al. (1992b); Maniv et al. (2001). On the other hand, near the accessible quasiclassical $H_{C2}$ the Landau-level indices are large and quantization effects are weak. As a consequence, in superconductors with large Fermi surfaces one can expect only very weak quantum oscillations of the temperature or angle dependence of $H_{C2}$ noticeable in extremely clean materials at very low temperatures.

A direct consequence of small electronic bands in FeSCs is very high upper critical fields in these materials, ranging from 15 to 100 tesla for different compounds and dopings Gurevich (2011). For compounds near the Lifshitz transition, the orbital effect is the weakest for the shallow band. Therefore one can expect that this band strongly influences the upper critical field. In contrast to single-band materials, the cyclotron frequency may be comparable with the Fermi energy of a small-size band at the upper critical field meaning that only few Landau levels may be occupied. In this case the quasiclassical approximation breaks down and the Landau quantization is essential. Furthermore, the spin-splitting effects are more pronounced in the shallow band and typically can not be neglected. The role of spin-splitting effects in multiple-band superconductors within quasiclassical approximation has been recently investigated in Ref. Gurevich (2010).

In this paper we investigate the upper critical field in a clean two-dimensional two-band superconductor in the vicinity of the Lifshitz transition when one of the bands is very shallow 222This model can be straightforwardly generalized to the case of several identical deep bands. For this band we take into account the Landau quantization precisely, while for the deep band we use the standard quasiclassical approximation. We will demonstrate that in such system the transition temperature $T_{C2}(H)$ has giant oscillations and is resonantly enhanced at the magnetic fields $H_{\ell}$ corresponding to the crossing of the $\ell$’s shallow-band Landau level and the chemical potential. This enhancement is most pronounced for the lowest Landau level, $\ell=0$, and rapidly decreases with the increasing Landau-level index. We mostly focus on the case when the highest field $H_{0}$ is close to the quasiclassical upper critical filed $H_{C2}^{\mathrm{qc}}$. In the case $H_{0}>H_{C2}^{\mathrm{qc}}$, the temperature-field phase diagram may acquire the reentrant superconducting region located at low temperatures around $H\sim H_{0}$, that is disconnected from the main low-field superconducting region. This reentrant piece merges with the main part when the pocket size diminishes. The specific behavior depends on the relative strength of the intraband and interband coupling constants and it is most pronounced when the interband coupling dominates. The Zeeman spin splitting strongly reduces the sizes of the reentrant regions and changes their location in the parameter space. However, the reentrant superconductivity reappears in resonance conditions, when the spin splitting energy $2\mu_{z}H$ exactly matches the separation between the Landau levels Maniv et al. (1992); Norman et al. (1992b).

The paper is organized as follows. In Secs. II and III, we describe the model two-band Hamiltonian and the corresponding Gor’kov equations. In Sec. IV, we discuss the transition temperature in zero magnetic field in the presence of a shallow band. In Sec. V we derive equations which determine the superconducting instability in the magnetic field. This instability is mostly determined by the field and temperature dependences of the pairing kernels. The behavior of the shallow-band kernel, strongly influenced by the Landau quantization, is discussed in Sec. V.3. In Sec. V.4, we investigate the influence of the Zeeman spin splitting on this kernel. The numerically-computed temperature-magnetic field phase diagrams are presented in Sec. VI. We also discuss in this section the dependence of the high-field transition temperatures on the strength of interband coupling. We conclude the paper in Sec. VII.

II Two-band model

To investigate the shallow-band effects in superconductors, we consider the simplest two-band model

[TABLE]

where $s\!=\uparrow,\downarrow=\!+\!,\!-$ is the spin index and $\alpha\!=\!e$ ($h$) represents the $e$-band ($h$-band) with the energy dispersion $\varepsilon_{\bm{k}+\bm{Q}}^{e}=\frac{\bm{k}^{2}}{2m_{e}}$ ($\varepsilon_{\bm{k}}^{h}=-\frac{\bm{k}^{2}}{2m_{h}}+\varepsilon_{0}$) with $m_{e}$ and $m_{h}$ being the band masses. In the $e$-band dispersion the momentum is measured from the nesting wave vector $\bm{Q}$. In the real-space operator $\varepsilon_{\bm{k}}^{\alpha}$ the wave vector $\bm{k}$ has to be replaced by the gauge-invariant gradient operator $\bm{k}\rightarrow-\mathrm{i}\bm{\nabla}_{\mathbf{r}}-\frac{e}{c}\bm{A}$. 333We use the natural system of units in which the Planck constant $\hbar$ and the Boltzmann constant $k_{B}$ are set to unity. The second term in the first line describes the Zeeman spin splitting, $\sigma^{z}=\text{diag}[1,-1]$ in the spin space, and, for simplicity, we have set the magnetic moments for both $e$- and $h$-band electrons to be $\mu_{z}$.

In normal state this model has two Lifshitz transition points at the chemical potential $\mu=0$ and $\mu=\varepsilon_{0}$ (see Fig. 1). For definiteness, we consider the system in the vicinity of $\mu=\varepsilon_{0}$ transition, i. e., we assume that the hole band is shallow and the electron band is deep. An equivalent model can also be used for description of the system with several identical deep bands.

III The Gor’kov equations

To tackle the many-body Hamiltonian in Eq. (1), we use the mean-field method to approximate the many-body quantum states as the Hartree-Fock states. In this approximation, the Hamiltonian becomes a one-body operator

[TABLE]

where we introduced the Nambu vector $\psi_{\alpha}^{T}(\bm{r})=\left[c_{\alpha,\downarrow}(\bm{r}),c_{\alpha,\uparrow}^{\dagger}(\bm{r})\right]$, $\hat{\tau}^{z}=\text{diag}[1,-1]$, $\bm{k}\rightarrow-\mathrm{i}\bm{\nabla}_{\mathbf{r}}-\frac{e}{c}\hat{\tau}^{z}\bm{A}$ in the operator $\hat{\varepsilon}_{\bm{k}}^{\alpha}$,

[TABLE]

with the gap parameters $\Delta^{\alpha}(\bm{r})\!=\!\sum_{\beta}U_{\alpha\beta}\langle c_{\beta,\downarrow}(\bm{r})c_{\beta,\uparrow}(\bm{r})\rangle$ and $\bar{\Delta}^{\alpha}(\bm{r})$ is its complex conjugate. The $2\times 2$ imaginary-time ($\tau$) Green’s function is defined as

[TABLE]

where $\psi_{\alpha}(\bm{r,}\tau)\!=\!\mathrm{e}^{-\tau(\mathcal{H}_{\text{HF}}-\mu\mathcal{N})}\psi_{\alpha}(\bm{r})\mathrm{e}^{\tau(\mathcal{H}_{\text{HF}}-\mu\mathcal{N})}$ with the chemical potential $\mu$ and the total-number operator $\mathcal{N}=\sum_{\alpha,\sigma}\int d^{2}\bm{r}c_{\alpha,\sigma}^{\dagger}(\bm{r})c_{\alpha,\sigma}(\bm{r})$. The Green’s function satisfies the matrix Gor’kov equation in the frequency representation,

[TABLE]

where $\omega_{n}=(2n+1)\pi T$ is the Matsubara frequency, $\hat{\tau}^{0}$ is the $2\times 2$ identity matrix, and $\hat{\xi}_{\bm{k}}^{\alpha}=\hat{\varepsilon}_{\bm{k}}^{\alpha}-\mu$. The gap parameter is expressed in terms of the anomalous Green’s function $F_{\omega_{n}}^{\alpha}(\bm{r},\bm{r}^{\prime})\equiv[\mathcal{G}_{\omega_{n}}^{\alpha}(\bm{r},\bm{r}^{\prime})]_{12}$ as,

[TABLE]

where $F_{\omega_{n}}^{\beta}(\bm{r})\equiv F_{\omega_{n}}^{\beta}(\bm{r},\bm{r})$. To analyze the behavior of the superconducting gap, one have to solve Eqs. (5) and (6) self-consistently.

In this paper, we are only interested in the region near the upper critical field ($H_{C2}$), where $\Delta^{\alpha}(\bm{r})\rightarrow 0$ and it is sufficient to keep only the anomalous Green’s function $F_{\omega_{n}}^{\alpha}(\bm{r},\bm{r}^{\prime})$ linear in $\Delta^{\alpha}(\bm{r})$. Iteration of Eq. (5) gives

[TABLE]

with the kernel

[TABLE]

where the normal-state Green’s function $G_{0,\omega_{n},\pm}^{\alpha}(\bm{r},\bm{r}^{\prime})$ satisfies $\left(\mathrm{i}\omega_{n}\!\mp\mu_{z}H+\hat{\xi}_{\bm{k}}^{\alpha}\right)G_{0,\omega_{n},\pm}^{\alpha}(\bm{r},\bm{r}^{\prime})\!=\!\delta\left(\bm{r}\!-\!\bm{r}^{\prime}\right)$.

In the following sections, we will utilize Eqs. (6), (7), and (8) to derive the conditions for superconducting instabilities in zero and finite magnetic fields.

IV Transition temperature at zero magnetic field

The influence of shallow bands on the transition temperature $T_{C}$ has been discussed in several recent papers Innocenti et al. (2010); Bang (2014); Chen et al. (2015); Chubukov et al. (2016); Valentinis et al. (2016). In this section, for completeness, we present the derivation of the transition temperature for our model. In the absence of the magnetic field, the band gap functions are homogeneous, $\Delta^{\alpha}(\bm{r})\!=\!\Delta^{\alpha}_{0}$. Substitution of these constant gaps into Eq. (7) gives equation which determines the transition temperature of the system with $\int d\bm{r}^{\prime}K_{\alpha}(\bm{r},\bm{r}^{\prime};\omega_{n})=\sum_{\mathbf{k}}\left[\omega_{n}^{2}+(\xi^{\alpha}_{\bm{k}})^{2}\right]^{-1}.$ The integration over the momentums in the $e$- ($h$-) band can be performed using the standard relations $\sum_{\mathbf{k}}\approx N_{e}\int^{\Omega}_{-\Omega}d\xi^{e}$ ($\approx N_{h}\int^{\mu_{h}}_{-\Omega}d\xi^{h}$) where $\Omega\gg T_{C}$ is the high-energy cutoff and $N_{\alpha}=m_{\alpha}/(2\pi)$ is the 2D density of states. The only difference from the standard BCS scheme is that for the shallow band the energy integration is limited by the band edge rather than by $\Omega$. The resulting gap equation can be presented as 444We replaced a smooth cutoff at $\omega_{n}\sim\Omega$ following from the $\xi^{\alpha}$-integration by the sharp cutoff, because the results weakly depend on the details of the high-energy behavior

[TABLE]

with the dimensionless coupling matrix $\hat{\Lambda}_{\alpha\beta}=U_{\alpha\beta}N_{\beta}$, and $\eta_{h}(\omega_{n})=\frac{1}{\pi}\tan^{-1}\frac{\mu_{h}}{\omega_{n}}$. The sign of the off-diagonal coupling constants $\Lambda_{eh}$ and $\Lambda_{he}$ determines relative sign of the order parameters in two bands. The case $\Lambda_{eh},\Lambda_{he}<0$ corresponds to $s_{\pm}$ superconducting state. In absence of interband scattering, this sign has no influence on the behavior of the upper critical field.

Introducing the following notations

[TABLE]

where

[TABLE]

and $\gamma_{E}\approx 0.5772$ is the Euler-Mascheroni constant, we can write the gap equation in a compact form as

[TABLE]

with $\hat{\Lambda}_{0}=\text{diag}[\Lambda_{0,h},\;\Lambda_{0,e}]$. In the limit $\mu_{h}\gg T_{C}$ the function $\Upsilon_{C}$ has the following asymptotics $\Upsilon_{C}\approx\tfrac{1}{2}\ln[2\mathrm{e}^{\gamma_{E}}\mu_{h}/(\pi T_{C})]$. In this limit both $\Lambda^{-1}_{0,h}$ and $\Lambda^{-1}_{0,e}$ have the same form $\ln(T_{0}/T_{C})$ but with different cutoff energies $T_{0}$. Defining the matrix $\hat{W}\equiv\hat{\Lambda}^{-1}-\hat{\Lambda}^{-1}_{0}$, we can present the equation for $T_{C}$ as the condition of degeneracy of this matrix,

[TABLE]

This is the instability condition for superconducting ground state. It leads to the explicit result for the effective coupling constant $\Lambda_{0,e}$, which directly determines $T_{C}$ by Eq. (10a), see Appendix A,

[TABLE]

with $\delta_{\Lambda}\!=\!-\mathrm{sign}\left[(1\!-\!\Upsilon_{C}\Lambda_{hh})\mathcal{D}_{\Lambda}\right]$ and $\mathcal{D}_{\Lambda}\!=\!\Lambda_{ee}\Lambda_{hh}\!-\!\Lambda_{eh}\Lambda_{he}$. The detailed investigation of the dependences of $T_{C}$ on $\mu_{h}$ for different pairing models has been performed in Ref. Chen et al. (2015). As our main goal is the investigation of the upper critical field, we only need Eqs. (11) and (13) as the zero-field references.

V Superconducting transitions in finite magnetic field

In the presence of the magnetic field the problem becomes nontrivial, since the superconducting states are not uniform. The upper critical field is mostly determined by the eigenvalues of the pairing kernels $K_{\alpha}(\bm{r},\bm{r}^{\prime};\omega_{n})$, (15). In this section, we describe evaluation of these kernels in magnetic field for the deep and shallow bands.

In the uniform magnetic field the Green’s function can be written as

[TABLE]

with $\phi_{\bm{A}}(\bm{r},\bm{r}^{\prime})=\frac{e}{c}\bm{A}\left(\frac{\bm{r}+\bm{r}^{\prime}}{2}\right)\cdot\left(\bm{r}-\bm{r}^{\prime}\right)$. In the symmetric gauge $\bm{A}(\bm{r})=\frac{1}{2}\bm{H}\times\bm{r}=\frac{H}{2}(-y,x,0)$, the phase factor in the exponent becomes $\phi_{\bm{A}}(\bm{r},\bm{r}^{\prime})=-\left[\bm{r\times r}^{\prime}\right]_{z}/(2l^{2})$, where $l=\sqrt{c/eH}$ is the magnetic length. This allows us to present the kernel in Eq. (8) as

[TABLE]

with $\bm{\rho}=\bm{r}^{\prime}-\bm{r}$.

In the isotropic case and for not too strong Zeeman spin splittingBuzdin and Brison (1996a); *Buzdin:PLettA218.1996, the shape of the gap function at the upper critical field is given by the ground-state eigenfunction of the particle with the charge $2e$ in the uniform magnetic field Helfand and Werthamer (1966); *Werthamer:PRev147.1966; Gruenberg and Gunther (1968); Rajagopal and Vasudevan (1966); Kogan and Prozorov (2012), i.e., $\Delta^{\alpha}(\mathbf{r})$ satisfies the equation $-\frac{1}{2}l^{2}\left(\nabla_{\mathbf{r}}-\mathrm{i}\frac{2e}{c}\bm{A}\right)^{2}\Delta^{\alpha}(\bm{r})=\Delta^{\alpha}(\bm{r})$. In the symmetric gauge, the shape of $\Delta^{\alpha}(\bm{r})$ is 555The complete ansatz gap function has to be multiplied by a arbitrary holomorphic function $f(x+iy)$. This function can only be fixed when the non-linear terms in $F^{\alpha}$ that account for the interactions between vortices are consideredTešanović et al. (1991); Maniv et al. (1992)

[TABLE]

For isotropic bands this $\Delta^{\alpha}(\bm{r})$ is an eigenfunction of the kernel $K_{\alpha}(\bm{r},\bm{r}^{\prime};\omega_{n})$, Eq. (15), for arbitrary function $g_{0}^{\alpha}(\rho,\omega_{n})$ Gruenberg and Gunther (1968). Indeed, substituting this ansatz into Eq. (7), we obtain

[TABLE]

with $\Phi_{\bm{r}}(\rho)=\int^{2\pi}_{0}\mathrm{d}\theta\exp\{-\frac{1}{l^{2}}(\mathrm{i}[\bm{r}\times\bm{\rho}]_{z}+\bm{r}\cdot\bm{\rho})\}$, where $\theta$ is the angle between the vectors $\bm{\rho}$ and $\bm{r}$. Noting that, $\mathrm{i}[\bm{r}\times\bm{\rho}]_{z}+\bm{r}\cdot\bm{\rho}=\mathrm{i}r\rho\sin\theta+r\rho\cos\theta$, we calculate the $\theta$-integral as

[TABLE]

This result can be easily obtained by expanding the exponential into the power series and noting that only the zeroth-order term survives after the $\theta$-integration. Therefore, presenting the anomalous Green’s function as

[TABLE]

where $N_{\alpha}=m_{\alpha}/2\pi$ is the density of states for the $\alpha$-band, we obtain the following general result for the dimensionless kernel eigenvalue

[TABLE]

Therefore, the problem of the upper critical field is mostly reduced to evaluation of these eigenvalues which in turn depend on approximations made for the normal-state Green’s functions $g_{0,\!\pm}^{\alpha}(\rho)$.

In general, depending on external conditions, the chemical potential may vary with the magnetic field, see, e.g., discussion in Refs. Maniv et al. (1992); Itskovsky et al. (2000). In the situation we consider here this variation can be neglected, because the deep band acts as a charge reservoir and fixes the chemical potential.

For the deep $e$-pocket the quasiclassical approximation can be employed on the Green’s function of the $e$-band. In contrast, for the small-size $h$-pocket the effects of Landau-level (LL) quantizations have to taken into account precisely. We will consider first the case of negligible spin splitting ($\mu_{z}\to 0$) when the quantization effects are most pronounced. After that, we will investigate in detail the role of spin-splitting effects.

V.1 Kernel eigenvalues without spin splitting

V.1.1 The deep e-band: quasiclassical approximation

The equation for the upper critical field in the quasiclassical approximation has been derived long time ago Helfand and Werthamer (1966), see also recent review Kogan and Prozorov (2012). Nevertheless, we include a minimum discussion of this well-known result in order to make a direct comparison with the latter calculation for the shallow $h$-band.

The essence of the quasiclassical approximation is to exploit the fact that the relevant length scale in the kernel $\rho$ is of the order of the coherence length $\xi$, which is much larger that the inverse Fermi wave vector $k_{F}^{-1}$. Also, typically the cyclotron frequency $\omega_{c}$ is much smaller than the Fermi energy. This allows us to neglect the Landau quantization and use the zero-field Green’s function $g_{0}^{e}(\rho,\omega_{n}),$

[TABLE]

As we neglect spin-splitting effects, we dropped the spin index in this function.

For the product of Green’s functions in Eq. (18), we obtain

[TABLE]

We introduce the variables $\bm{k}=\bm{\bar{k}}+\bm{q}/2,\bm{k}^{\prime}=\bm{\bar{k}}-\bm{q}/2$ and expand $\xi_{\bm{k}}^{e}\approx\xi_{\bm{\bar{k}}}^{e}+\bm{v}_{e}\bm{q}/2,\xi_{\bm{k}^{\prime}}^{e}\approx\xi_{\bm{\bar{k}}}^{e}-\bm{v}_{e}\bm{q}/2$. In quasiclassical regime $\bar{k}\!\sim\!k_{F}\gg q\!\sim\!1/\xi$. This allows us to approximately perform the integration over $\bm{\bar{k}}$ by using the standard transformation $\int\frac{d\bm{\bar{k}}}{(2\pi)^{2}}\approx\int\frac{\mathrm{d}k_{F}^{e}}{4\pi^{2}v_{e}}\int_{-\infty}^{\infty}d\xi^{e}$ and neglecting $\bm{\bar{k}}$ dependence of $\bm{v}_{e}$, which leads to the following result

[TABLE]

where $N_{e}=m_{e}/2\pi$ is the e-band density of states and $\langle\dots\rangle_{e}$ means averaging over the electron Fermi surface, $\langle\dots\rangle_{e}=\int\dots\frac{\mathrm{d}k_{F}^{e}}{4\pi^{2}v_{e}N_{e}}$.

Substituting this presentation into Eq. (18), we obtain

[TABLE]

We can further transform this results using the transformation $A^{-1}=\int_{0}^{\infty}\mathrm{d}s\mathrm{e}^{-sA}$, which leads to

[TABLE]

This gives the following well-known result

[TABLE]

Remark that, for the sake of simplicity, we consider here the case of an isotropic band meaning that the averaging $\left\langle\ldots\right\rangle_{e}$ can be omitted. In the case of a single band, generalization to elliptic anisotropy is straightforward. However, there is no accurate analytical description of multiple bands with different anisotropies. Without elaborated numerical calculations, this case can only be treated approximately Kogan and Prozorov (2012).

V.1.2 The shallow h-band: Landau-level quantization

In the shallow hole band the typical length scale of the kernel may be comparable with $k_{F}^{-1}$ and, in magnetic field, the cyclotron frequency may be comparable with the Fermi energy. This means that the quasiclassical approximation is not applicable and we have to use the exact normal-state Green’s functions $g_{0}^{h}(|\bm{r}-\bm{r}^{\prime}|,\omega_{n})$ in the kernel $K_{h}(\bm{r},\bm{r}^{\prime};\omega_{n})$, Eq. (15). In this case the shape of the kernel is influenced by the Landau quantization. For single-band materials, such exact presentation of the kernel was derived in several theoretical works Gruenberg and Gunther (1968); Vavilov and Mineev (1997); Mineev (2000); *Champel:PhilMagB81.2001; Zhuravlev and Maniv (2012); Ueta and Hioki (2013); Rajagopal and Ryan (1991). The normal-state Green’s function for the hole band is determined by the equation

[TABLE]

with $\bm{D}_{\mathbf{r}}=\nabla_{\mathbf{r}}-\mathrm{i}\frac{e}{c}\mathbf{A}(\bm{r})$. The solution is given by Eq. (14) with

[TABLE]

where $\rho=|\bm{r}-\bm{r}^{\prime}|$, $\omega_{c}=eH/(cm_{h})$, and $L_{\ell}(x)$ are the Laguerre polynomials.

There are several routes to transform and simplify the kernel eigenvalue $\lambda^{h}_{\omega_{n}}$ in Eq. (18). Using the integral representation $\{\omega_{n}\!\pm\mathrm{i}[\omega_{c}(\ell+\tfrac{1}{2})-\mu_{h}]\}^{-1}\!=\!\zeta_{\omega}\int^{\infty}_{0}\mathrm{d}s\exp\left(-\zeta_{\omega}s\left\{\omega_{n}\!\pm\!\mathrm{i}\left[\omega_{c}(\ell\!+\tfrac{1}{2})\!-\mu_{h}\right]\right\}\right)$ with $\zeta_{\omega}\equiv\mathrm{sign}(\omega_{n})$ and the generating function of Laguerre polynomials

[TABLE]

we can carry out the summation over the Landau levels Zhuravlev and Maniv (2012). After that, the integration over $\rho$ can be done exactly (see Appendix C) leading to

[TABLE]

where we introduced the dimensionless variables $\bar{\omega}_{n}\!=\!\omega_{n}/\omega_{c}$, and $\bar{\mu}_{h}\!=\!\mu_{h}/\omega_{c}$. We can see that the replacement $\zeta_{\omega}\to-\zeta_{\omega}$ is equivalent to the interchange $s_{1}\leftrightarrow s_{2}$ and, therefore, the factor $\zeta_{\omega}$ can be dropped meaning that $\lambda^{h}_{\omega_{n}}$ is even function of $\omega_{n}$.

The above presentation of $\lambda^{h}_{\omega_{n}}$ can be further transformed by breaking the $\bar{s}_{1}$- and $\bar{s}_{2}$-integrations into infinite sums, $\int^{\infty}_{0}\mathrm{d}\bar{s}_{i}=\sum^{\infty}_{n=0}\int^{2(n+1)\pi}_{2n\pi}\mathrm{d}\bar{s}_{i}$, and changing of variables $\bar{s}=\frac{1}{2}(\bar{s}_{1}+\bar{s}_{2})$ and $\bar{u}=\frac{1}{2}(\bar{s}_{1}-\bar{s}_{2})$. This gives us the following result (see Appendix C for details)

[TABLE]

where

[TABLE]

We can observe that the denominator in Eq. (24) oscillates with $\omega_{c}$ and has minimums at $\bar{\mu}_{h}=\text{integer}+\frac{1}{2}$ corresponding to matching of the chemical potential with the Landau levels. In high-field and low-temperature regime, since $\cosh(2\pi\bar{\omega}_{n})\sim 1$, the denominator produces strong peaks in $\lambda^{h}_{\omega_{n}}$ at $\bar{\mu}_{h}=\ell+\frac{1}{2}$, which diverge at zero temperature. The identical oscillating factor also appears in the quasiclassical result for the kernel eigenvalue Maniv et al. (1992).

V.2 Equation for the upper critical field

To study the superconducting state near $H_{C2}$, we can just substitute the result for $F_{\omega_{n}}^{\alpha}(\bm{r})$ in Eq. (17) into the gap equation (6) which leads to

[TABLE]

However, similar to the zero-field case, this gap equation contains logarithmic divergences as $\Omega\to\infty$ (UV divergences) which has to be cut at $\omega_{n}\sim\Omega$. These logarithmic UV divergences in $\sum_{\omega_{n}}\lambda^{\alpha}_{\omega_{n}}$ can be compensated by explicitly subtracting

[TABLE]

from both side of the gap equation. Using definitions in Eqs. (10a) and (10b), this leads to the following regularized gap equation (see Appendix B)

[TABLE]

where $\mathcal{A}_{1}=\frac{1}{2}\ln t-\Upsilon_{T}+\Upsilon_{C}$, $\mathcal{A}_{2}=\ln t$, $t=T/T_{C}$,

[TABLE]

$\Upsilon_{C}\equiv\Upsilon_{T_{C}}$, and

[TABLE]

Now the right hand sides of the above equations remain finite with $\Omega\to\infty$, since the logarithmic divergences are canceled by the $1/\omega_{n}$ terms. Assuming $T_{C},\omega_{c}\ll\Omega$, we took the limit $\Omega\to\infty$ in the frequency sums. All the information about UV cutoff is absorbed by the parameter $T_{C}$. The functions $\mathcal{A}_{\alpha}(T)$ and $\mathcal{J}_{\alpha}(H,T)$ are defined in such a way that $\mathcal{A}_{\alpha}(T)\to 0$ for $T\to T_{C}$ and $\mathcal{J}_{\alpha}(H,T)\to 0$ for $H\to 0$. In Eq. (28b) the summation over the Matsubara frequencies can be carried out leading to the following well-known presentation Kogan and Prozorov (2012)

[TABLE]

The upper critical field is the magnetic field at which a nontrivial solution of the linear gap equation, Eq. (27), appears. This corresponds to the condition

[TABLE]

As the matrix $\hat{W}$ is degenerate, this leads to the concise equation

[TABLE]

which determines superconducting instability in the magnetic field. The constants $W_{\alpha\alpha}$ can be directly connected with the coupling constants as

[TABLE]

with $\mathcal{D}_{\Lambda}=\Lambda_{ee}\Lambda_{hh}-\Lambda_{eh}\Lambda_{he}$, $\delta_{W}=\mathrm{sign}[\mathcal{D}_{\Lambda}(1-\Upsilon_{C}\Lambda_{hh})]$, and

[TABLE]

All information about the coupling matrix is contained in these two constants, $W_{11}$ and $W_{22}$, which also weakly depend on the ratio $\mu_{h}/T_{C}$. These parameters are typically large in absolute values because they scale as $\Lambda_{\alpha\beta}^{-1}$, but they can be either positive or negative depending on the sign of the determinant $\mathcal{D}_{\Lambda}$. The relative contribution of the band $\alpha$ to the superconducting instability is inversely proportional to $|W_{\alpha\alpha}|$.

The behavior of the upper critical field is mostly depends on the field and temperature dependences of the functions $\mathcal{J}_{\alpha}$ which determine the field-induced contributions to the pairing kernels. The quasiclassical kernel $\mathcal{J}_{2}$ has monotonic field and temperature dependences. If only the deep band is present, the conventional monotonic upper critical field is determined by the equation $\mathcal{J}_{2}(H,T)=\ln t$. In contrast, due to the Landau-level quantization, $\mathcal{J}_{1}(H,T)$ is an oscillating function of the magnetic field at low temperatures and this leads to the anomalous behavior of the upper critical field. In the next section we discuss in details the behavior of the kernel $\mathcal{J}_{1}$.

V.3 Shape of the quantum field-dependent pairing kernel $\mathcal{J}_{1}(H,T,\mu_{h})$ without spin splitting

In this section, we examine in detail behavior of the function $\mathcal{J}_{1}$ in Eq. (28a). This function can be evaluated numerically for any temperature except $T\!=\!0$. We only present the results that are relevant to the discussions, and leave the mathematical details in Appendix D.

The direct numerical evaluation of $\mathcal{J}_{1}$ from $\lambda^{h}_{\omega_{n}}$ in Eq. (24) is doable but not very efficient. To derive presentation better suitable for numerical evaluation, one can trade the slowly-converging frequency sum to another rapidly-convergent series sum, see Appendix D,

[TABLE]

where $\bar{\tau}=\pi^{2}T/\omega_{c}$, $z_{j}^{\pm}=j\!\pm\!\left(1\!-\!\frac{\bar{s}}{\pi}\right)$, and derivative of the function $\mathcal{I}(\bar{s})$, Eq. (25), can be transformed to the following form

[TABLE]

In Fig. 2(a), we present the numerically calculated $\mathcal{J}_{1}$ within the range $\omega_{c}/\mu_{h}\in[0.05,3.5]$ for three temperatures, $T/\mu_{h}\!=\!0.02$, $0.1$, and $0.5$. We can see that at low temperatures $\mathcal{J}_{1}$ is strongly oscillating function with the peaks at $\omega_{c}/\mu_{h}\!=\!1/(\ell\!+\!1/2)$. The strongest peak is realized at the lowest Landau level, $\ell\!=\!0$. The peak amplitudes rapidly decrease with increasing temperature so that at $T=0.5$ the function $\mathcal{J}_{1}$ is already monotonic.

At low temperatures and $\bar{\mu}_{h}$ not close to half-integers, the $\mathcal{J}_{1}$ can be approximated as (see Appendix D.1)

[TABLE]

with $\Upsilon_{T}\approx\tfrac{1}{2}\ln[2\mathrm{e}^{\gamma_{E}}\mu_{h}/(\pi T)]$. Note that the first two terms are logarithmically divergent as $T\to 0$. They exactly cancel with corresponding divergent terms in left-hand-side of the gap equation (27) so that $\mathcal{A}_{1}-\mathcal{J}_{1}$ approaches a finite value at $T\to 0$. This expansion breaks down for the values of $\mu_{h}$ close to the Landau-levels, $\mu_{h}=\omega_{c}(\ell+1/2)$. In the vicinity of the Landau levels the applicability condition of this asymptotics becomes $T\ll|\mu_{h}-\omega_{c}(\ell+1/2)|$. In particular, near the lowest Landau level $\mu_{h}\sim\omega_{c}/2$, we derive in Appendix D.1 the following presentation

[TABLE]

We can see indeed that the low-temperature asymptotic is realized for $T\ll|2\mu_{h}\!-\!\omega_{c}|$. At $\mu_{h}\!=\!\omega_{c}/2$ the function $\mathcal{J}_{1}$ diverges as $\omega_{c}/(8T)$ for $T\to 0$. Similar behavior is realized at higher Landau levels, for $\mu_{h}\!=\!\omega_{c}(\ell+1/2)$ the function $\mathcal{J}_{1}$ diverges as $[(2\ell)!/(2^{\ell}\ell!)^{2}]\omega_{c}/(8T)$, see Appendix D.4. These divergencies were pointed out, e.g., in Ref. Gruenberg and Gunther (1968). They reflect enhanced Cooper pairing due to $\delta$-function singularities of the density of states at the Landau levels.

For better exposition of the LL-quantization effects, we derive in Appendix D.2 an approximate result for $\mathcal{J}_{1}$ in which these effects are completely neglected,

[TABLE]

This result describes behavior of $\mathcal{J}_{1}$ in the limit $\omega_{c}\ll\pi T,\mu_{h}$ and is similar to quasiclassical approximation except that it is valid for arbitrary relation between $\mu_{h}$ and $T$. It corresponds to the form of the Green’s function given by Eq. (14) in which the magnetic field is only taken into account in the phase factor $\phi_{\bm{A}}(\bm{r},\bm{r}^{\prime})$ and for $g_{0}^{h}(|\bm{r}-\bm{r}^{\prime}|,\omega_{n})$ the exact zero-field Green’s function is substituted. The function $\mathcal{J}_{1}^{\mathrm{qc}}$ reduces to the standard quasiclassical result similar to Eq. (29) in the limit $\mu_{h}\gg T$. In the Fig. 2(a) we plot the function $\mathcal{J}_{1}^{\mathrm{qc}}$ together with exact results and we see that this monotonic function well reproduces the exact shape of $\mathcal{J}_{1}$ whenever the quantum oscillations become small due to temperature smearing. It is important to note that the condition $\omega_{c}\ll\pi T,\mu_{h}$ does not yet imply that $\mathcal{J}_{1}\propto\omega_{c}$. This linear low-field asymptotics formally requires the condition $\omega_{c}\ll T^{2}/\mu_{h}$ and in the case $T\ll\mu_{h}$ the parameter $T^{2}/\mu_{h}$ is much smaller than both $T$ and $\mu_{h}$. We can indeed see in Fig. 2(a), that the approximation in Eq. (36) well reproduces the exact result at small $\omega_{c}$ even in the region where $\mathcal{J}_{1}(\omega_{c})$ is strongly nonlinear.

V.4 The kernel eigenvalue and the functions $\mathcal{J}_{\alpha}$ with finite spin splitting

As discussed in the introduction, the relative role of spin-splitting effects on suppression of superconductivity is characterized by the Maki parameter $\alpha_{M}$ which in clean case scales inversely proportional to the Fermi energy. Therefore, one can expect that these effects may be essential for shallow bands. The kernel eigenvalues $\lambda^{\alpha}_{\omega_{n}}$ can be straightforwardly generalized to the case with finite Zeeman effects by replacing $|\omega_{n}|\to\zeta_{\omega}(\omega_{n}+\mathrm{i}\mu_{z}H)$ in Eqs. (20) and (24). Therefore, we have

[TABLE]

Here we introduced the parameter $\gamma_{z}\!=\!\mu_{z}H/\omega_{c}\!=\!\mu_{z}m_{h}c/e\!=\!gm_{h}/4m_{0}$ characterizing the relation between the spin-splitting energy and Landau-level separation. Here $m_{0}$ is the free-electron mass and $g$ is the spin g-factor. For free electrons $\gamma_{z}\approx 0.5$. As the cyclotron frequency is determined by the z-axis component of the magnetic field and the spin-splitting energy is determined by the total field, the effective spin-splitting factor can be enlarged by tilting the magnetic field away from the z axisWosnitza (1996), for field tilted at the angle $\theta$ with respect to the z axis, $\gamma_{z}(\theta)=\gamma_{z}(0)/\cos\theta$.

With finite Zeeman splitting, the eigenvalues $\lambda^{\alpha}_{\omega_{n}}$ become complex and one has to take the real part of the right-hand sides in the definitions of the functions $\mathcal{J}_{i}$, Eqs. (28). Similarly to zero spin-splitting case, we can trade the Matsubara-frequency sum to the fast convergent series. Derivations presented in Appendix D.3 give the following presentations

[TABLE]

where the functions $\mathcal{I}(\bar{s})$, $\mathcal{I}^{\prime}(\bar{s})$ are given by Eqs. (25), (33) and $z_{j}^{\varsigma}=j+\varsigma(1-\bar{s}/\pi)$. Alternatively, one can derive a presentation for $\mathcal{J}_{1}$, in which the summations over the Landau levels are preserved, see Appendix D.4,

[TABLE]

This presentation is similar to one derived and used in Refs. Rajagopal and Vasudevan (1966); Gruenberg and Gunther (1968); Maniv et al. (1992). Even though the presentations in Eqs. (38a) and (39) look very different, they do describe the same function and can be used for studying different properties of this function.

The derivations of the low-temperature asymptotics of $\mathcal{J}_{1}$ for different cases are presented in Appendix D.4. For noninteger $2\gamma_{z}$, the function $\mathcal{J}_{1}\!-\!\mathcal{A}_{1}$ approaches finite limits at $T\!\to\!0$ for any value of $\mu_{h}$. For small $\gamma_{z}$ these limiting values are large at the shifted LLs, $\mathcal{J}_{1}\!-\!\mathcal{A}_{1}\!\approx\!\!\frac{(2\ell_{0})!}{2^{2\ell_{0}}(\ell_{0}!)^{2}}\frac{1}{4\gamma_{z}}$ for $\mu_{h}\!=\omega_{c}(\ell_{0}+\tfrac{1}{2}\pm\gamma_{z}\pm 0)$. When $2\gamma_{z}$ equals integer $j_{z}$ the spin-splitting energy $2\gamma_{z}\omega_{c}$ exactly matches the LL spacing, see example in the inset of Fig. 2(c) for the free-electron spin-splitting, $j_{z}\!=\!1$. In these resonance cases, the function $\mathcal{J}_{1}-\mathcal{A}_{1}$ again becomes divergent at low temperatures for $\mu_{h}=\omega_{c}\left(\ell_{0}+j_{z}/2+1/2\right)$,

[TABLE]

However, the numerical coefficient in this asymptotics rapidly decreases with $j_{z}$. In particular, for the lowest Landau level, $\ell_{0}=0$, $\mathcal{J}_{1}-\mathcal{A}_{1}\simeq 2^{-j_{z}}\omega_{c}/(8T)$.

We plot the functions $\mathcal{J}_{1}(\omega_{c})$ at $T=0.02\mu_{h}$ for different $\gamma_{z}$ in Fig. 2(b). The spin-splitting effects effectively suppress the spin-singlet pairing in each LLs already at small values of $\gamma_{z}$ leading to rapid suppression of the $\mathcal{J}_{1}$ peaks. The peaks are replaced by the downward and upward steps at $\omega_{c}/\mu_{h}\!=\!(\ell\!+\!\gamma_{z}\!+\!\tfrac{1}{2})^{-1}$ and $\omega_{c}/\mu_{h}\!=\!(\ell\!-\!\gamma_{z}\!+\!\tfrac{1}{2})^{-1}$ respectively. However, the pair-breaking effect of spin-splitting is somewhat reduced for integer values of $2\gamma_{z}=j_{z}$, see inset in Fig. 2(c). For these special values, the peaks in $\mathcal{J}_{1}(\omega_{c})$ reappear at $\mu_{h}=\omega_{c}\left(\ell_{0}+j_{z}/2+1/2\right)$. The upper critical field is directly determined by the difference $\mathcal{J}_{1}\!-\!\mathcal{A}_{1}$, which has finite limit at $T\to 0$ for all parameters except the resonance values of $\gamma_{z}$ and $\omega_{c}/\mu_{h}$. The contour plot of this function in the plane $\omega_{c}/\mu_{h}$–$\gamma_{z}$ is shown in Fig. 2(c) and provides somewhat clearer illustration of the general behavior with increasing $\gamma_{z}$. One can again see that this function has steps when the chemical potential crosses the Zeeman-shifted Landau levels and very sharp peaks at the resonance parameters $\gamma_{z}=2j_{z}$ and $\mu_{h}/\omega_{c}\!=\!\ell_{0}\!+\!j_{z}/2\!+\!1/2$ reflecting enhancement of the pairing strength. We can also see that away from these steps and peaks the dependence on $\gamma_{z}$ is very weak. In the next section we will use the derived formulas for the kernels $\mathcal{J}_{\alpha}$ to compute the upper critical fields for different coupling matrices and spin-splitting parameters.

VI Temperature-field phase diagrams: reentrant Landau-level regions

In previous sections we derived general relations which determine the superconducting instabilities in the magnetic field for clean two-dimensional superconductors with two bands, deep and shallow. At this stage we have all the ingredients to determine the upper critical field in such a system. In this section, we discuss the shapes of the magnetic field-temperature phase diagrams for several representative cases. First, we present simple analytical results for different limits.

For $T\to T_{C}$ and $\omega_{c}\to 0$, we can keep only linear terms in $\mathcal{J}_{\alpha}$ with respect to $H$ (see Appendix D.2), $\mathcal{J}_{\alpha}\approx-\,H\mathcal{Y}_{\alpha}$ with

[TABLE]

where $\zeta(x)$ is the Riemann zeta-function, $\zeta(3)\approx 1.202$, and expand $\mathcal{A}_{\alpha}$ with respect to $1-t$,

[TABLE]

We remark that the shallow-band results in Eqs. (40a) and (41a) are somewhat different from the WH approach. As near the transition temperature $|(\mathcal{A}_{\alpha}-\mathcal{J}_{\alpha})/W_{\alpha\alpha}|\ll 1$, the equation for $H_{C2}$, Eq. (30), becomes

[TABLE]

Substituting the linear expansions for $\mathcal{A}_{\alpha}$ and $\mathcal{J}_{\alpha}$, we obtain

[TABLE]

Furthermore, for the shallow band, $\mu_{h}\ll\mu$, and this implies $\mathcal{Y}_{1}\ll\mathcal{Y}_{2}$. This allows us to simplify $H_{C2}$ near $T_{C}$ as

[TABLE]

Near $T_{C}$, the quantum and spin-splitting effects are negligible and the shallow band gives a relatively small correction to $H_{C2}$.

As demonstrated in Sec. V.3, without the spin-splitting effects the function $\mathcal{J}_{1}(H,T)$ diverges for $T\!\to 0$ as $1/T$ at $\omega_{c}\!=\!\mu_{h}/(\ell+1/2)$. As a consequence, the transition temperature is usually finite at these field values. For the lowest LL ($2\mu_{h}\!=\!\omega_{c}$), this transition temperature $T^{(0)}_{C2}$ can be calculated from Eq. (30) in the case $T^{(0)}_{C2}\ll T_{C},\mu_{h}$ by using the low-temperature asymptotics of $\mathcal{J}_{1}$, Eq. (35), and $\mathcal{J}_{2}$, Eq. (29), see Appendix E.1, which yields

[TABLE]

with

[TABLE]

where $H_{c2}^{e}=(2\pi^{2}/\mathrm{e}^{\gamma_{E}})cT_{C}^{2}/ev_{e}^{2}$ is the orbital upper critical field of the deep band. We focus on the regime $H>H_{c2}^{e}$ meaning that $r_{C}>1$.

For special values of spin-splitting parameters $2\gamma_{z}=j_{z}$, the transition temperature may be also finite at $\omega_{c}\!=\!\mu_{h}/(\ell\!+\!j_{z}/2\!+\!1/2)$. In particular, we derive in appendix E.2 the transition temperature $T_{C2}^{(1)}$ for the important particular case of free-electron spin splitting, $j_{z}\!=\!1$, and the lowest resonance field, $\omega_{c}\!=\!\mu_{h}$,

[TABLE]

where $\tilde{r}_{C}=r_{C}[1\!+\!2m_{e}\omega_{c}\gamma_{z}^{2}/(m_{h}\mu)]$ accounts for weak Zeeman correction in the deep band. The result for $T_{C2}^{(1)}$ is similar to $T_{C2}^{(0)}$ but contains a smaller numerical factor.

The transition temperatures $T_{C2}^{(i)}$ emerge as a result of the interplay between the pairing strengths in two bands which is accounted for by the first term within the square brackets in Eqs. (45) and (47). These temperatures are finite if the expressions inside the square brackets are positive which is true for most parameter sets 666Strictly speaking, the transition temperatures $T_{C2}^{(i)}$ are not always finite, in contrast to the single-band case. As follows from Eqs. (45) and (47), $T_{C2}^{(i)}$ vanish at $r_{C}=\exp(-2W_{22})$ which corresponds to the case $\Lambda_{eh}\Lambda_{he}>\Lambda_{ee}\Lambda_{hh}$ and to either very deep electron band or very strong interband coupling. The values of $T_{C2}^{(i)}$ are determined not only by overall strength of the Cooper pairing but also by the relative weights with which two bands contribute to superconducting instability. Therefore, they are very sensitive to the coupling-matrix structure.

In particular, for the dominating deep-band coupling, $\Lambda_{ee}\!>\!\Lambda_{hh},|\Lambda_{eh}|,|\Lambda_{he}|$, a noticeable reentrant $T_{C2}$ only appears for sufficiently strong interband couplings. Indeed, in this scenario, the constants $W_{\alpha\alpha}$ can be estimated as $W_{11}\approx\Lambda_{ee}/\mathcal{D}_{\Lambda}$ and $W_{22}\approx\Lambda_{eh}\Lambda_{he}/(\Lambda_{ee}\mathcal{D}_{\Lambda})$ with $|W_{22}|\ll|W_{11}|$. In the case $\ln r_{C}\ll|W_{22}|$ and $W_{11}/W_{22}\!\approx\!\Lambda_{ee}^{2}/\Lambda_{eh}\Lambda_{he}\!\gg\!\Upsilon_{C},\ln(\mu_{h}/T_{c})$, we obtain a simple estimate, $T_{C2}^{(0)}\approx\mu_{h}\Lambda_{eh}\Lambda_{he}/(2\ln r_{C}\Lambda_{ee}^{2})$ showing that $T_{C2}^{(i)}$ indeed vanish for $\Lambda_{eh},\Lambda_{he}\!\to\!0$. We also see that in this case $T_{C2}^{(i)}$ decrease with the increasing deep-band coupling constant $\Lambda_{ee}$. Such counterintuitive behavior is caused by the reduction of the shallow-band weight at the superconducting instability.

In the opposite limit of the dominating interband coupling $|\Lambda_{he}|,|\Lambda_{eh}|\!\gg\!\Lambda_{ee},\Lambda_{hh}$, assuming that $\Upsilon_{C}\!\ll\!1/\sqrt{\Lambda_{eh}\Lambda_{he}}$, we obtain $W_{11}\!\approx\!W_{22}/2\!\approx\!-1/\sqrt{2\Lambda_{eh}\Lambda_{he}}$. In the limit $\ln r_{C}\!\ll\!1/\sqrt{\Lambda_{eh}\Lambda_{he}}$ we obtain a simple estimate for the transition temperature

[TABLE]

which does not depend on coupling constants at all.

For further understanding the relative role of the deep-band and intraband coupling strengths, we analyze in more detail the case of vanishing pairing in the shallow band $\Lambda_{hh}\!=\!0$. We consider the evolution of $T^{(0)}_{C2}$ with the increasing interband coupling, assuming that $T_{C}$ is fixed, meaning that the effective coupling $\Lambda_{0,e}$ in Eq. (13) remains unchanged. In this case $\Lambda_{ee}=\Lambda_{0,e}$ for $\Lambda_{eh}\Lambda_{he}=0$ and $\Lambda_{eh}\Lambda_{he}=\Lambda_{0,e}^{2}/(\frac{1}{2}+\Lambda_{0,e}\Upsilon_{C})$ for $\Lambda_{ee}=0$. In the case $\Lambda_{hh}\!=\!0$, we can strongly simplify presentations for the parameters $W_{\alpha\alpha}$ in Eqs. (31) by relating them with $\Lambda_{0,e}$,

[TABLE]

This allows us to rewrite the result for $T^{(0)}_{C2}$, Eq. (45), more transparently as

[TABLE]

Similar presentation can be obtained for $T^{(1)}_{C2}$. We can see that at fixed $T_{C}$ the temperatures $T^{(i)}_{C2}$ monotonically increase with the interband couplings and have some tendency to saturation when these couplings become large.

The dependence of $T^{(i)}_{C2}$ on $\mu_{h}$ and, correspondingly, on $\omega_{c}\!\propto\!H$ is characterized by the three typical scales: the transition temperature $T_{C}$, the value at which the Landau-level magnetic field matches the deep-band upper critical field $\mu_{C2}^{(i)}=c_{i}m_{e}T_{C}^{2}/m_{h}\mu$ with $c_{0}\approx c_{1}/2\approx 2.77$, and the large scale $\mu_{\Lambda}^{(i)}=\mu_{C2}^{(i)}\exp(2/\Lambda_{0,e})$. In the range $T_{C},\mu_{C2}^{(i)}\ll\mu_{h}\ll\mu_{\Lambda}^{(i)}$ the temperatures $T^{(i)}_{C2}$ increase with $\mu_{h}$, similar to the prediction for the single-band caseGruenberg and Gunther (1968); Tešanović et al. (1991); Maniv et al. (2001). For $\Lambda_{eh}\Lambda_{he}\ll\Lambda_{0,e}^{2}$, the function $T_{C2}^{(i)}(\mu_{h})$ has minimum at $\mu_{h}\!=\!\mathrm{e}\mu_{C2}^{(i)}$. At larger interband couplings the minimum is displaced to larger values which are determined by interplay between $T_{C}$ and $\mu_{C2}^{(i)}$. In the case $\mu_{\Lambda}^{(i)}\ll\epsilon_{0}$, which may realize only for very deep band, $T_{C2}^{(i)}(\mu_{h})$ reaches maximum for $\mu_{h}=\mu_{\Lambda}^{(i)}/\mathrm{e}$ and vanishes at $\mu_{h}=\mu_{\Lambda}^{(i)}$. The latter behavior, however, corresponds to very large magnetic fields and, probably, it is of only academic interest.

Figure 3 shows the dependences of $T^{(i)}_{C2}$ on the interband coupling $\Lambda_{he}$ for three values of $\mu_{h}/T_{C}$ and the representative parameters $\Lambda_{0,e}=0.2$, $m_{h}=m_{e}$, and $\epsilon_{0}=10T_{C}$. For this choice of parameters $\mathrm{e}\mu_{C2}^{(0)}\approx 0.75T_{C}$ and $\mathrm{e}\mu_{C2}^{(1)}\approx 1.5T_{C}$. Consistent with above estimates, $T^{(0)}_{C2}$ increases with $\mu_{h}$ for all $\Lambda_{he}$, while $T^{(1)}_{C2}$ weakly depends on $\mu_{h}$ at small $\Lambda_{he}$ and decreases with $\mu_{h}$ at large $\Lambda_{he}$. The maximum $T^{(1)}_{C2}\approx 0.17T_{C}$ realized for purely interband coupling case is roughly two times smaller than the maximum $T^{(0)}_{C2}\approx 0.34T_{C}$. The values of $T^{(0)}_{C2}$ and $T^{(1)}_{C2}$ provide natural measures for the strength of the high-field reentrant behavior which we discuss below.

In the whole temperature-magnetic field region, we computed the superconducting instability boundaries from Eq. (30) for several parameter sets. We consider first the case of zero spin splitting. Figure 4 shows the typical phase diagram in this case for representative parameters. The field scale in this and other plots $cm_{h}T_{C}/e$ is around 37 T for $T_{C}\!=50$K and $m_{h}\!=$ free-electron mass. The most remarkable feature is the existence of the reentrant superconducting regions at high magnetic fields whenever the highest occupied LL crosses $\mu_{h}$. These regions appear due to sharp enhancement of the density of states at these magnetic fields. At higher temperatures the thermal smearing of the Landau levels erases the quantization effects. As the result, the reentrant states disappear and the $H_{C2}$ curve approaches the quasi-classical result. The reentrance effect is most pronounced for the lowest Landau level and in the following consideration we mostly concentrate on this case.

The specific behavior is sensitive to the structure of the coupling matrix. In particular, it is quantitatively different for two cases discussed in the introduction, interband-coupling scenario and induced superconductivity in the shallow band. Figure 5 shows evolution of the temperature-field diagram with $\mu_{h}/T_{C}$ without spin splitting for these two cases. The qualitative behavior is similar in both cases, with increasing the chemical potential the strong bump appears at low temperatures and then it separates from the main domain and becomes a separate high-field superconducting region. The size of this reentrance region is much larger for the interband-coupling case, in which maximum $T_{C2}$ almost reaches $T_{C}/3$. These numerically computed $T_{C2}$ are in perfect agreement with the analytical result, Eq. (45).

The spin-splitting effects rapidly suppress the high-field reentrant regions. This can be seen in Fig. 6, in which we plot the dependence of superconducting boundaries at low temperatures on the spin-splitting factor $\gamma_{z}$ for the case of induced superconductivity in the shallow band (the same parameters as in the right plot of Fig. 5). A very small value $\gamma_{z}\approx 0.05$ is already sufficient to eliminate the separated region. In the interband-coupling scenario this value is somewhat larger, $\gamma_{z}\approx 0.1$. Another noticeable feature in Fig. 6 is a significant suppression of $T_{C2}$ at the magnetic fields for which the chemical potential falls in between the spin-up and spin-down Landau levels. This leads to stepwise behavior of the main boundary with steps corresponding to crossing the spin-down Landau levels and may cause the appearance of normal regions inside the superconducting domain. We also see in Fig. 6 that the reentrant regions reappear near the integer values of $2\gamma_{z}$, 1 and 2, corresponding to crossing of Landau levels with different spin orientations shown by the dashed lines. The reentrance is well developed for $\gamma_{z}\!\approx\!0.5$ when the chemical potential matches coinciding the spin-up 0th and spin-down 1st Landau levels. Figure 7 (left) illustrates the low-temperature part of phase diagrams for resonant values of spin splitting, $\gamma_{z}=0$, $0.5$, and $1$, for the same parameters as in Fig. 6. We can see that the size of reentrant region decreases with increasing $\gamma_{z}$. When $\gamma_{z}$ deviates from the resonance values, the reentrance rapidly disappears. Before disappearance, two small reentrant domains typically exist at fields corresponding to matching of the chemical potential with LL for two spin orientations, as illustrated in Fig. 7 (right).

Figure 8 illustrates evolution of the temperature-magnetic field phase diagrams with decreasing chemical potential for $\gamma_{z}=0.5$ in the interband-coupling case. We used the same parameters as in the left plot of Fig. 5. We can see that the behavior is similar to the case $\gamma_{z}=0$ except that the reentrant regions are smaller and their location is shifted to different values of $\mu_{h}/T_{C}$. In the case of induced superconductivity into the shallow band corresponding to the right plot of Fig. 5 the behavior is similar but the reentrant regions are even smaller.

VII Summary and discussion

In multiple-band superconductors, the shallow band can play an important role in spite of its low carriers concentration. In the presence of high magnetic field, the highest occupied LL has a very low quantum number. As a consequence, the Landau quantization causes the reentrant high-field superconducting regions at the low temperatures. The quantitative behavior depends on the relative strength of intra and interpocket pairing interactions. The reentrance is most pronounced when the interpocket coupling dominates.

The Zeeman spin splitting rapidly suppresses the high-field reentrant regions. However, such regions reappear in the special cases when spin-splitting energy exactly matches the LL spacing. The magnitude of the Zeeman term is determined by two factors, the $g$ factor and the band effective mass. In real materials both these factors may significantly differ from the free-electron values. In particular, the relative role of spin-splitting is reduced for light quasiparticles due to higher Landau-level energies, see, e. g., Ref. Zhuravlev et al. (2017).

In this paper we limited ourselves to the case when at the superconducting instability the lowest-Landau-level gap solution, Eq. (16), realizes. It was demonstrated in Refs. Buzdin and Brison (1996a); *Buzdin:PLettA218.1996 that in the case of small-size single 2D band and large spin splitting, the gap shape at the superconducting instability may be given by the higher-Landau-level wave functions. We verified that this does not happen for the parameter range we considered. Such scenario typically requires smaller Fermi energy for the deep band, $\epsilon_{0}\lesssim 5T_{C}$.

We only considered clean superconductor with very small scattering rate. In general, impurities are expected to have the same effect as the temperature which broadens the LLs leading to diminishing of the reentrant behavior. Similar to other quantum oscillations, we expect that this behavior persists until $\tau\omega_{c}\gg 1$, where $\tau$ is the scattering time.

We mention that a different orbital mechanism for the reappearance of superconductivity at high magnetic fields was predicted for quasi-one- and quasi-two-dimensional metals when the field is applied along the high-conductivity directions Lebed (1986); *Lebed:PhysRevLett1998. In this case the restoration of superconductivity is caused by the interplay between the quantum orbital motion of quasiparticles in the direction perpendicular to the conducting chains or planes and interchain/interplane periodicity.

The reentrant superconductivity in high magnetic field similar to one predicted here has been observed in Eu-doped Chevrel phases, EuxSn1-xMo6S8, with $T_{C}\approx 4$K Meul et al. (1984); *Rossel:JAP1985. This material has a wide isolated semi-elliptical superconducting region for $T<1$K and 4T$<H<$22.5T. The quasiclassical model used for the interpretation of this behavior assumed very weak orbital effect of magnetic field with Maki parameter $\alpha_{M}\approx 4.8$ and the compensation of the Zeeman-splitting effects due to interaction of quasiparticles with local magnetic moments, Jaccarino-Peter effectJaccarino and Peter (1962). While the second assumption looks very reasonable due to the presence of the magnetic Eu ions, the reason for extreme weakness of the orbital effects in this material is not very clear. We can not exclude that quantization effects play a role in the formation of the reentrant region in this material.

The presence of tunable shallow bands, as well as high values of the transitions temperatures and upper critical fields make FeSCs natural candidates for the reentrant behavior. An essential requirement is a sufficiently strong pairing interaction in between deep and shallow bands. Observation of the very large superconducting gap in the shallow hole band of LiAsFe Miao et al. (2015) suggests that such strong interband coupling is indeed present at least in some FeSC compounds. In this paper we limited ourselves to the two-dimensional case for which the quantization effects are the strongest. At the qualitative level, we expect that our results are applicable to the FeSe monolayer on SrTiO3 for which the Lifshitz transition has been reported recently Shi et al. (2016) or for intercalated FeSe compounds. The reentrant behavior is expected when the chemical potential of the shallow band is tuned to transition temperature. For $T_{C}\sim 50$K this corresponds to $\mu_{h}\sim 4.3$ meV, which is about 10 times smaller than the Fermi energy of the deep electron band at M point. Experimental probes of the electronic spectrum in the bulk FeSe by quantum oscillations Terashima et al. (2014); Watson et al. (2015); Audouard et al. (2015) and ARPES Nakayama et al. (2014); Shimojima et al. (2014); Watson et al. (2015), show that the quasiparticles in this material have heavy effective masses exceeding 3 $\sim$ 4 times the free-electron mass, probably due to correlation effects. The FeSe monolayer has similarly large effective masses Liu et al. (2012); Shi et al. (2016). This factor should enhance Zeeman effects in the shallow bands. On the other hand, we are not aware of direct measurements of $g$ factors in iron-based superconductors. In addition, quantitative consideration requires knowledge of the coupling matrix. A challenging practical requirement is the fabrication of a clean monolayer with very small scattering rate.

In the bulk FeSC materials one has to consider three-dimensional electronic spectrum. We expect that the main qualitative features preserve even though the quantum effects are weaker in the 3D case. An additional complication is the possibility of the Fulde-Ferrell-Larkin-Ovchinnikov modulation of the order parameter along the direction of magnetic field for strong spin splitting which has been considered recently in Refs. Gurevich (2010); Takahashi et al. (2014) within the quasiclassical approximation.

Even though our consideration has been motivated by physics of FeSCs, it may be applicable to other multicomponent superconducting systems. Recently, another promising possibility to realize similar LL quantization effects has been discussed for ultra-cold system of two different fermionic atoms with the artificial magnetic fieldÜnal and Oktel (2016). The mathematical description of superconducting instability for this system is very close to multiple-band metals. The reentrant behavior is always expected in this case, since, in contrast to multiple-band superconductors, the Zeeman spin-splitting and disorder effects are absent.

Acknowledgements.

The authors would like to thank Konstantin Matveev and Michael Norman for useful discussions. The work was supported by the U.S. Department of Energy, Office of Science, Materials Sciences and Engineering Division. K. W. S. was supported by the Center for Emergent Superconductivity, an Energy Frontier Research Center funded by the US DOE, Office of Science, under Award No. DEAC0298CH1088

Appendix A Formula for the effective coupling constant $\Lambda_{0,e}$

According to Eq. (10a), $T_{C}$ is directly determined by the effective coupling constants $\Lambda_{0,e}$,

[TABLE]

with $C=2\mathrm{e}^{\gamma_{E}}/\pi\approx 1.134$. The constants $\Lambda_{0,e}$ and $\Lambda_{0,h}$ obey the instability condition $\det(\hat{\Lambda}^{-1}-\hat{\Lambda}^{-1}_{0})=0$, which gives

[TABLE]

Using Eq. (10b) we can exclude $\Lambda_{0,h}$, $\Lambda^{-1}_{0,h}=\Lambda^{-1}_{0,e}/2+\Upsilon_{C}$. This gives the quadratic equation for $\Lambda_{0,e}^{-1}$, which we can solve as

[TABLE]

where, $\delta_{\Lambda}=\pm 1$, $b\!=\!(\Lambda_{ee}\!+\!\tfrac{1}{2}\Lambda_{hh})/\mathcal{D}_{\Lambda}\!-\!\Upsilon_{C}$ and $\mathcal{D}_{\Lambda}=\det\hat{\Lambda}$. The sign $\delta_{\Lambda}$ giving the largest $T_{C}$ (i.e., the largest positive $\Lambda_{0,e}$) has to be selected. Analyzing different cases, we derive $\delta_{\Lambda}=-\mathrm{sign}[(1-\Upsilon_{C}\Lambda_{hh})/\mathcal{D}_{\Lambda}]$. Note that $\Upsilon_{C}$ depends on the ratio $\mu_{h}/T_{C}$ and therefore, formally, Eq. (52) is an implicit equation for $T_{C}$. It is convenient however to treat the ratio $\mu_{h}/T_{C}$ as an independent parameter. In this case Eqs. (50) and (52) determine $T_{C}$ as a function of this parameter, the pairing energy scale $\Omega$, and the coupling matrix.

Appendix B The regularization of the gap equation at finite magnetic field

In this appendix, we describe regularization of logarithmic (UV) divergence in the gap equation for finite magnetic field, Eq. (26), which allows us to absorb all the information about the energy cut-off into $T_{C}$. In the above equation the sums $2\pi T\sum_{\omega_{n}>0}\lambda_{\omega_{n}}^{\alpha}(H,T)$ are logarithmically diverging and have to be cut at the high-energy scale $\Omega$, similar to the zero-field case, Eq. (9). To regularize Eq. (26), we make the standard decomposition

[TABLE]

where, according to Eq. (9), $\lambda_{\omega_{n}}^{e}(0,T)\!=\!1/\omega_{n}$ and $\lambda_{\omega_{n}}^{h}(0,T)\!=\!\left[\tfrac{1}{2}+\eta_{h}(\omega_{n})\right]/\omega_{n}$. As follows from the definitions (10), $2\pi T\sum_{\omega_{n}>0}\lambda_{\omega_{n}}^{\alpha}(0,T_{C})=\Lambda^{-1}_{0,\alpha}$, and this is the only sum containing logarithmic divergence. In other two terms the summation over $\omega_{n}$ can be extended to infinity. In particular, we have

[TABLE]

Using also the definitions of $\mathcal{J}_{\alpha}(H,T)$ in Eq. (28), we obtain the relation $2\pi T\sum_{\omega_{n}>0}\lambda_{\omega_{n}}^{\alpha}(H,T)=\Lambda^{-1}_{0,\alpha}-\mathcal{A}_{\alpha}(T)+\mathcal{J}_{\alpha}(H,T)$ which leads to the regularized gap equation, Eq. (27).

Appendix C Derivation of the kernel eigenvalue $\lambda^{h}_{\omega_{n}}$

In this appendix we present derivation of Eqs. (23) and (24) for the shallow-band eigenvalue $\lambda^{h}_{\omega_{n}}$. Using definition in Eq. (18) with the Green’s functions in Eq. (21) we present $\lambda^{h}_{\omega_{n}}$ as

[TABLE]

where $E_{\ell}=\omega_{c}(\ell+1/2)$, $\zeta_{\omega}=\text{sign}(\omega_{n})$, $x=\rho^{2}/(2l^{2})$, and we used the integral representation $(\omega_{n}\pm\mathrm{i}E_{\ell})=\zeta_{\omega}\int^{\infty}_{0}\mathrm{d}s\exp[-\zeta_{\omega}s(\omega_{n}\pm\mathrm{i}E_{\ell})]$.Zhuravlev and Maniv (2012) The summation of Laguerre polynomials can now be done by using the generating function in Eq. (22) which gives

[TABLE]

where $\vartheta_{1}=\exp(-\mathrm{i}\zeta_{\omega}\omega_{c}s_{1})$ and $\vartheta_{2}=\exp(\mathrm{i}\zeta_{\omega}\omega_{c}s_{2})$. We can now integrate out the variable $x$, and this leads to

[TABLE]

We have introduced the dimensionless quantities $\bar{s}_{1,2}=s_{1,2}\omega_{c}$, $\bar{\omega}_{n}=\omega_{n}/\omega_{c}$, and $\bar{\mu}_{h}=\mu_{h}/\omega_{c}$. This result is equivalent to Eq. (23).

To derive presentation better suited for numerical evaluation, we eliminate the infinite integral using splitting $\int^{\infty}_{0}\mathrm{d}\bar{s}_{1}=\sum^{\infty}_{m=0}\int^{2(m+1)\pi}_{2m\pi}\mathrm{d}\bar{s}_{1}$ and $\int^{\infty}_{0}\mathrm{d}\bar{s}_{2}=\sum^{\infty}_{n=0}\int^{2(n+1)\pi}_{2n\pi}\mathrm{d}\bar{s}_{2}$, and translating the $\bar{s}_{1}\to\bar{s}_{1}+2m\pi$ and $\bar{s}_{2}\to\bar{s}_{2}+2n\pi$ in each term of the sum. Noting that $\vartheta_{1,2}$ remain unchanged by shifting $\bar{s}_{1,2}\to\bar{s}_{1,2}+2n\pi$, this leads to

[TABLE]

in which $Q(\bar{\omega}_{n})$ is determined by the double sum

[TABLE]

Evaluation of this sum gives the quantum oscillating factor

[TABLE]

Next, we change the integration variables in Eq. (54), $\bar{s}=\frac{1}{2}(\bar{s}_{1}+\bar{s}_{2})$ and $\bar{u}=\frac{1}{2}(\bar{s}_{1}-\bar{s}_{2})$, as illustrated in Fig. 9. This gives us the following presentation

[TABLE]

where

[TABLE]

with $\lim_{\bar{s}\rightarrow+0}\mathcal{I}(\bar{s})=1$. The first and second term in Eq. (56) correspond to the integration over the domains $R_{1}$ and $R_{2}$ in Fig. 9. In the new coordinates, the integrand exponentially decays with $\bar{s}$ and oscillates in the $\bar{u}$ direction. Making substitution $\bar{s}\rightarrow 2\pi-\bar{s}$ in the second term and using the relation

[TABLE]

we obtain Eq. (24) of the main text.

Appendix D Representations and asympotics of the function $\mathcal{J}_{1}(H,T)$

In this Appendix we derive presentation for the function $\mathcal{J}_{1}$, Eq. (32), used for numerical calculations. The starting point is the presentation (28a) in which $\lambda^{h}_{\omega_{n}}$ is defined by Eq. (24). Integrating this presentation by parts we obtain

[TABLE]

where the derivative $\mathcal{I}^{\prime}(\bar{s})\equiv d\mathcal{I}(\bar{s})/d\bar{s}$ is explicitly given by Eq. (33). Substituting this result into Eq. (28a), we represent $\mathcal{J}_{1}$ as

[TABLE]

To transform this equation further, we will expand expression under the sum with respect to $\exp\left(-2\pi\bar{\omega}_{n}\right)$, which will allow us to carry out the $\omega_{n}$ summation. Using the relation

[TABLE]

with $\varphi_{\omega_{n}}=\exp[-2\pi(\bar{\omega}_{n}+i\bar{\mu}_{h})]$ and expanding it with respect to $\varphi_{\omega_{n}}$ gives

[TABLE]

Substituting this expansion in Eq. (59) and using Matsubara-frequency summation formula

[TABLE]

we derive

[TABLE]

with $\bar{\tau}\!=\!\pi^{2}T/\omega_{c}$ and $z_{j}^{\pm}=j\!\pm\!\left(1\!-\!\frac{\bar{s}}{\pi}\right)$. Changing the summation variable, $j\rightarrow j-1$, we can transform the first term in the sum as

[TABLE]

which allows us to further simplify $\mathcal{J}_{1}$,

[TABLE]

This presentation is used in numerical calculations.

D.1 Low-temperature limit

To derive the low-temperature asympotics of $\mathcal{J}_{1}$, we start with the intermediate result (61) which we rewrite as

[TABLE]

with

[TABLE]

We can extract the $\ln T$ divergent terms in $P_{T}$ by further rearranging the $j$-summation as follows.

[TABLE]

In the $T\to 0$ limit for $\bar{\mu}_{h}$ not close to half-integers, we may expand $\tanh x\approx x$ and this yields

[TABLE]

Substituting this result into Eq. (63), we obtain Eq. (34).

We note again that the above low-temperature approximation is valid for any $\bar{\mu}_{h}$, except near the half-integers, where the $j$-sum diverges. We can derive more accurate result in the vicinity of $\bar{\mu}_{h}=1/2$ using presentation given by Eq. (62) as a starting point. The last term in this presentation containing the double integration vanishes at $\bar{\mu}_{h}=1/2$ and can be neglected. Defining $\delta_{\mu}=2\bar{\mu}_{h}-1$ with $\delta_{\mu}\ll 1$ we represent $\mathcal{J}_{1}$ as $\mathcal{J}_{1}\approx-\!\Upsilon_{T}+\mathcal{T}$

[TABLE]

For computation of the sum at low temperatures we introduce the intermediate scale $N,$ $1\ll N\ll\omega_{c}/(\pi^{2}T)$ and approximate the sum for $j>N$ by the integral

[TABLE]

where we have used $-\int_{0}^{\infty}dx\cos\left(\pi\delta_{\mu}x\right)\ln\tanh\left(\bar{\tau}x\right)=\frac{1}{2\delta_{\mu}}\tanh\left(\frac{\pi^{2}\delta_{\mu}}{4\bar{\tau}}\right)$. As the sum of the first two terms is not singular at $\delta_{\mu}\rightarrow 0$, we can set $\delta_{\mu}=0$ in them. Evaluating the sum and integral, we obtain

[TABLE]

where $\Gamma(z)$ is the Gamma function, $\Gamma(N+1)=N!$. Using the Stirling’s formula, $\Gamma(z)\approx\mathrm{e}^{-z}z^{z-1/2}\sqrt{2\pi}$ for $z\gg 1$, we obtain

[TABLE]

As expected, the intermediate scale dropped out from the final result. This corresponds to the equation (35) for $\mathcal{J}_{1}$ in the real variables.

D.2 Modified quasiclassical approximation for $\mathcal{J}_{1}$ in the

limit $\omega_{c}\ll\pi T,\mu_{h}$

In the limit $\omega_{c}\ll\pi T,\mu_{h}$ the Landau-quantization effects are very weak allowing for significant simplification of $\mathcal{J}_{1}$. As we consider the case when $\mu_{h}$ is small, we can not simply use the conventional quasiclassical approach but will derive the approximation which also accounts for the case $\mu_{h}\sim T$. We take the presentation for $\mathcal{J}_{1}$ in Eq. (62) as a starting point. In the limit $\bar{\tau}=\pi^{2}T/\omega_{c}\gg 1$ all terms in the sum over $j$ are exponentially small except $z_{j}^{-}$ term for $j=1$ in the $\bar{s}$ integral, $z_{1}^{-}=\bar{s}/\pi$. Also, the $\bar{s}$ integration converges at $\bar{s}\ll 1$ and, therefore, it can be extended to infinity, giving the following approximation

[TABLE]

The function $\mathcal{I}^{\prime}(\bar{s})$, Eq. (33), for $\bar{s}\ll 1$ and $\omega_{c}\ll\mu_{h}$ can be approximated as

[TABLE]

Using also the following presentation for the function $\Upsilon_{T}$

[TABLE]

we derive

[TABLE]

Returning back to the real coordinates, we obtain Eq. (36) of the main text. This result is valid for arbitrary relation between $\mu_{h}$ and $T$ provided $\omega_{c}\ll\pi T,\mu_{h}$. In the limit $\mu_{h}\gg T$ it reproduces the conventional quasiclassical approximation.

In the limit $\omega_{c}\rightarrow 0$ the function $\mathcal{J}_{1}$ vanishes linearly with $\omega_{c}$. In general, the condition $\omega_{c}\ll\pi T,\mu_{h}$ does not yet imply that we can use this linear asympotics. It requires the condition $\omega_{c}\ll T^{2}/\mu_{h}$, which may be stronger because in the limit $T\ll\mu_{h}$ the parameter $T^{2}/\mu_{h}$ is much smaller than both $T$ and $\mu_{h}$.

To compute the ratio $\mathcal{J}_{1}/\omega_{c}$ in the limit $\omega_{c}\rightarrow 0$, we evaluate the $u$ integral as

[TABLE]

Using the integrals $\int_{0}^{\infty}\frac{\sin u}{u}du\!=\!\frac{\pi}{2}$ and $\int_{0}^{\infty}dss\ln\tanh s\!=\!-\frac{7\zeta(3)}{16}$, we transform $\mathcal{J}_{1}$ to the following form

[TABLE]

Here $\zeta(x)$ is the Riemann zeta-function, $\zeta(3)\approx 1.202$. Using substitutions $u=svT/\mu_{h}$ and $s=\tilde{s}/\pi T$, we transform the double integral $\mathcal{L}$ into a single integral as

[TABLE]

The integral in this formula can be transformed as

[TABLE]

giving

[TABLE]

Substituting this result into Eq. (65), we finally obtain asymptotics of $\mathcal{J}_{1}$ for $\omega_{c}\rightarrow 0$

[TABLE]

This result will be used for evaluation of the $H_{C2}$ slope for $T\rightarrow T_{C}$. In the limit $\mu_{h}\gg T$ the second term vanishes and we reproduce well known classical result.

D.3 $\mathcal{J}_{\alpha}(H,T)$ with finite Zeeman splitting

With finite Zeeman splitting Eq. (24) becomes

[TABLE]

Integrating by parts, we can transform it as

[TABLE]

where $\omega_{n}>0$ and for $\omega_{n}<0$, $\lambda^{h}_{\omega_{n}}$ is just the complex conjugate of the above equation. With such $\lambda_{\omega_{n}}^{h}$ the function

[TABLE]

takes the following form

[TABLE]

Using the expansion

[TABLE]

and performing summation over the Matsubara frequencies, we obtain

[TABLE]

with $z_{j}^{\varsigma}=j\!+\!\varsigma\left(1\!-\!\frac{\bar{s}}{\pi}\right)$. Using the relation

[TABLE]

this result can be simplified as

[TABLE]

Taking explicitly the $\bar{s}$ derivative, we obtain Eq. (38a) of the main text.

For $\mathcal{J}_{2}$ with Zeeman effects, we substitute Eq. (37a) into Eq. (28b). To eliminate the $1/\omega_{n}$ in the latter equation, we integrate it by parts which gives

[TABLE]

Using the identity in Eq. (60) and taking the derivative, this yields

[TABLE]

Making change of variable $\bar{s}=\omega_{c}s$, and taking the real part of the above equation, we obtain Eq. (38b) of the main text.

D.4 Presentation for $\mathcal{J}_{1}$ with Landau-level summation

In this appendix we derive an alternative presentation for the function $\mathcal{J}_{1}$ in which the summation over Landau levels is preserved. This presentation is used, e. g., in Refs. Rajagopal and Vasudevan (1966); Gruenberg and Gunther (1968); Maniv et al. (1992). As follows from Eqs. (21) and (18), with finite spin splitting the kernel eigenvalue is

[TABLE]

Using the result

[TABLE]

we can perform integration over the coordinate which yields

[TABLE]

The Matsubara-frequency sum of $\lambda_{\omega_{n}}^{h}$ can be computed using the relation $T\sum_{\omega_{n}}\frac{1}{i\omega_{n}-z}=-\frac{1}{2}\tanh\frac{\beta z}{2}$ and we obtain

[TABLE]

This sum diverges at large $\ell,\ell^{\prime}$ and has to be cut at $\omega_{c}(\ell^{\prime}+\ell)<2\Omega$. To obtain the converging function $\mathcal{J}_{1}$, we have to subtract the zero-field limit of this sum

[TABLE]

which yields

[TABLE]

Introducing new summation variable $m=\ell+\ell^{\prime}$ and making variable change $x=\omega_{c}z/2$ in the integral, we obtain

[TABLE]

Note that there is no divergence when the denominator $\omega_{c}(m\!+\!1)\!-\!2\mu_{h}$ approaches zero, because the nominator $\tanh\frac{\omega_{c}(\ell+\gamma_{z}+\frac{1}{2})-\mu_{h}}{2T}\!+\!\tanh\frac{\omega_{c}(m-\ell-\gamma_{z}+\frac{1}{2})-\mu_{h}}{2T}$ also vanishes. Splitting the integral into the sum of integrals over the intervals $m+\tfrac{1}{2}<z<m+\tfrac{3}{2}$, after some rearrangements we arrive at the presentation given by Eq. (39) of the main text, in which $\mathcal{J}_{1}$ is separated into the converging parts allowing us to take the limit $\Omega\to\infty$.

The presentation (75) can be used to derive the low-temperature behavior of $\mathcal{J}_{1}$ for several interesting particular cases. Without spin splitting and for $\mu_{h}$ matching the Landau level with the index $\ell_{0}$, $\mu_{h}=\omega_{c}(\ell_{0}+1/2)$, the main diverging term at $T\to 0$ in the sum in Eq. (75) is coming from $m=2\ell_{0}$ and $\ell=\ell_{0}$ giving $\mathcal{J}_{1}\simeq(\omega_{c}/8T)(2\ell_{0})!/\left[2^{2\ell_{0}}(\ell_{0}!)^{2}\right]$.

At finite $\gamma_{z}$ such that $2\gamma_{z}$ is not integer, $\mathcal{J}_{1}$ has only logarithmic divergence for $T\to 0$ identical to one in $\mathcal{A}_{1}(T)$. In this case the function $\mathcal{J}_{1}(\omega_{c})$ at $T\rightarrow 0$ has steps when the Zeeman-shifted Landau levels cross the chemical potential at $\omega_{c}=\omega_{c,\ell_{0},\pm}\equiv\mu_{h}/(\ell_{0}\pm\gamma_{z}+\frac{1}{2})$,

[TABLE]

where $+$ ($-$) corresponds to the spin up (spin down) state. For $\gamma_{z}<0.5$ the largest term in this sum is at $m=2\ell_{0}$

[TABLE]

meaning that $\mathcal{J}_{1}(\omega_{c})$ steps down (up) when $\omega_{c}$ crosses $\omega_{c,\ell_{0},+}$ ($\omega_{c,\ell_{0},-}$), see Figs. 2(b) and 2(c). At small $\gamma_{z}$ the zero-temperature value of $\mathcal{J}_{1}-\mathcal{A}_{1}$ is small when $\mu_{h}$ is located between the Zeeman-shifted Landau levels with opposite spins, $|\mu_{h}\!-\omega_{c}(\ell_{0}\!+\tfrac{1}{2})|\!<\!\gamma_{z}\omega_{c}$ reflecting strong pair breaking. It jumps to large values approximately given by $|\Delta\mathcal{J}_{1\pm}|$ when $\mu_{h}$ crosses these levels.

In the case $2\gamma_{z}$ equals integer $j_{z}$ and $\mu_{h}\!=\omega_{c}(\ell_{0}\!+\frac{j_{z}}{2}\!+\tfrac{1}{2})$, the function $\mathcal{J}_{1}(H,T)$ again diverges $\propto 1/T$ at $T\to 0$. The diverging term at $m=2\ell_{0}+j_{z}$ and $\ell=\ell_{0}$ is $\mathcal{J}_{1}\simeq(\omega_{c}/8T)\left(2\ell_{0}+j_{z}\right)!/(2^{2\ell_{0}+j_{z}}\ell_{0}!\left(\ell_{0}+j_{z}\right)!)$. In the important particular case $j_{z}\!=1$ and $\ell_{0}\!=0$, we derived more accurate asymptotics from Eq. (39),

[TABLE]

This result allows us to evaluate the transition temperature for $\omega_{c}=\mu_{h}$.

Appendix E $T_{C2}$ for particular cases

E.1 The lowest Landau-level without spin splitting, $\gamma_{z}=0$ and $\mu_{h}=2\omega_{c}$

The $T_{C2}$ for the lowest LL with zero spin splitting can be derived analytically as follows. First, we calculate $\mathcal{J}_{1}$ at $\omega_{c}=2\mu_{h}$. From Eq. (35), we have

[TABLE]

For $\mathcal{J}_{2}$ in the limit $T\ll T_{C}$, we can expand $\ln\tanh(\pi Ts)\simeq\ln t+\ln(\pi T_{C}s)$ in Eq. (29). This yields

[TABLE]

where $\omega^{e}_{c}=eH/m_{e}c=(m_{h}/m_{e})\omega_{c}$ and

[TABLE]

Therefore, Eq. (30) in the limit of small temperatures becomes

[TABLE]

Solving this equation for $T\equiv T^{(0)}_{C2}$, we obtain Eq. (45) in the main text.

E.2 The case $\gamma_{z}=0.5$ and $\mu_{h}=\omega_{c}$

In this appendix we derive the value of the transition temperature $T_{C2}^{(1)}$ for the free-electron spin splitting $\gamma_{z}=0.5$ when the chemical potential is located at the coinciding spin-up lowest LL and spin-down first LL, $\mu_{h}=\omega_{c}$. Using the result for $\mathcal{J}_{1}$ for these parameters in Eq. (77), we obtain

[TABLE]

For the quasiclassical function $\mathcal{J}_{2}$, Eq. (38b), the low-temperature result in Eq. (79b) has to be modified to account for the finite spin splitting. In the typical case $(\omega_{c}^{e}/\mu)\gamma_{z}^{2}\ll 1$ the paramagnetic effect influences weakly the quasiclassical pairing kernel and can be taken into account perturbatively. Expansion of Eq. (38b) with respect to $\gamma_{z}$ gives $\mathcal{J}_{2}(\gamma_{z})\approx\mathcal{J}_{2}(0)-(\omega_{c}^{e}/\mu)\gamma_{z}^{2}$ meaning that we have to replace $r_{C}$ in Eq. (79b) with $\tilde{r}_{C}=r_{C}[1+2(\omega_{c}^{e}/\mu)\gamma_{z}^{2}]$. Combining the above results, we transform Eq. (30) to the following form

[TABLE]

Solution of this equation for $T=T_{C2}$ gives Eq. (47).