Quasiclassical approach to interacting many-body systems has proved to be a powerful tool in describing their transport and thermodynamic properties. Within this method, the quantum mechanical averages of an operator corresponding to a physical quantity are replaced with the averages of its classical counterpart over all classical trajectories. Alternatively, one can formulate the quasiclassical theory by using the quasiclassical functions which are obtained from the quantum mechanical single-particle propagators by integrating them over all single particle energies. Qualitatively, for a superconductor with pairing gap Δ and quasiparticles with Fermi momentum p_F and Fermi velocity v_F, this procedure corresponds to averaging over the short length scales of the problem ~pF-1 and retaining the physics at long scales ~ v_F/Δ. Quasiclassical theory was particularly useful in the comparatively recent analysis of the problem of far-from-equilibrium order parameter dynamics in charge-neutral superfluids [1–5].

Most recently, several non-trivial phenomena have been observed in a family of iron-based superconductors and their alloys [6]. One example of such phenomena is an observation of the peak in the penetration depth in BaFe₂(As_1−xP_x)₂ as a function of phosphorus concentration [7–10], in Ba_1−xK_xFe₂As₂ as a function of potassium concentration [11] and, most recently in Ba(Fe_1−xCo_x)₂As₂ as a function of cobalt concentration [12]. Another example is the experimental observation of the spin-density-wave order which is characterized by two magnetic ordering vectors, Q→1 and Q→2, in various iron-based superconducting alloys [13–20].

Due to the fact that in iron-based superconductors the superconductivity is often observed near magnetic instability, quasiclassical approaches initially developed for the purely superconducting states have been re-formulated to specifically include the effects of competition between superconducting and magnetic phases as well as the effects of disorder [21–24]. The experimental observations of the peak in the London penetration depth remains only partially understood [25, 26] which provides an additional motivation to look for possible explanations of this effect.

In turn, the experimental discovery of the double-Q→ magnetic state in iron-based superconductors has lead to an appearance of many theoretical works discussing the emergence of this state and its various properties as well as its relation with other magnetic states [27–35]. Most recently, the effects of disorder on the stability of the single- and double-Q→ states have been discussed [36]. In particular, it was found that disorder leads to suppression of the single-Q→ state in favor the double-Q→ one.

Inspired by the earlier work on this problem, in this paper we use a slightly simplified version of the model introduced in reference [36] to formulate a quasiclassical theory of the double-Q→ state in iron-based superconductors. Specifically, we consider the disordered model which incorporates both interband and intraband disorder. In agreement with the earlier results [36], we find that when the interband disorder can be ignored, the intraband disorder promotes the emergence of the double-Q→ state.

This paper is organized as follows. In the next section II introduce the model Hamiltonian. Section III is devoted to the formulation of the quasiclassical approach with the derivation of the quasiclassical equations. In section IV contains the results of the Landau expansion for the free energy using the quasiclassical equations. Section V contains the discussion of the results and comments related to the further development of the presented formalism in the context of the physics of iron-based superconductors. Sections with acknowledgments and Appendix with some technical details conclude the paper.

In what follows we first introduce the model Hamiltonian, which consists of three terms:

The first term on the right hand side of this expression is a single-particle Hamiltonian which describes the band-structure consisting of three bands: hole-like band at the Γ point and two electron-like bands centered at Q→X=(π,0), Q→Y=(0,π) of the two-dimensional Brillouin zone. We use the compact notations to write down H^₀ using the six-component spinor Ψ^k†=(c^k↑†, c^k↓†, d^k↑†, d^k↓†, f^k↑†, f^k↓†):

where σ^0 is a unit 2 × 2 matrix and single particle energy spectra are given by ε_Γ(k) = ξ_k, ξk=ϵ0-k2/2, ε_X(k) = −ξ_k + δ₀ + δ₂ cos 2ϕ_k, ε_Y(k) = −ξ_k + δ₀ − δ₂ cos 2ϕ_k, ϵ₀ is the energy which amounts to the off-set between the bands and k = (k cos ϕ_k, k sin ϕ_k). Here δ₀ is a parameter which is defined relative to the chemical potential μ and describes the deviation from nesting: the bands are perfectly nested when δ₀ = 0. Lastly, δ₂ is an anisotropy parameter which accounts for the ellipticity of the corresponding Fermi pockets [36].

The second term, H^_sdw, appearing in (1) accounts for the spin-density-wave order within the mean-field approximation:

Here m→X, m→Y are the magnetizations corresponding to two structure vectors Q→X and Q→Y. In what follows, we will assume that magnetic state has Ising-like anisotropy, so we replace m→X,Y·σ→→mX,Yσ^3. Within the mean-field approach we have adopted here, the order parameters m_X,Y must be computed self-consistently.

Finally, the last term on the r.h.s. side of Equation (1) introduces the disorder potential in a system. In principle, the disorder should scatter quasiparticles within each band (intraband scattering) as well as between the bands (interband scattering). The disorder unavoidably leads to the suppression of itinerant magnetism. In this paper we will limit ourselves to the case of an intraband disorder only, for an interband disorder scattering only plays a crucial role in the problem of co-existence of magnetism and superconductivity [21–23], while for the problem at hand it will only lead the faster suppression of the magnetic order. Thus, we write for the last term in (1)

and the summation is performed over the impurity sites.

In order to formulate the quasiclassical theory, we first introduce a single-particle correlation function

in the Matsubara representation, Ψ^α(x)=Ψ^α(r,τ), and the averaging is performed over the ground state of the Hamiltonian (Equation 1). Next step consists in employing the equations of motion for the propagator (5):

Here H^_r acts on r, the self-energy part Σ^ is generated by the disorder potential and its action on the propagator is

The summation over the repeated indices is assumed. Next, we perform the Wigner transformation

In the presence of the quenched disorder, propagators will be dependent on R = (r + r′)/2. In what follows we assume that the disorder in uncorrelated and will average the propagator over the disorder distribution which corresponds to self-consistent Born approximation. Lastly, we introduce the following matrices:

Quasiclassical equations can now be derived after we multiply the first equation (6) from the left and the second equation from the right by M^3. Subtracting the resulting first equation from the second one we find

where we introduced the quasiclassical function, [f^,g^] implies the usual commutation relation and

where G^(iω_n, k) is the Matsubara transform of the matrix function G^(τ, k). The self-energy part is determined by the quasiclassical function and disorder scattering rate Γ=πνF|u|2 (ν_F is the density of states at the Fermi level per valley per spin):

In the absence of disorder, quasiclassical equation (10) is linear in G^ and therefore is not sufficient to find G^ unambiguously. In order to define the problem completely, one has to complement (11) with a certain constraint. To derive this constraint, we introduce a new (matrix) function [37]

Equation for this matrix function can be easily derived from (10). It then follows that quasiclassical functions must satisfy the following normalization condition:

As we will demonstrate below, inclusion of disorder potential does not violate this condition. In order to solve the quasiclassical equations (10) self-consistently, we need to specify the matrix structure of the function G^.

To derive an expression for the free energy in terms of the quasiclassical functions, we can employ an expression for the effective action corresponding to the model Hamiltonian (1). Omitting the disorder potential for now, we have [36] F(mX,mY)=(mX2+mY2)/gsdw-S(λ=1) with

Here G^₀(iω_n, k) is the single-particle propagator for the non-interacting system, W^=-mXS^X-mYS^Y and

The expression for the free energy in terms of the quasiclassical functions can be derived by following the steps in the calculation of reference [38]. First, we note

where G^λ(iωn,ϕk) is found from solving the quasiclassical equations (10) in which order parameters have been rescaled by parameter λ, m_X,Y → λm_X,Y. The resulting expression for the free energy reads

This expression can also be employed for the case of non-zero disorder by using the solution of equations (25) with the rescaled magnetizations.

It is a hopeless task to evaluate the free energy (29) analytically, but it is amenable to the numerical analysis. However, our numerical computation of the free energy for the single-Q and double-Q states ran into an unexpected problem: the difference between the free energies of the corresponding states fall within the numerical error of the calculation. Thus, in order to determine which one of the two magnetic states will be energetically favorable, below we derive the Landau expansion.

As we have already pointed out in the Introduction, our main goal was to demonstrate how the quasi-classical method can be applied to analyse the competition between magnetic states in multiband metals in the presence of disorder. Having accomplished that goal, we can now generalize it to investigate the problem of an interplay between superconductivity and magnetism. It is already well-established that by including the interband disorder scattering Anderson-Abrikosov-Gor'kov theorem makes it possible for superconductivity and magnetism to co-exist in a certain region of the phase diagram, which size is determined by the ratio of the intra- and inter-band scattering rates [22, 23]. The question is then would be to check if superconducting order may provide an additional contribution in determining which of the two competing magnetic states would be energetically favorable. These results may be employed to provide a qualitative understanding as to why nematicity has been observed in stoichiometric iron selenide in contrast to electron-doped iron selenide.

Lastly, we would like to mention that the inclusion of the interband disorder scattering would not affect our results in any substantial way. Indeed, compared to the case of intraband disorder, the inclusion of the interband scattering leads primarily to the faster suppression of the critical temperature, without affecting the ground state energies of the single- and double-Q→ states significantly.

To summarize, in this paper we have formulated the quasi-classical approach to analyze the relative stability of the single- and double-Q→ spin-density-wave states with respect to band and effective mass anisotropy as well as disorder scattering. Generally, we find that with an increase in intraband disorder scattering rate, the system favors the single-Q→ for moderately high values of the Fermi surface anisotropy parameter, δ₂.

All datasets generated for this study are included in the article/supplementary material.

MD has performed numerical and analytical calculations. MK has performed the analytical calculations. All authors participated in the discussions of the ideas and the ways to solve the problem, and participated in writing the manuscript. All authors contributed to the article and approved the submitted version.

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.