# Peculiarities of the orbital effect in the Fulde-Ferrell-Larkin-Ovchinnikov state in quasi-one-dimensional superconductors

Mihail D. Croitoru, Alexandre I. Buzdin

► To cite this version:

Mihail D. Croitoru, Alexandre I. Buzdin. Peculiarities of the orbital effect in the Fulde-Ferrell-Larkin-Ovchinnikov state in quasi-one-dimensional superconductors. Physical Review B: Condensed Matter and Materials Physics (1998-2015), 2014, 89 (22), pp.224506 (1-9). 10.1103/PhysRevB.89.224506. hal-01058631

HAL Id: hal-01058631

https://hal.science/hal-01058631

Submitted on 23 Feb 2016

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
# Peculiarities of the orbital effect in the Fulde-Ferrell-Larkin-Ovchinnikov state in quasi-one-dimensional superconductors

M. D. Croitoru¹,² and A. I. Buzdin¹,³

¹Université Bordeaux I, LOMA, UMR 5798, F-33400 Talence, France

²Departement Fysica, Universiteit Antwerpen, B-2020 Antwerpen, Belgium

³Institut Universitaire de France, Paris, France

(Received 23 April 2014; published 9 June 2014)

Using the quasiclassical formalism, we determine the low-temperature phase diagram of a quasi-one-dimensional superconductor, taking into account the interchain Josephson coupling and the paramagnetic spin splitting. We show that the anisotropy of the onset of superconductivity changes in the FFLO state as compared with the conventional superconducting phase. It can result in anomalous peaks in the field-direction dependence of the upper critical field when the magnetic field length equals to the FFLO period. This regime is characterized by the lock-in effect of the FFLO modulation wave vector, which is governed by the magnetic length. Furthermore, in the FFLO phase, the anisotropy of the upper critical field is inverted at $T_c^{**} = 0.5T_{c0}$, where the orbital anisotropy disappears. We suggest that an experimental study of the anisotropy of the upper critical field can provide very reach information about the parameters of the FFLO phase in quasi-1D samples.

DOI: 10.1103/PhysRevB.89.224506

PACS number(s): 74.70.Kn, 74.78.Fk

## I. INTRODUCTION

Among the numerous studies of the physics of superconductivity in low-dimensional systems the effect of magnetic field has been (stirred) one of the main research subjects (interests). This is because a series of new superconducting compounds with unique (remarkable) properties in magnetic field and pronounced reduced dimensionality has been discovered. Notable examples include the quasi-1D organic Bechgaard salts (TMTSF)₂X, where anion X is PF₆, ClO₄, etc. [1–8], polysulfur nitride (SNₓ) [9], the metal-chalcogenide-based compounds [10–15], transition-metal carbides [16], the quasi-1D M₂Mo₆Se₆ compounds (relatives of the quasi-3D Chevrel phases) [17–19], lithium purple bronze [20], and arrays of 4-Å superconducting carbon nanotubes embedded in the linear pores channels of AFI zeolite single crystals [21,22]. Moreover, the analysis of the behavior of the upper critical field, especially for a field oriented in the highly conducting plane formed by a sheet of highly conducting 1D chains can provide information about the nature of superconductivity in these compounds.

In particular, at low temperature in the compound (TMTSF)₂PF₆, the upper critical field exhibits a pronounced upturn with no sign of saturation. Its magnitude exceeds the Pauli paramagnetic limit, $H_P$, for a magnetic field aligned parallel to their most conducting plane and reaches values $H_{c2} = 90$ kOe [6], which is more than four times larger than $H_P \simeq 22$ kOe. An enhancement of almost two times over $H_P \simeq 27$ kOe is observed in the compound (TMTSF)₂ClO₄, $H_{c2} \simeq 50$ kOe [3]. An array of superconducting 4-Å carbon nanotubes grown in zeolite crystals shows a Pauli limit in a field of $H_P \sim 18$ T, where the resistivity tends to saturate. However, a small finite slope indicates the persistence of superconducting correlations up to 28 T [22]. The upper critical field of purple bronze Li₀.₉Mo₆O₁₇ is also found to exhibit a large anisotropy. The most striking feature is that in the magnetic field aligned along the most conducting axis, $H_{c2}$ increases monotonically with decreasing temperature to a value five times larger than the estimated paramagnetic pair-breaking field ($H_{c2} \lesssim 15$ T, $H_P \simeq 3.1$ T) [20]. A new

transition metal-chalcogenide compound Nb₂Pd₀.₈₁Se₅ with $T_c \sim 6.5$ K displays a remarkably large superconducting upper critical field for fields applied along the most conductive or needle axis. $H_{c2} \simeq 37$ T in this multiband superconductor with nearly cylindrical and quasi-one-dimensional Fermi surface sheets having hole and electron character, significantly exceeds the Pauli limit imposed by quantum statistics on spin singlet superconductors, $H_P \simeq 12$ T, making a very large Maki parameter $\alpha_M = \sqrt{2H_{c2}/H_P} \simeq 4.4$ [13]. The superconducting state in quasi-one-dimensional compounds Ti₂Mo₆Se₆ and In₂Mo₆Se₆ undergoes a normal metal state transition at critical fields approaching the Clogston limit $H_P \sim 4.9$ T and the upper critical field does not saturate down to 50 mK [17,18]. This inorganic compound has great uniaxial anisotropy, comparable to the most anisotropic organic superconductors.

The high crystallographic quasi-one-dimensional structure of these compounds and high values of the Maki parameter convert them into prime candidates for the search of additional superconducting phases such as the Fulde-Farrel-Larkin-Ovchinnikov [23,24] state [25–28] and triplet state [5,6,29,30]. The unusual high critical fields, as well as NMR relaxation and Knight shift results, suggest the possibility that equal spin triplet pairing may be responsible for the superconductivity in the organic compound (TMTSF)₂PF₆, although this statement is not completely settled by the community. The NMR measurements as well as the anomalous in-plane anisotropy of the upper critical field [3,31] in the organic compound (TMTSF)₂ClO₄ favor, for a magnetic field parallel to the conducting chains, the existence of FFLO state with a spin singlet pairing. Moreover, it has been shown within the layered superconductor model under the applied in-plane magnetic field in the FFLO phase that the experimentally observed anomalous field-direction dependence of the upper critical field [3] as well as the anomalous cusps in this dependence due to resonances between the modulation wave vector and the vector potential can serve as a direct evidence for the appearance of the FFLO phase in layered superconductors [32–35]. This peculiar behavior of the in-plane anisotropy of the FFLO phase completes the anomalous oscillatory
out-of-plane angular dependence of the critical field in FFLO state [36–40]. The authors of Ref. [41] have explained the experimentally observed destruction of superconductivity in these compounds in a magnetic field parallel to the most conducting axis in the framework of a triplet superconductivity scenario, provided that the interplane distance is less than the corresponding coherence length. Nevertheless, as indicated in Ref. [20], the experimentally observed highly anisotropic superconducting state in the Q1D purple bronze with a very large Maki parameter can also be attributed to the formation of the FFLO phase.

Motivated by the recent experimental findings of clean quasi-one-dimensional compounds with the orbital pair-breaking effect sufficiently weaker than the Pauli paramagnetic limit, we investigate in this work the possibilities of the spatially modulated superconducting phase formation in quasi-1D superconductors with s-wave pairing. We show that the vector potential of magnetic field modulates the coupling between the conducting chains. The period of this modulation may be in resonance with the FFLO modulation leading to the anomalies in the upper critical field behavior and to a lock-in effect of the FFLO. This study provides a valuable insight into this unconventional state in quasi-1D samples along with a few of its benchmarks (keys findings). Their experimental observation may serve as a direct evidence for the FFLO phase and hence inspire experimentalists for further research with such superconductors.

## II. THE MODEL HIGHLIGHTS

We consider a quasi-one-dimensional (quasi-1D) conductor with the following electron spectrum:

$$E_{\mathbf{p}} = \frac{p_x^2}{2m_x} + 2t_y \cos(p_y d_y) + 2t_z \cos(p_z d_z), \quad (1)$$

where $d_y$ and $d_z$ are the interchain distances along the y and z axes, respectively. We assume that the couplings between chains in both directions are small [see Fig. 1], i.e., $t_z \ll T_{c0}$ and $t_y \ll T_{c0}$, but sufficiently large to suppress the CDW and SDW transitions, to stabilize the superconducting long-range order and to make the mean-field treatment justified, $T_{c0}^2/E_F \ll t_z$, $T_{c0}^2/E_F \ll t_y$. Here, $T_{c0}$ is the critical temperature of the system at $H = 0$.

We choose the gauge for which the vector potential of the external magnetic field is $A = [\mathbf{H} \times \mathbf{r}]$ [$\mathbf{r} = (x,0,0)$] is a coordinate along the x-chain], i.e., $A_z = -xH_y = -xH \sin\alpha \sin\theta$, $A_y = xH_z = xH \cos\alpha \sin\theta$, $A_x = 0$, where H is the amplitude of the magnetic field and $\alpha$ accounts for its direction from the positive z axis in the y-z plane, while $\theta$ is the

<figure><img src="image_2.png" /><figcaption>FIG. 1. (Color online) Scheme of the quasi-1D conductor in a parallel magnetic field.</figcaption></figure>

angle accounting the field direction from the positive x axis. Taking into account that the system is near the second-order phase transition, we can employ the linearized Eilenberger equation for a quasi-1D superconductor in the presence of the magnetic field (in the momentum representation with respect to the coordinate y and z) [42]:

$$(\Omega_n + \hat{\Pi}) f_{\omega}(x, k_{\perp}, p_{\perp}) = \left[ \Delta(x, k_{\perp}) + \frac{\langle f_{\omega} \rangle}{2\tau} \right] \text{sign}(\omega_n) \quad (2)$$

with

$$\hat{\Pi} \equiv \frac{\hbar}{2} v_F \frac{\partial}{\partial x} + 2it_y \sin(p_y d_y) \sin \left( Q_y x - \frac{k_y d_y}{2} \right) + 2it_z \sin(p_z d_z) \sin \left( Q_x x - \frac{k_z d_z}{2} \right),$$

where $k_{\perp} = (k_y, k_z)$, $p_{\perp} = (p_y, p_z)$ and we have assumed that the vector potential varies slowly at the interchain distances (this assumption means that we neglect the diamagnetic screening currents and take the magnetic field as uniform and given by the external field, H), and $Q_y = \pi d_y H_z / \phi_0 = \pi d_y H \cos\alpha \sin\theta / \phi_0$, $Q_x = -\pi d_z H_y / \phi_0 = -\pi d_z H \sin\alpha \sin\theta / \phi_0$ with $\phi_0 = \pi \hbar c / e$, $h = \mu_B H$ is the Zeeman energy, $v_F = v_{F,i}$ is the Fermi velocity along the x axis, and $\Omega_n \equiv \omega_n + 1/2\tau - i\text{hsign}(\omega_n)$. Usually a quasi-1D superconductor can be considered as a system in the clean limit, meaning that the mean free path is much larger than the corresponding intrachain coherence length, $\xi_0^x = v_{F,i}/(2\pi T_{c0})$. The order parameter is defined self-consistently as

$$\frac{1}{\lambda} \Delta(x, k_y, k_z) = 2\pi T \text{Re} \sum_{\omega>0} \langle f_{\omega}(x, p_y, k_y, p_z, k_z) \rangle, \quad (3)$$

where $\lambda$ is the pairing constant and the brackets denote averaging over $p_y$ and $p_z$ and $v_{F,i}$,

$$\langle \dots \rangle \equiv N^{(1)} \int_{-\frac{\pi}{d_y}}^{\frac{\pi}{d_y}} \frac{d_y dp_y}{2\pi} \int_{-\frac{\pi}{d_z}}^{\frac{\pi}{d_z}} \frac{d_z dp_z}{2\pi} \frac{1}{2} \sum_{\pm v_F} (\dots) \quad (4)$$

with $N^{(1)} \equiv m/2\pi p_x d_y d_z$. The temperature unit is so that the Boltzmann constant $k_B = 1$. The maximum critical temperature corresponds to $k_y = 0$ and we have freedom in choosing $k_z$, so we make $k_z = 0$. For the purpose of simplicity, choosing the external magnetic field in the x-y plane, finally we obtain

$$\hat{\Pi} \equiv \frac{\hbar}{2} v_F \frac{\partial}{\partial x} + 2it_z \sin(p_z d_z) \sin(Q_x x). \quad (5)$$

If the angle between the external magnetic field H and the x axis, which is along the chains, is not too small, we can proceed as follows (quantitative analysis will be given later). We seek the solution of the linearized Eilenberger equation, Eq. (2), in the form

$$f_{\omega}(x, y, p_z) = e^{iqr} \sum_{m} e^{imQ \cdot r} f_{m}(\omega_{n}, p_z). \quad (6)$$

Because of the form for $f_{\omega}(x, y, p_z)$ of Eq. (6), one can write $\Delta(x, y)$ as

$$\Delta(x, y) = e^{iqr} \sum_{m} e^{im2Q \cdot r} \Delta_{2m}(k_z). \quad (7)$$

We reserve here the possibility to fix the modulation vector q not only along the x axis. Therefore we keep in the derivation
also its y component. From symmetry considerations, it follows that $\Delta_{-2m} = \Delta_{2m}$. Adopting a second-order approximation in the small parameter $t_z/T_{c0}$ to the solution of Eq. (2) and substituting Eqs. (6) and (7) back into Eq. (2) one gets

$$L_n(\mathbf{q})f_0 + \tilde{t}f_{-1} - \tilde{t}f_1 = \Delta_0, \quad (8)$$

$$L_n(\mathbf{q} \pm 1\mathbf{Q})f_{\pm 1} \pm \tilde{t}f_0 \mp \tilde{t}f_{\pm 2} = 0, \quad (9)$$

$$L_n(\mathbf{q} \pm 2\mathbf{Q})f_{\pm 2} \pm \tilde{t}f_{\pm 1} \mp \tilde{t}f_{\pm 3} = \Delta_{\pm 2}, \quad (10)$$

$$L_n(\mathbf{q} \pm 3\mathbf{Q})f_{\pm 3} \pm \tilde{t}f_{\pm 2} \mp \tilde{t}f_{\pm 4} = 0, \quad (11)$$

$$L_n(\mathbf{q} \pm 4\mathbf{Q})f_{\pm 4} \pm \tilde{t}f_{\pm 3} \mp \tilde{t}f_{\pm 5} = \Delta_{\pm 4}, \quad (12)$$

$$L_n(\mathbf{q} \pm 5\mathbf{Q})f_{\pm 5} \pm \tilde{t}f_{\pm 4} = 0, \quad (13)$$

where $f_m \equiv f_m(\omega_n, p_z)$, $L_n(\mathbf{s}) = \Omega_n \pm i\hbar v_F s_x/2$ [$L_n(\mathbf{q} \pm \mathbf{Q}) = \Omega_n \pm i\hbar v_F(\mathbf{q} \pm \mathbf{Q})_x/2$] and $\tilde{t} = t_z \sin(p_z d_z)$. Here we took into account that $\Delta_{\pm(2m+1)} = 0$ and also $f_{\pm 4}$ and $f_{\pm 5}$ harmonics in order to make description symmetric for the case of resonances (see further discussion).

Solving the system of coupled equations (8)–(11) in a way similar to that described in Ref. [32] and making use of the definition for the order parameter function, Eqs. (3) and (4), we obtain two systems of coupled equations:

$$\begin{aligned}
\Delta_0[P + t_z^2 b_0] &= t_z^2 [g_{\mp 0} \Delta_{\mp 2} + g_{\pm 0} \Delta_{\pm 2}], \\
\Delta_{+2}[P + \delta_{\pm 2} + t_z^2 b_{\pm 2}] &= t_z^2 [g_{\pm 0} \Delta_0 + g_{\pm 2} \Delta_{\pm 4}], \\
\Delta_{-2}[P + \delta_{\mp 2} + t_z^2 b_{\mp 2}] &= t_z^2 g_{\mp 0} \Delta_0, \\
\Delta_{+4}[P + \delta_{\pm 4} + t_z^2 b_{\pm 4}] &= t_z^2 g_{\pm 2} \Delta_{\pm 2},
\end{aligned} \quad (14)$$

where

$$\begin{aligned}
\delta_{\xi n} &= \pi T \sum_{\omega} \left\langle \frac{1}{L_n(q_x)} - \frac{1}{L_n(q_x + \xi n Q_x)} \right\rangle \Big|_{T_{cP}} \\
&= \frac{\pi T}{2} \sum_{\omega, \pm} \frac{1}{B_0^{\pm}} - \frac{1}{B_{4\xi}^{\pm}},
\end{aligned} \quad (15)$$

$$\begin{aligned}
b_{\xi n} &= \pi T \sum_{\omega, \xi} \langle L_n^{-2}(q_x + \xi n Q_x) L_n^{-1}(q_x + \xi n Q_x + \zeta Q_x) \rangle \Big|_{T_{cP}} \\
&= \frac{\pi T}{2} \sum_{\omega, \xi, \pm} \frac{1}{B_{n\xi}^{\pm} [B_{n\xi}^{\pm} + i \hbar v_F \zeta Q_x / 2]},
\end{aligned} \quad (16)$$

$$\begin{aligned}
g_{\xi n} &= \pi T \sum_{\omega} \langle L_n^{-1}(s_x)L_n^{-1}(s_x+\xi Q_x)L_n^{-1}(s_x+\xi 2Q_x) \rangle \Big|_{T_{cP}} \\
&= \frac{\pi T}{2} \sum_{\omega, \pm} \frac{1}{B_0^{\pm} B_{1\xi}^{\pm} B_{2\xi}^{\pm}},
\end{aligned} \quad (17)$$

where $s_x = q_x + \xi n Q_x$, $B_{n\xi}^{\pm} = \Omega_n \pm i\hbar v_F(q_x + \xi n Q_x)/2$ with $\xi = 0, \pm 1$. In Eq. (14), the following notations are introduced: $P = (T_c - T_{cP})/A_c T_c$ and $A_c = 1 - \frac{\hbar}{T_{cP}} \frac{\partial T_{cP}}{\partial \hbar}$. Here, $T_{cP}$ is the superconducting onset temperature in the pure Pauli limit determined by (if takes into account influence of the nonmagnetic impurities)

$$\ln \left( \frac{T_{c0}}{T_{cP}} \right) = \pi T_{cP} \sum_{n} \frac{1}{\omega_n} - \frac{\langle L_n^{-1}(\mathbf{q}) \rangle}{1 - \langle L_n^{-1}(\mathbf{q}) \rangle / 2\tau}, \quad (18)$$

where $\langle L_n^{-1}(\mathbf{q}) \rangle = \frac{1}{2} \sum_{\pm} (\Omega_n \pm i\hbar v_F q_x/2)^{-1}$ and the FFLO modulation vector $\mathbf{q}$ corresponds to the maximum of $T_{cP}(\mathbf{q})$ and its direction is along x axis (see Appendix A). $\Delta T_{cP} = T_c - T_{cP}$ is the shift of the critical temperature due to the orbital effect. The solution of the system of Eqs. (14) is given by

$$T_c = T_{cP} [1 - A_c(S_O + S_R)] \quad (19)$$

with the “orbital” term $S_O \equiv t_z^2 b_0(\mathbf{q}, \mathbf{Q})$ and the “resonance” term $S_R(\mathbf{q}, \mathbf{Q}) = \min_{\xi} S_R^{\xi}(\mathbf{q}, \mathbf{Q})$ with

$$S_R^{\xi}(\mathbf{q}, \mathbf{Q}) = \frac{(\bar{a}_{\xi} - \bar{b}_{\xi})t_z^2 - \delta_{\xi}}{2} - \frac{t_z^2}{2} \sqrt{(\bar{a}_{\xi} - \bar{b}_{\xi} - \delta_{\xi}/t_z^2)^2 + 4g_{1\xi}^2},$$

where

$$\bar{a}_{\xi} = b_0 - \frac{g_{-\xi}^2 t_z^2}{\delta_{-\xi}^2}, \quad \bar{b}_{\xi} = b_{\xi} - \frac{g_{\pm 2}^2 t_z^2}{\delta_{\xi}^4}$$

and those values of $\xi = \pm$ are chosen to maximize the critical temperature. When the system is out of resonance [32], the second harmonic of the order parameter, $\Delta_{\pm 2}$, can be neglected and the solution is just $T_c = T_{cP}/[1 + A_c S_O]$ [34]. However, if $\delta_{\xi} = 0$ then the term $(P + t_z^2 b_{\xi=\pm} + \delta_{\xi})$ in the left-hand side (l.h.s.) of the first line of Eq. (14) is the same as the corresponding term in the second or third line of Eq. (14) and the precise resonance is established with $S_R^{\xi}(\mathbf{q}, \mathbf{Q}) = -c_{\xi}t_z^2$. In this expression, where the FFLO modulation vector is fixed along the x axis, the resonance occurs when $|\mathbf{q} \pm 2\mathbf{Q}| = |\mathbf{q}|$, i.e., $q = \pm \mathbf{Q}$.

Let us consider the limits of validity of our calculations. As was indicated in Ref. [32], our method is valid for $t_z \ll \hbar v_F Q$ or $H \gg t_z \phi_0 / \pi \hbar d v_F$. Since for the magnetic field orientation we have chosen the hopping $t_y$ does not play any role, we get $t_z \ll \hbar v_F Q_x = -\hbar v_F (\pi d H / \phi_0) \sin \theta$ and hence $\sin \theta \gg t_z \phi_0 / \pi \hbar d H v_F$. Therefore the developed formalism is valid only for the field orientation satisfying the following relation:

$$\sin \theta \gg \frac{t_z \phi_0}{\pi \hbar d_z H v_F}. \quad (20)$$

In consequence, all calculations were performed for the external magnetic fields satisfying this condition.

### III. RESULTS AND DISCUSSIONS

In our numerical calculations, we restrict ourselves to the parameters mostly similar to that of Nb$_2$Pd$_{0.67}$Se$_5$ with $d_x = 0.339$ nm, $d_y = 1.283$ nm, and $d_z = 1.539$ nm [12]. The Fermi velocity has been varied in the range $v_F = 0.5 \div 3 \times 10^5$ m/sec, the superconducting temperature was chosen to be $T_{c0} = 6.5$ K [13], and the band couplings are $t_z = 1.6$ K $\sim$ $t_y$. Introducing the dimensionless Fermi velocity parameter, $\eta = \pi \hbar d_z v_F / \phi_0 \mu_B = \hbar v_F Q_x / \mu_B$, this value of $v_F$ corresponds to $\eta = 1.3 \div 7.8$. The summation over the Matsubara frequencies was performed numerically.

Figure 2 illustrates the temperature dependence of the absolute value of the FFLO wave vectors $\mathbf{q}$ (red lines) and the upper critical field in the Pauli limit (green lines) calculated in the clean limit (solid lines) and when moderate disorder is
<figure><img src="image_3.png" /><figcaption>FIG. 2. (Color online) The upper critical field $H_{c2}$ and the absolute value of the FFLO modulation vector $q$ as a function of $T_{cP}/T_{c0}$ for an s-wave quasi-1D superconductor in the pure Pauli limit. Solid lines: $q \neq 0$, $1/\tau = 0$; short dashed lines: $q \neq 0$, $1/\tau \neq 0$; dashed line: $q = 0$.</figcaption></figure>

present in the sample with $\tau k_B T_c / \hbar \simeq 1$ (short-dashed line). The magnitude of the $q$ vector monotonically increases from zero at the tricritical point until infinity in the clean sample or until some large values depending on the mean-free path while accounting for the disorder [43,44]. The optimal direction of the modulation vector is along the conductive chains for $0 < T < T^*$ with $T^* \simeq 0.56 T_{c0}$ for the clean sample. The upper critical field in the inhomogeneous superconducting state diverges at $T \to 0$ in the clean limit and it is appreciably higher than the paramagnetic limit already for $T_{cP}/T_{c0} < 0.3$. At $\tau \neq 0$, $H_{c2}^P$ approaches a finite value. As it was shown in Ref. [45], if the corrugation of the Fermi surface due to the interchain coupling is taken into account, the system looses its 1D character at $T \to 0$ resulting in the finite value of the $H_{c2}(T \to 0)$.

Figure 3 presents the normalized correction of $\Delta T_c = T_c - T_{cP}$ due to the orbital effect as a function of $T_{cP}/T_{c0}$ for different directions of the applied field. In sequel, the tilting of the external field is characterized by $\theta$, the angle the applied field $\mathbf{H}$ makes from the principal $x$ axis of the quasi-1D compound. This figure shows the results for $\eta = 3.16$. The dashed lines illustrate the solution $P = -t_z^2 b_0$, which is justified when the second harmonics of the order parameter is negligible, $\Delta_{\pm 2} = 0$. The orbital effect generates the higher harmonics $q \pm Q$, $q \pm 2Q$, ... as well. For $\hbar v_F Q_x \gg t_z$, we may treat them perturbatively anywhere except for the resonance case $q \approx Q$, when the order parameter harmonics $\Delta_q$ and $\Delta_{q-2Q}$ should be treated on equal footing. The solid lines exhibit the results obtained with the full Eq. (19). At some values of $T_{cP}/T_{c0}$, we see an essential discrepancy between the solid and dashed lines. This discrepancy is induced by the resonance effect discussed above (see Fig. 4) [32]. The physics of the resonance effect in a quasi-1D sample is of the same origin as in the quasi-2D case. As seen from Eq. (5), the

<figure><img src="image_4.png" /><figcaption>FIG. 3. (Color online) Contribution of the orbital effect as a function of $T_{cP}/T_{c0}$ for several magnetic field tilting angles $\theta$ between $\mathbf{H}$ and $q$, when the Fermi velocity is $\eta = 3.16$ and $t/T_{c0} = 0.25$.</figcaption></figure>

projection of the vector potential of the magnetic field on the $x$ axis results in a modulation of the interchain coupling as $t_z \sin(Q_x x)$. The period of this modulation, $\lambda_H = 2\phi_0/d_z H_x$, which depends on the magnitude and direction of the external field, may interfere with the period of the FFLO modulation vector $\lambda_{FFLO}$, leading to the anomalies in $\Delta T_c$, when $\lambda_{FFLO} = \lambda_H$. Figures 5 and 6 show the influence of the interchain coupling on the orbital-motion induced normalized correction of $\Delta T_c = T_c - T_{cP}$. In addition to an overall decrease of the orbital corrections with a decrease of interchain coupling strength $t_z$, one may note also that the resonant cusp becomes narrower.

<figure><img src="image_5.png" /><figcaption>FIG. 4. (Color online) Temperature variation of the magnetic wave vector and its relation to the FFLO modulation wave vector as a function of $T_{cP}/T_{c0}$ for several magnetic field tilting angles $\theta$, when the Fermi velocity is $\eta = 3.16$ and $t/T_{c0} = 0.25$.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>FIG. 5. (Color online) Contribution of the orbital effect as a function of $T_{cP}/T_{c0}$ for several magnetic field tilting angles $\theta$ between **H** and **q**, when the Fermi velocity is $\eta = 3.16$ and $t/T_{c0} = 0.35$.</figcaption></figure>

Figure 7 displays the orbital correction of $\Delta T_c$ as a function of reduced temperature for different magnitudes of the Fermi velocity parameter, $\eta$, when the external field is fixed at $\theta = 45^\circ$. We see that the orbital correction exhibits strong anomalies in the upper critical field in a wide range of the Fermi velocities. As seen from Fig. 2, the curves of the FFLO wave vector **q** and those of the magnetic wave vectors **Q** are almost parallel in a wide range of the reduced temperature. This makes the resonance condition easier to fulfill in practice. For instance, our calculations show that the anomaly in the curve $\Delta T_c$ as a function of $T_{cP}/T_{c0}$ for $\theta = 90^\circ$ occurs in the range $\eta = 1.4 \div 2.2$, for $\theta = 60^\circ$, it takes place in the range $\eta = 1.7 \div 2.5$, and for $\theta = 35^\circ$, in the range $\eta = 2.9 \div 3.7$.

<figure><img src="image_7.png" /><figcaption>FIG. 6. (Color online) Contribution of the orbital effect as a function of $T_{cP}/T_{c0}$ for several magnetic field tilting angles $\theta$ between **H** and **q**, when the Fermi velocity is $\eta = 3.16$ and $t/T_{c0} = 0.1$.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>FIG. 7. (Color online) Contribution of the orbital effect as a function of $T_{cP}/T_{c0}$ for different Fermi velocity parameters in the range $\eta = 1.4–6.7$ and for $\theta = 45^\circ$; $t_z = 1.625$ K, $t/T_{c0} = 0.25$.</figcaption></figure>

The orbital wave vector $Q_x$ is determined by the paramagnetic effect and turns out rather large for $\hbar v_F Q_x \gg \sqrt{t_z T_{c0}}$ for all angles $\theta$ except a narrow vicinity of $\theta = 0$. In such a situation, the angular ($Q$) dependence of the critical temperature is given by the $Q$ expansion of $b_0$ in Eq. (B8) [see Appendix B]. Interestingly, at $T_1^{**} = 0.5 T_{c0}$, the coefficient in front of $Q^2$ in this expansion disappears, making the orbital correction isotropic. Exactly this we observe in Fig. 7, where the curves for different field orientations intersect at the same point $T_1^{**}$. Note that for $T < T_1^{**}$, the orbital effect increases the critical field—the Lebed’s effect ([41]) and a similar situation is realized at $T > T_2^{**} = 0.8 T_{c0}$. However, in the range at $T_1^{**} < T < T_1^{**}$, we have the opposite role of the orbital effect. Therefore we may expect a qualitative change of the critical field anisotropy at $T_2^{**}$ and $T_1^{**}$. Figure 7 illustrates how the anisotropy of the upper critical field becomes two times inverted at these points while sweeping the external conditions of the system from $T_{cP}/T_{c0} > 0.8$ till deep into the FFLO phase. These effects are based on the properties of the FFLO state and hence the experimental observation of the anisotropy inversion may also serve as an additional evidence for the FFLO state.

The existence of the temperature and the Fermi velocity-dependent anomalies in the orbital-motion-induced corrections result in particular features of the anisotropy of the upper critical field in the spatially modulated FFLO phase. Figure 8 shows the magnetic field angular dependence of the normalized superconducting transition temperature, $T_c(\theta)/T_{cP}$, calculated at several values of $T_{cP}/T_{c0}$ for $\eta = 3.16$, when the modulation wave vector is along the **x** axis. For the sake of clarity, the inset of this figure shows the same result but in the polar system. In the polar plot, the direction of each point seen from the origin corresponds to the magnetic field direction, and the distance from the origin corresponds to the normalized critical temperature, when the orbital destructive effect is taken into account. The regions in the close
<figure><img src="image_9.png" /><figcaption>FIG. 8. (Color online) Normalized transition temperature $T_c(\theta)/T_{cP}$ as a function of $\theta$ for several values of $T_{cP}/T_{c0}$ and $\eta = 3.16$; $t_z = 1.625$ K, $t/T_{c0} = 0.25$. Solid lines: $\Delta_{\pm 2} \neq 0$; Dashed lines: $\Delta_{\pm 2} = 0$.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>FIG. 9. (Color online) The absolute value of the FFLO modulation vector $q$ vs $T_{cP}/T_{c0}$ for an s-wave quasi-1D superconductor when accounting for the orbital effects within the first iteration. Here, $t_c/T_c = 0.2$. Solid line takes into account the influence of resonance. Short dashed line shows the field vector $Q$.</figcaption></figure>

vicinity of field orientation angles $\theta = 0, 180$ were excluded, $\sin\theta < t_z\phi_0/\pi\hbar d_z H v_F = t_z/\eta\mu_B H$, because the developed formalism is not applicable, see Eq. (20). The short-dashed lines are extrapolations to correct values. Indeed, for $\theta = 0$, as it is seen from Eq. (3), the orbital influence of the magnetic field is absent and the upper critical field should reach the value of the paramagnetic limit, $T_{cP}$. The solid lines are the results accounting for the contribution of the second harmonics of the order parameter, $\Delta_{\pm 2} \neq 0$, while the dashed ones are obtained when approximating it by $\Delta_{\pm 2} = 0$. In this case, the solution (19) simplifies to $T_c = T_{cP}[1 - A_c S_O]$ and such a solution is valid for $\sqrt{t_z T_{c0}} \ll \hbar v_F Q_x$, or $\sin\theta < \sqrt{t_z T_{c0}}/\eta\mu_B H$, which is the beginning of the superconductivity re-entrant regime [41]. For the intermediate tilting angles of the applied field $t_z/\eta\mu_B H < \sin\theta < \sqrt{t_z T_{c0}}/\eta\mu_B H$, one has to account for the contributions from the higher-order perturbations, $S_R(q, Q)$ in Eq. (19). The second harmonics of the order parameter generate the Lowerence-Doniach term in the final solution [34]. We see that the difference between both solutions is negligible for large tilting angles of the field. However, it is noticeable already for angles smaller than $20^\circ$. So we can infer that the observed knobs or upturn at low angle are due to the Lowerence-Doniach term [34]. Again, we would like to stress one more time that at $\theta = 0$ the upper critical field should become equal to $T_{cP}$ (short-dashed lines). Because the second-harmonics term becomes less important with the field, the knobs are absent for large tilting angles. Nevertheless, in Fig. (8), one can see cusps at some large angles. The deviation between the solutions this time is due to the resonance between the FFLO modulation wave vector and the interchain coupling modulated by the vector potential. Thus, in addition to the overall anisotropy induced by the FFLO modulation, additional cusps develop for certain directions of the applied field, when the resonance conditions are realized. To describe the resonances, we have to account for the second harmonics, $\Delta_{\pm 2}$.

In the previous discussions, we assumed, while investigating the orbital contributions to the upper critical field, that the absolute value of the FFLO wave vector is determined by Eq. (18). However, the orbital contribution, in its turn, may influence the absolute value of the FFLO modulation vector. To investigate this issue, we optimized the full Eq. (19) with respect to the vector $q$. The result of the calculations depicted in Fig. 9, which illustrates the modulus of the FFLO wave vector versus the reduced temperature, shows that the orbital effect on the modulation vector is weak except for the region in the close proximity to the resonance. In the vicinity of the resonance, where the unperturbed $q$ and $Q$ curves intersect, an interesting lock-in effect appears; while sweeping the $T_{cP}/T_{c0}$ across the resonance, the FFLO vector jumps from its almost unperturbed value $q_0$, determined without the orbital contribution, to the magnetic field vector $Q$, crosses the resonance, and then jumps back on the original $q_0$ curve. Further iterations on $q$ and $Q$ will make the resonance observed in Fig. 3 wider due to the lock-in effect.

In conclusion, we have described the behavior of the quasi-1D superconductors in high magnetic field with different tilting with respect to the principal axis of the compound. We demonstrated that (i) the experimentally observed in quasi-1D compounds, such as purple bronze $Li_{0.9}Mo_6O_{17}$, the transition-metal-chalcogenide compound $Nb_2Pd_{0.81}Se_5$, or 4-Å carbon nanotubes grown in zeolite crystals, large values of the upper critical field in the fields parallel to the compound principal axis and at low temperature may be attributed to the stabilization of the FFLO superconducting phase. The low-temperature FFLO phase is characterized by an essential reduction of the magnetic field induced orbital effect as compared to the high-temperature FFLO phase. This makes the $T_c(\theta)/T_{cP}$ dependence much closer to the paramagnetic FFLO limit. The finite magnitude of the paramagnetic limit
may result, for example, from the impurities scattering. (ii) We found that the upper critical field is maximal for the modulation vector **q** parallel to the most conducting direction. (iii) We demonstrate that in quasi-1D compounds one can explore experimentally all field regimes, from the weak field up to the re-entrant one at high magnetic field, by rotating the magnetic field direction [46]. The high-field superconductivity survives because the quasi-1D Fermi surface stabilizes the FFLO state that can exist above the Pauli limited fields. (iv) The resonance between the modulation vector of the FFLO phase and the vector potential of the magnetic field may lead to anomalous cusps in the field-direction dependence of the upper critical field analogous to those calculated previously in quasi-2D conductors. At the resonance, the interplay between the orbital and paramagnetic effects may result in a lock-in effect. (v) We discussed that the inverse effects pertinent to the $T_1^{**} = 0.5T_{c0}$ and $T_2^{**} = 0.8T_{c0}$ may become an additional tool to evidence the FFLO state. We suggest that observation of these effects may serve as a direct proof for the appearance the FFLO phase in quasi-1D superconductors. The possibility to experimentally observe the orbital effect in quasi-1D superconductors should provide very reach information about the parameters of the FFLO phase. In our study, we have assumed s-wave pairing, however, the model can be easily extended to the d-wave pairing [35]. We expect that all the obtained results will remain qualitatively similar.

## ACKNOWLEDGMENTS

We thank D. Jerome for useful discussions. We acknowledge the support by the French ANR program "ElectroVortex" and European NanoSC COST Action MP1201. M.D.C. acknowledges the support by the BELSPO Back to Belgium Grant.

## APPENDIX A: DIRECTION OF THE FFLO WAVE VECTOR

In the above discussion, we explicitly assumed that the FFLO modulation vector is fixed along the **x** axis due to the symmetry of the crystal field [17]. The symmetry of the order parameter can provide additional source of pinning for the modulation vector [35]. The direction of the FFLO vector can be modified by the orbital effects too. Let us investigate the latter case in more details and consider the situation near the tricritical point, where we may expand the Eilenberger equation with respect to $g_x = g_x v_F$, $g_y = 2t_y \sin(\frac{g_y d_y}{2})$ and $g_z = 2t_z \sin(\frac{g_z d_z}{2})$ up to the fourth order, where $g = q + \xi Q$. The result can be written as

$$ \ln \left( \frac{T_{P0}}{T_P} \right) = -\frac{1}{8(2\pi T)^2} \left( g_x^2 + \frac{g_y^2 + g_z^2}{2} \right) \chi^{(2)} \left( \frac{ih}{2\pi T} \right) + \frac{1}{384(2\pi T)^4} \left[ g_x^4 + 3g_x^2(g_y^2 + g_z^2) + \frac{3g_y^2g_z^2}{2} + \frac{3(g_y^4 + g_z^4)}{8} \right] \chi^{(4)} \left( \frac{ih}{2\pi T} \right), $$

where $\chi^{(n)}(ix) \equiv 1/2[\psi^{(n)}(1/2 + ix) + \psi^{(n)}(1/2 - ix)]$. The FFLO modulation vector **q** corresponds to the maximum of $T_{cP}(q)$ and its direction is either along **x** axis or along **y** axis.

If we neglect the orbital contribution then $g = q$. Let us first suppose that $g_x = 0$ and $g_z = 0$ than this equation gives

$$ \ln \left( \frac{T_{P0}}{T_P} \right) = - \frac{\frac{g_y^2}{2}}{8(2\pi T)^2} \chi^{(2)} \left( \frac{ih}{2\pi T} \right) + \frac{\frac{3g_x^4}{8}}{384(2\pi T)^4} \chi^{(4)} \left( \frac{ih}{2\pi T} \right). $$

Extremum of this expression results in

$$ \left. \frac{T_{P0} - T_P}{T_P} \right|_{Q_x=0} = - \frac{\chi^{(2)} \left( \frac{ih}{2\pi T} \right)^2}{\chi^{(4)} \left( \frac{ih}{2\pi T} \right)}. \quad (A1) $$

Next, we suppose that $g_y = 0$ and $g_z = 0$ and this constrain gives

$$ \ln \left( \frac{T_{P0}}{T_P} \right) = - \frac{g_x^2}{8(2\pi T)^2} \chi^{(2)} \left( \frac{ih}{2\pi T} \right) + \frac{g_x^4}{384(2\pi T)^4} \chi^{(4)} \left( \frac{ih}{2\pi T} \right). $$

The extremum of this expression gives rise to

$$ \left. \frac{T_{P0} - T_P}{T_P} \right|_{Q_y=0} = - \frac{3}{2} \frac{\chi^{(2)} \left( \frac{ih}{2\pi T} \right)^2}{\chi^{(4)} \left( \frac{ih}{2\pi T} \right)}. \quad (A2) $$

Near the tricritical point, $h/2\pi T \sim 0.3$. In this case, the right-hand sides (r.h.s.) of these expressions are negative and we can conclude that the direction is fixed along the **x** axis. If the in-plane orbital effects become strong then we recover the result of Ref. [47].

## APPENDIX B: ORBITAL TERM EXPANSION

Below the tricritical point, $q = 0$, if $(\hbar v_F Q_x) \ll T_{c0}$, or $H_x \ll \frac{\phi_0}{\pi \hbar d_x v_F} T_{c0}$, then $P + t^2 b_{\pm} \ll \delta_{\pm}$ and we find from Eq. (14) that

$$ \Delta_{\pm 2} \approx \frac{t^2 c_{\pm}}{\delta_{\pm}} \Delta_0, \quad (B1) $$

with [see Eqs. (15) and (17)]

$$ \delta_{\pm} = \pi T_{cP} \sum_{n} \frac{1}{\Omega_{n}} \left( 1 - \frac{1/2}{1 + iQ} - \frac{1/2}{1 - iQ} \right), \quad (B2) $$

$$ c_{\pm} = \pi T_{cP} \sum_{n, \pm} \frac{1}{2\Omega_n^3} \frac{1}{(1 \pm iQ)(1 \pm iQ/2)}, \quad (B3) $$

where $Q \equiv \hbar v_F Q_x / \Omega_n$. Expansion of these expressions with respect to $g \ll 1$ gives

$$ \delta_{\pm} \approx \pi T_{cP} \sum_{n} \frac{(\hbar v_F Q_x)^2}{\Omega_n^3}, \quad (B4) $$

$$ c_{\pm} \approx \pi T_{cP} \sum_{n} \frac{1}{\Omega_n^3} \left[ 1 - \frac{7(\hbar v_F Q_x)^2}{4\Omega_n^2} \right], \quad (B5) $$

and from Eq. (B1), we find that $\Delta_{\pm 2}$ reads as

$$ \Delta_{\pm 2} \approx \frac{t^2}{(\hbar v_F Q)^2} \Delta_0. \quad (B6) $$
Substitution of $\Delta_{\pm 2}$ back into Eq. (14) leads to the following equation, determining the temperature $T_c$ of the onset of the superconducting state, when the orbital effects of the applied magnetic field are accounted for within the second-order approximation in parameter $t/T_{c0}$,

$$P + t^2 b_0 = \frac{t^4}{(\hbar v_F Q)^2} \sum_{\pm} c_{\pm}, \qquad (B7)$$

where

$$b_0 = \pi T_{cP} \sum_n \frac{1}{2\Omega_n^3} \left( \frac{1}{1+iQ/2} + \frac{1}{1-iQ/2} \right). \qquad (B8)$$

Making use of the expansion of $b_0$ into a series

$$b_0 \approx \pi T_{cP} \sum_n \frac{1}{\Omega_n^3} \left[ 1 - \frac{1}{4} \frac{(\hbar v_F Q)^2}{\Omega_n^2} \right], \qquad (B9)$$

[1] D. Jerome and H. J. Schulz, Adv. Phys. 31, 299 (1982).

[2] A. I. Buzdin and L. N. Bulaevskii, Usp. Fiz. Nauk 144, 415 (1984) [Sov. Phys. Usp. 27, 830 (1984)].

[3] S. Yonezawa, S. Kusaba, Y. Maeno, P. Auban-Senzier, C. Pasquier, K. Bechgaard, and D. Jérome, Phys. Rev. Lett. 100, 117002 (2008).

[4] J. Shinagawa, Y. Kurosaki, F. Zhang, C. Parker, S. E. Brown, D. Jérome, K. Bechgaard, and J. B. Christensen, Phys. Rev. Lett. 98, 147002 (2007).

[5] I. J. Lee, M. J. Naughton, G. M. Danner, and P. M. Chaikin, Phys. Rev. Lett. 78, 3555 (1997).

[6] I. J. Lee, P. M. Chaikin, and M. J. Naughton, Phys. Rev. B 62, R14669 (2000).

[7] I. J. Lee, P. M. Chaikin, and M. J. Naughton, Phys. Rev. B 65, 180502(R) (2002).

[8] J. I. Oh and M. J. Naughton, Phys. Rev. Lett. 92, 067001 (2004).

[9] L. F. Lou, J. Appl. Phys. 66, 979 (1989).

[10] S. Srivastava and B. N. Avasthi, J. Mater. Sci. 27, 3693 (1992).

[11] Y. F. Lu, T. Takayama, A. F. Bangura, Y. Katsura, D. Hashizume, and H. Takagi, J. Phys. Soc. Jpn. 83, 023702 (2014).

[12] Q. R. Zhang, D. Rhodes, B. Zeng, T. Besara, T. Siegrist, M. D. Johannes, and L. Balicas, Phys. Rev. B 88, 024508 (2013).

[13] Q. Zhang, G. Li, D. Rhodes, A. Kiswandhi, T. Besara, B. Zeng, J. Sun, T. Siegrist, M. D. Johannes, and L. Balicas, Sci. Rep. 3, 1446 (2013).

[14] C. Q. Niu, J. H. Yang, Y. K. Li, B. Chen, N. Zhou, J. Chen, L. L. Jiang, B. Chen, X. X. Yang, C. Cao, J. Dai, and X. Xu, Phys. Rev. B 88, 104507 (2013).

[15] D. J. Singh, Phys. Rev. B 88, 174508 (2013).

[16] W. Scherer, Ch. Hauf, M. Presnitz, E.-W. Scheidt, G. Eickerling, V. Eyert, R.-D. Hoffmann, U. C. Rodewald, A. Hammerschmidt, Ch. Vogt, and R. Pöttgen, Angew. Chem. Int. Ed. 49, 1578 (2010); G. Eickerling, Ch. Hauf, E.-W. Scheidt, L. Reichardt, Ch. Schneider, A. Muñoz, S. Lopez-Moreno, A. H. Romero, F. Porcher, G. André, R. Pöttgen, and W. Scherer, Z. Anorg. Allg. Chem. 639, 1985 (2013).

[17] R. Brusetti, A. Briggs, O. Laborde, M. Potel, and P. Gougeon, Phys. Rev. B 49, 8931 (1994).

[18] A. P. Petrović, R. Lortz, G. Santi, M. Decroux, H. Monnard, O. Fischer, L. Boeri, O. K. Andersen, J. Kortus, D. Salloum, P. Gougeon, and M. Potel, Phys. Rev. B 82, 235128 (2010).

we obtain an equation for $T_c$:

$$P = -\pi T_{cP} \sum_n \frac{t^2}{\Omega_n^3} \left[ 1 - \frac{1}{4} \frac{(\hbar v_F Q)^2}{\Omega_n^2} - \frac{t^2}{(\hbar v_F Q)^2} + \frac{7t^2}{4\Omega_n^2} \right]. \qquad (B10)$$

Neglecting the last term, which is beyond our approximation, we get

$$P = \pi T_{cP} \sum_n \frac{t^2}{\Omega_n^3} \left[ -1 + \frac{1}{4} \frac{(\hbar v_F Q)^2}{\Omega_n^2} + \frac{t^2}{(\hbar v_F Q)^2} \right]. \qquad (B11)$$

[19] B. Bergk, A. P. Petrović, Z. Wang, Y. Wang, D. Salloum, P. Gougeon, M. Potel, and R. Lortz, New J. Phys. 13, 103018 (2011).

[20] J.-F. Mercure, A. F. Bangura, Xiaofeng Xu, N. Wakeham, A. Carrington, P. Walmsley, M. Greenblatt, and N. E. Hussey, Phys. Rev. Lett. 108, 187003 (2012).

[21] Z. K. Tang, L. Y. Zhang, N. Wang, X. X. Zhang, G. H. Wen, G. D. Li, J. N. Wang, C. T. Chan, and P. Sheng, Science 292, 2462 (2001).

[22] Z. Wang, W. Shi, H. Xie, T. Zhang, N. Wang, Z. K. Tang, X. X. Zhang, R. Lortz, P. Sheng, I. Sheikin, and A. Demuer, Phys. Rev. B 81, 174530 (2010).

[23] A. I. Larkin and Yu. N. Ovchinnikov, Zh. Eksp. Teor. Phys. 47, 1136 (1964) [Sov. Phys. JETP 20, 762 (1965)].

[24] P. Fulde and R. A. Ferrell, Phys. Rev. 135, A550 (1964).

[25] G. Zwicknagl and J. Wosnitza, in BCS: 50 Years, edited by L. N. Cooper and D. Feldman (World Scientific, Singapore, 2011).

[26] R. Beyer and J. Wosnitza, Low Temp. Physics 39, 225 (2013).

[27] Y. Matsuda and H. Shimahara, J. Phys. Soc. Jpn. 76, 051005 (2007).

[28] A. Buzdin, Physica B: Condensed Matter 407, 1912 (2012).

[29] L. I. Burlachkov, L. P. Gor'kov, and A. G. Lebed, Europhys. Lett. 4, 941 (1987).

[30] A. G. Lebed, Phys. Rev. B 59, R721 (1999); A. G. Lebed, K. Machida, and M. Ozaki, ibid. 62, R795 (2000).

[31] M. D. Croitoru, M. Houzet, and A. I. Buzdin, Phys. Rev. Lett. 108, 207005 (2012).

[32] M. D. Croitoru and A. I. Buzdin, Phys. Rev. B 86, 064507 (2012).

[33] M. D. Croitoru, M. Houzet, and A. I. Buzdin, J. Supercond. Nov. Magn. 25, 1283 (2012).

[34] M. D. Croitoru and A. I. Buzdin, Phys. Rev. B 86, 224508 (2012).

[35] M. D. Croitoru and A. I. Buzdin, J. Phys.: Condens. Matter 25, 125702 (2013).

[36] M. Houzet and A. Buzdin, EuroPhys. Lett. 50, 375 (2000).

[37] A. I. Buzdin and M. Kulic, J. Low Temp. Phys. 54, 203 (1984).

[38] A. I. Buzdin and J. P. Brison, EuroPhys. Lett. 35, 707 (1996).
[39] A. I. Buzdin and J. P. Brison, Phys. Lett. A **218**, 359 (1996).

[40] J. P. Brison, N. Keller, P. Lejay, A. D. Huxley, L. Schmidt, S. Buzdin, N. Bernhoeft, V. P. Mineev, A. N. Stepanov, J. Flouquet, D. Jaccard, S. Julian, and G. Lonzarich, Physica B **199**, 70 (1994).

[41] A. G. Lebed and O. Sepper, Phys. Rev. B **87**, 100511(R) (2013).

[42] N. B. Kopnin, *Theory of Nonequilibrium Superconductivity* (Clarendon Press, Oxford, 2001).

[43] A. I. Buzdin and L. N. Bulaevskii, Pis'ma Zh. Eksp. Teor. Fiz. **31**, 708 (1980) [JETP Letters **31**, 750 (1980)].

[44] A. I. Buzdin and S. V. Polonskii, Zh. Eksp. Teor. Fiz. **93**, 747 (1987) [JETP **66**, 422 (1987)].

[45] N. Miyawaki and H. Shimahara, J. Phys. Scoc. Japan. **83**, 024703 (2014).

[46] N. Dupuis and G. Montambaux, Phys. Rev. B **49**, 8993 (1994).

[47] M. D. Croitoru and A. I. Buzdin, J. Supercond. Nov. Magn. **26**, 1657 (2013).