Eur. Phys. J. B 72, 25–39 (2009) DOI: 10.1140/epjb/e2009-00321-2

THE EUROPEAN PHYSICAL JOURNAL B

Regular Article

The electromagnetic response of a fermion system at zero temperature ´ atek2 P. Tarasewicz1,2, A. Jasionowski2,a , and A. Swi¸ 1 2

Department of Physics, University of Stellenbosch Private Bag X1, Matieland 7602, Stellenbosch, South Africa Faculty of Pharmacy, Collegium Medicum in Bydgoszcz Nicolaus Copernicus University in Toru´ n, ul. Jagiello´ nska 13, 85-067 Bydgoszcz, Poland Received 16 June 2008 / Received in final form 7 February 2009 c EDP Sciences, Societ` Published online 26 September 2009 –  a Italiana di Fisica, Springer-Verlag 2009 Abstract. An investigation is shown whether the Meissner-Ochsenfeld effect exists in a gas of spin 1/2 fermions in which an attraction (denoted by W ) between particles having the same momenta and opposite spins as well as a BCS type four-fermion interaction between pairs are present in a system. The former potential is equivalent to the effective kinetic energy of free fermion pairs. The latter interaction (denoted by V4 ) is responsible for the presence of fermion quadruples in the system. The Meissner effect proves to be weaker than in BCS theory, implying a larger penetration depth λ of the external magnetic field in the regime of the weak attraction W ; however, this effect turns out to exist even when the gap parameter of quadruples vanishes. In this regime the strength of the Meissner effect is the same as in the BCS case. PACS. 74.20.-z Theories and models of superconducting state – 74.20.fg BCS theory and its development

1 Introduction

where T =

Many discoveries of new exotic phases in solids have become the fact in the contemporary physics. Hightemperature superconductivity, heavy-fermion systems, strange normal phases possessing completely different properties from the Fermi-liquid behavior are good examples of such phenomena. Searching for the mechanisms standing behind such phenomena poses a challenge for the theoretical as well as experimental solid state physics. One of such intriguing questions is how many particles are involved in creating some clusters responsible for occurrence of those amazing puzzles of nature. A good example of such a cluster can be a Cooper pair in superconductivity. It is known from [1] that the classical superconductivity is caused by pairs of electrons with opposite momenta and spins. However, we cannot be sure if this concept suffices and a possibility that three or four-particle clusters occur should be taken into account and investigated. Moreover, there can be different pairing channels, not necessary the BCS one. In this paper we investigate a gas of spin 1/2 fermions described in terms of the following Hamiltonian H4 = T + W + V4 , a

e-mail: [email protected]

(1)



ξk nkσ ,

kσ





γk b∗k bk

(2)

and finally a BCS type interaction  V4 = −Λ−1 gkk b∗k b∗−k b−k bk .

(3)

W =

γk nk+ nk− =

k

k

k,k

The operators are defined in the following way: bk = ak+ ak− and nkσ = a∗kσ akσ , operators a∗k± , ak± denote the creation and annihilation operators. gkk and γk represent couplings of the four-fermion interaction and two-particle interaction W respectively, which are assumed to be as follows gkk = −gχ(k)χ(k ), γk = 2γχ(k),

with

with

g > 0,

γ < 0.

The symbol Λ denotes the volume of the system. The factor 2 is adopted for convenience. The function χ(k) denotes the characteristic function of the set S := {k : −δ ≤ 2 2 k ξk ≤ δ} with ξk = 2m − μ standing for the one-electron energy measured from the chemical potential μ. δ is the

26

The European Physical Journal B

typical cutoff of the energy of phonons. The incorporation of the function χ(k) enables us to restrict the range of both interactions to the thin shell of width δ around the chemical potential. The potential W binds fermions with the same momenta and opposite spins into preformed pairs. One can look at the potential W as a kinetic operator of these fermion pairs as well. Such a pair has a quasi-boson character, namely [bk , bk ] = [b∗k , b∗k ] = 0,

for k = k [bk , b∗k ] = 0 ∗ [bk , bk ] = 1 − nk+ − nk− .

The quasi-boson nature of bk manifests itself in the subspace spanned by states unoccupied by particles with k+ and k−. Outside this subspace the nature of these operators deviates from that of real bosons. The interaction W can come from the presence of localized distinguishable spin 1/2 magnetic impurities interacting with conduction electrons via the reduced s–d potential which involves only the z-components of spin operators M Γ2   Vr = − (nk− − nk+ )σαz , N α=1

 σαz

=

k

1 0 0 −1

 ,

where N denotes the number of atoms of the host metal while M is the number of magnetic impurities. Γ 2 is the coupling constant [2]. This interaction is a reduced form of the following interaction Vs−dz = − ×

M Γ2   exp [(k − k )Rα ] N  α=1

k,k ∗ (ak − ak−

− a∗k + ak+ )σαz ,

where the vectors Rα describe the locations of M magnetic impurities in an alloy. The reduction relies on taking into account only the diagonal part of Vs−dz , (k = k ). In paper [2] the thermodynamics of a system with V4 and Vr was investigated. To this end standard Bogolyubov’s methods were used to derive an effective mean field Hamiltonian that contains the term W with the negative coupling constant and other ones. Apart from the quadruple gap there is another parameter describing the ordering of impurities. The mean field Hamiltonian (the effective impurity-electron Hamiltonian) reads [2] h0 (x) =

 k>0

hk (x) − Γ |S|x

M 

1 σαz + |S|d(Γ 2 + N 2 x2 ), 2 α=1

1 hk (x) = (ξk + Γ 2 )(nk+ + nk− + n−k+ + n−k− ) 2 − Γ 2 (nk+ nk− + n−k+ n−k− ) + −2Δk (Bk + Bk∗ ) − Γ xN (nk− − nk+ + n−k− − n−k+ ) + Δk (Bk h0 + Bk∗ h0 ),

where Bk = b−k bk ,

Δk = Λ−1



gkk Bk h0 = Δ,

k

N x = Γ nk− − nk+ + σαz h0 . Δ is the quadruple gap parameter and x describes the impurity subsystem. |S| denotes the number of points in the set S while d = M/N stands for the concentration of impurities and is assumed to be small in the thermodynamic limit. In the grand canonical ensemble the average Ah0 = Tr Ae−βh0 , where A is an arbitrary operator. As one can Tr e−βh0 see the interesting term −Γ 2 (nk+ nk− + n−k+ n−k− ) appears in the Hamiltonian. This term is attractive (in our paper −Γ 2 corresponds to 2γ and its range is restricted to the set S). One of the results was to find a phase with coexistence of superconductivity and ferromagnetism. The calculations concerned diluted alloys and weak coupling between magnetic impurities and conduction electrons. In general, ferromagnetism suppresses the critical temperature and destroys superconductivity but there exist exceptions. Examples of coexistence of superconductivity and ferromagnetism are known in literature. For instance, Rickayzen mentions in [3] about the second-row transition metals doped with small amount of iron. In that case even the enhancement of critical temperature is observed. Of course, in [2] the scenario is different. Conduction electrons are coupled to magnetic impurities and this results in pairing between electrons with the same momenta and opposite spins. However, the resulting interaction competes with the four-fermion superconductivity and finally suppresses quadruples. One can ask if this mechanism could lead to superconductivity resulting from this exotic pairing. This paper provides an answer to this question. But we do not consider the ferromagnetic correlations here. As is widely known the existence of ferromagnetic correlations between impurity spins in a system is destructive for any kind of superconductivity. The Kondo effect is rather more common in these situations. As was already mentioned superconductivity and ferromagnetism can coexist for some range of parameters but after all ferromagnetism wins [2]. It is worthwhile to add that the coexistence of superconductivity and ferromagnetism in a few pseudoternary compounds are reported in [4]. This leads to a very important and still open question regarding the contribution of magnetic phenomena to the pairing. One is certain; the interaction W cannot be derived from the Fr¨ ohlich Hamiltonian. A coupling constant of the effective two-particle interaction depends on a transfer momentum q [5]. In order to obtain an interaction of the density-density form from that potential one has to put q = 0 but at cost of the vanishing of the coupling constant. This shows that the phonon exchange mechanism cannot be a source of W . In this paper we restricted our attention to the interactions W and V4 and neglected the other terms involving the order parameter of impurity subsystem. Thus, we disregarded the ferromagnetic correlations between localized spins of impurities at this stage. Our first goal was to find out what effect of the potential W on quadruples could be. The next step is going to be

P. Tarasewicz et al.: The electromagnetic response of a fermion system at zero temperature

the investigation of the influence of the neglected terms. In our opinion this can significantly affect the electromagnetic kernel. The interaction between magnetic impurities and conduction electrons effectively results in existence of electron pairs with the same momenta and opposite spins. The pairs propagating through a crystal enter the interaction V4 . This interaction, it is V4 , is mediated by the phonon field. At higher temperatures there is a mixture of free fermions and some fraction of free fermion pairs which are responsible for a gap above the critical temperature for quadruples. The fermion pairs enter the four-fermion structures to make a condensate at sufficiently low temperatures. This obviously should lead to superconductivity. However, as one can convince oneself the vanishing of the potential V4 does not imply the vanishing of superconductivity. The potential W has a strong influence on the quadruple order parameter manifesting itself by the attenuation of one. The stronger this potential is the smaller the quadruple gap is. The largest quadruple gap occurs for W = 0. But on the other hand for sufficiently large W and the zero order parameter we still obtain the fulfillment of the Schafroth’s criterion and in consequence the Meissner effect, what is a surprise. Particularly interesting is that the Meissner effect is weaker in the case of W = 0 [6] than that in the case of W strong enough to destroy the quadruple order parameter. As was already mentioned V4 is caused by the phonon field. Such a four-fermion attraction was derived by us by making use of two Fr¨ ohlich type transformations to the Fr¨ ohlich electron-phonon Hamiltonian [5]. The best scenario would take a place if the BCS channel was introduced to Hamiltonian (1). The full system would consist of Cooper pairs, quadruples and pairs with the same momenta and opposite spins. But such a system is hard to investigate and one has to do it gradually. Here, we assumed that the BCS term is somehow suppressed. In the papers [7,8] some results of experiments on 3He in aerogel were interpreted in terms of quadruples. Impurities in such a system might be regarded to be a cause of the destabilization of two-particle correlations like the Cooper correlations but they could strengthen higher order correlations. The idea of preformed fermion pairs was introduced in nineties of the last century, e.g. in [9]. These local fermion pairs were constructed from electrons hopping between sites belonging to a lattice. At sufficiently low temperatures and filling the pairs condensed and constituted a Bose-Einstein condensate. The natural frame of this investigation was the attractive Hubbard model. Here, in this paper, we introduced the k-space counterpart of such pairs but as one can convince oneself these pairs will lead to some interesting phenomenon in an external electromagnetic field. The Hamiltonian (1) was investigated by us in a series of papers, e.g. [10], in which we dealt with the thermodynamics of that Hamiltonian as was already mentioned here. For instance, the existence of pseudogap was shown. As is known the presence of uncorrelated fermion pairs leads to the occurrence of pseudogap in the excitation spectrum. Regarding the form of W we restricted

27

the range of that potential , similarly as V4 , to thin shell around the chemical potential. Moreover, it is assumed to be constant. As was mentioned before the interaction W can effectively be regarded as the kinetic energy of fermion pairs moving through a crystal, that is clear when looking at the latter form in equation (2). Recently, a few authors [11–17] considered the Hubbard model with pairhopping where singlet electron pairs can hop from one lattice site to the other one. Those papers showed an influence of that kind of hopping on the critical temperatures of some high-Tc materials. Generally speaking, the critical temperatures were significantly higher in the presence of that mechanism. Of course, the mechanisms underlying the pair-hopping are rather unknown though the Coulomb repulsion is believed to have something to do with. The former form in equation (2), totally equivalent to the first one, appeared in the papers [18, 19] as a concept leading to the so called statistical spin liquid (SSL) originally introduced in [20]. That state is obtained when positive γk becomes infinite and acts in a whole band. Narrow band systems pose usually a good frame for this scenario. In such a case, the exclusion of states occupied by electrons with the same momenta and opposite spins from momentum space occurs. Because of that effect the pairs with the same momenta and opposite spins cannot exist and the half of electrons is pushed out above the Fermi level of the free electron gas that leads to swelling the momentum space occupied by electrons. It resembles the relation between the Hubbard model and t − J model where the exclusion concerns real space. The underlying mechanism of SSL is the forward scattering process with the singular coupling in which the momentum k of a particle is almost unchanged during scattering process. Such a system differs considerably from the Fermi-liquid state. A derivation of the interaction W with positive γk comes from the Landau’s considerations on the Fermi liquid concept [21]. Among the other things the system undergoes the metalinsulator transition at half filling [18, 19]. It should be added that in the paper [19] authors dealt with both cases: U > 0 and U < 0. There have been some attempts to investigate a system with the BCS interaction in which the normal state was described in terms of SSL [22,23]. The effect of presence of quadruples in an investigated system is rather very weak and can usually be neglected; however, in materials with strong coupling between electrons and phonons such an effect could be visible and recognizable especially from surveys of magnetic flux quanta. The existence of such quanta was reported in [24]. It is widely accepted that if half-h/2e magnetic flux quanta appear among usual ones it points to the existence of quadruples in the investigated system. The model is somewhat unrealistic because the range of both interactions is infinite. In order to adjust this model to reality one should deal with a narrow band system. Thus, one would have to use the cosine dispersion relation instead of the parabolic one. As was mentioned before if the quadruples are under consideration then the BCS term should not be suppressed because this interaction leads indirectly to appearing quadruples in the

28

The European Physical Journal B

system. However, the reasons for doing that have already been given before. The pair hopping should rather be of finite range. Those all simplifications were used to get easier a picture of a four-fermion superconductor’s Meissner phase in the presence of pairs with the same momenta and opposite spins which strongly affect quadruples and such superconductivity. The idea of four-fermion interactions was used to deal with some problems in physics before. On the one hand Kleinert [25], R¨ opke et al. [26] addressed a four-fermion interaction in the nuclear physics. On the other hand Schneider and Keller [27], Kleinert [28] as well as Aligia at al. [29] used this concept in the theory of superconductivity as a possible mechanism of this phenomenon.

2 Interaction with external electromagnetic field at zero temperature The considerations of this section are partly based on perturbation theory used in the book by Abrikosov [30]. We shall focus our interest on the Hamiltonian H asymptotically equivalent to H4 , viz.,   (ξk (nkσ +n−kσ )+γk (nk+ nk− + n−k+ n−k−) H= k>0

σ

− 2Δk (Bk + Bk∗ )) + C =



Hk + C,

and H2 =

 where Bk = b−k bk , Δk = Λ−1 k gkk Bk H = Δ (the gap equation) and C is a constant. The average 2 2 e−βH k Bk H = TrTrBke−βH . Of course ξk = 2m − μ. Since H appears in what follows only in the exponential exp (−βH) of the thermal average, the constant C can be discarded. We assume that γk has the form γk = 2γχ(k). Moreover, γ is assumed to be negative. The function χ(k) allows us to impose the restriction on the range of interactions in the reciprocal space and the factor 2 is introduced for convenience. The sum in the Hamiltonian H is over {k : k > 0} denoting the set of all 1-fermion momenta restricted to a definite half-space of R 3 . Our objective is to demonstrate the existence of the Meissner-Ochsenfeld effect in the system described by the Hamiltonian Hf = H + H  , where H  represents the perturbation due to a weak static external electromagnetic field described by the vector potential A(r). Thus   1 H = d3 rψ ∗ (r) (−i∇ − eA/c)2 2m  − (−i∇)2 ψ(r). (5) After making expansion of the field operators ψ ∗ (r), ψ(r) one obtains H  = H1 + H2 with



a(q)a(q )a∗kσ ak−q−q σ ,

(7)

k,q,q ,σ

where a(k − k) = 1/Λ



d3 rA(r) exp [−i(k − k) · r].

Now we would like to calculate corrections to the ground state energy of the system EG up to second order in a(q). EG is connected with the eigenenergy 2ξk + γk − 2Ek taken from Table A.1 in Appendix A, where Ek = (ξk + 12 γk )2 + Δ2 . To this end, let us replace q with  k − k in H1 . In terms of new creation and annihilation operators a∗ki , aki with i = 1, 2, 3, 4 (defined as a∗k1 := a∗k+ , a∗k2 := a∗k− , a∗k3 := a∗−k+ , a∗k4 := a∗−k− ) H1 assumes the form e   H1 = − a(k −k, t)·(k+k)(a∗k 1 ak1 −a∗k4 ak 4 ) 2mc k>0 k >0

+ a(−k − k, t) · (−k + k)(a∗k 3 ak1 − a∗k4 ak 2 )

+ a(k + k, t) · (k − k)(a∗k 1 ak3 − a∗k2 ak 4 )

 + a(−k + k, t) · (−k − k)(a∗k 3 ak3 − a∗k2 ak 2 ) ,

(4)

k>0

e2 2mc2

where {k : k > 0} denotes the set of all 1-fermion momenta restricted to a definite half-space of R 3 . The expression above is a result of a rearrangement of terms entering (6) and switching k and k in some terms. Thus, we are able to calculate the first order corrections, viz.,

0, k = k; (1) E1 = G|H1 |G = (8) 0, k = k, and (1)

E2

= G|H2 |G e2  = a(q)a(q )G|a∗kσ ak−q−q σ |G 2mc2  kqq σ

e2  2 = 2vk a(q)a(−q), 2mc2

(9)

kq

where the ground state vector |G =

k>0

|Gk , with |Gk 

being the vector 15 in Table A.1 in Appendix A and vk2 = ξk 1 2 (1 − Ek ). In (9) operators are in the original form, it is, with momentum indices belonging to whole momentum (2) space. Now the calculation of the second order term E1 can be performed, viz., (2)

E1

=

 |G|H  |m|2 1 , EG − Em

m=G

e  H1 = − a(q) · (q + 2k)a∗k+qσ akσ , 2mc k,q,σ

(6)

where |m and Em denotes vectors representing excited states constructed from the vectors placed in Appendix A

P. Tarasewicz et al.: The electromagnetic response of a fermion system at zero temperature

and their energies, respectively. EG is the ground state energy. To get this correction the eigenstructure in that Appendix for momenta k and k has to be exploited. It is found that the only nonzero contributions are as follows: |k 1110 |k 0001 |a∗k 1 ak1 |Gk |Gk |2 = (uk vk )2 , |k 1000 |k 0111 | − a∗k4 ak 4 |Gk |Gk |2 = (uk vk )2 , |k 1110 |k 0100 |a∗k 3 ak1 |Gk |Gk |2 = (uk vk )2 , |k 1000 |k 1101 | − a∗k4 ak 2 |Gk |Gk |2 = (uk vk )2 , |k 1011 |k 0001 |a∗k 1 ak3 |Gk |Gk |2 = (uk vk )2 ,

to a vector potential. It yields  e2  2 vk a(q) δE = 2 mc2 kq 

e 2  (uk+q vk )2 + 2mc ξk+q − ξk − 2(Ek + Ek+q ) + γk+q kq  (vk+q uk )2 + ξk − ξk+q − 2(Ek + Ek+q ) + γk  × [a(q) · (2k + q)](2k + q) δa(−q).

|k 0010 |k 0111 | − a∗k2 ak 4 |Gk |Gk |2 = (uk vk )2 , |k 1011 |k 0100 |a∗k 3 ak3 |Gk |Gk |2 = (uk vk )2 , |k 0010 |k 1101 | − a∗k2 ak 2 |Gk |Gk |2 = (uk vk )2 , where u2k = 12 (1+ Eξkk ) and vk2 = 12 (1− Eξkk ). After summing all nonzero contributions we obtain

e 2  = 2mc

Having the expression above for δE the current density can now be obtained, e2  2 vk a(q) mcΛ k  2 e 2  (uk+q vk )2 − cΛ 2m ξk+q − ξk − 2(Ek + Ek+q ) + γk+q k  (vk+q uk )2 + ξk − ξk+q − 2(Ek + Ek+q ) + γk

jq = −2



(uk+q vk )2 ξk+q − ξk − 2(Ek + Ek+q ) + γk+q kq  (vk+q uk )2 + [a(q) · (2k + q)] ξk − ξk+q − 2(Ek + Ek+q ) + γk

(2) E1

× [a(−q) · (2k + q)]. (10)

× [a(q) · (2k + q)](2k + q).

kq

1 c



jδAd3 r,

where j is the current density operator. Averaging this relation over a given state , we obtain 1 δE = − c



jδAd3 r.

(11)

Next we obtain expression (11) in Fourier representation, viz., δE = −

Λ jq δa(−q), c q

with jq =

1 Λ



j0 = −K(0, 0)a(0),

(13)

with

It is well known that the Hamiltonian of a system in an external magnetic field fulfils the following relation [30]: δH = −

(12)

For q tending to zero the relation for current density takes the following form

∗ We  the sum  exploited a (q) = a(−q) here. Furthermore was replaced with the unrestricted sum . k>0 k >0

29

d3 rjeiq·r .

In order to get the expression for jq the corrections to the ground state energy have to be differentiated with respect

e2  2 8 1 K(0, 0) = 2 vk − mcΛ 3 cΛ



2



2(uk vk )2 . 4Ek − γk k k (14) The Meissner effect holds if K(0, 0) = const. > 0 which is known as Schafroth’s criterion [3]. Two zeros standing as the arguments in K(0, 0) correspond to static external field (ω = 0) and q = 0. e 2m

k2

3 The gauge invariance At this stage, one should ask whether the theory is gaugeinvariant. Unfortunately, the present approach does not guarantee this. To show this difficulty it suffices to resort to arguments used by Rickayzen in [3,33]. As is known, only a transverse external field induces the electric current inside a superconductor and Schafroth’s criterion is fulfilled in such a case. A longitudinal field, of course, cannot bring it about and j = 0 for such a field. Here, similarly as in BCS, there is a problem with the equality j = 0 for the longitudinal field. The consequence is that the continuity equation is not satisfied. In BCS theory the problem has been solved but so far we have not been able to cope with it in satisfactory manner due to the great complexity of the issue. As an example let us consider a method

30

The European Physical Journal B

used in [31]. This approach is efficient in the case of the BCS model. The electromagnetic kernel is correctly obtained both at zero and finite temperatures. This method  e is based on the concept of the superflow vs = m ∇ ϕ2 − cm A ϕ that is gauge-invariant. 2 is the phase of the Cooper pairs wavefunction and A is the vector potential. The incorporation of the superflow into a Hamiltonian is done by the following transformation k → k + mvs with the assumption that the lengths of the vectors fulfil the inequality |k|  m|vs |. As was already mentioned the method is successful in the BCS theory but unfortunately this is not so regarding the quadruple superconductivity. Let us shortly demonstrate this. For simplicity we will consider the model with γk = 0. The Hamiltonian (4) reads now  H= [ξk+mvs (nk+mvs + + nk+mvs− ) k>0

+ ξ−k+mvs (n−k+mvs + + n−k+mvs − )  ∗ − 2Δk (Bk+mvs + Bk+mv )] = Hk+mvs . s k>0

The constant C was discarded because plays no role in the subsequent steps. Moreover, we decided to include   e ∇ ϕ4 − cm A. in the indices. Let us add that now vs = m Here ϕ is divided by 4 unlike in the case of the Cooper pairs where the division is by 2. Next we can resort to the linear approximation ξ±k+mvs ≈ ξk ± kvs . The diagonalization of the Hamiltonian leads to the following eigenvalues: ξk + kvs (twofold), ξk − kvs (twofold), 2ξk (fourfold), 2(ξk + kvs ), 2(ξk − kvs ), 3ξk + kvs (twofold), 3ξk − kvs (twofold), 2(ξk + Ek ), 2(ξk − Ek ). Now having found the spectrum one can write the statistical sum down  Z=2 e−2βξk (2 + cosh (2βkvs ) k>0

+ 4 cosh (βkvs ) cosh (βξk ) + cosh (2βEk )). (15) Next we employ the relation for the variation of the free energy density with respect to the external magnetic field δF = jn δA, where jn is the normal particles current density or in other words the paramagnetic contribution to the superconducting current density. This can be calculated from e ∂ ln Z jn = . (16) mβΛ ∂vs The superconducting current density is given by the following expression js = envs − jn ,

(17)

where n is the concentration of electrons in the system. Finally, making use of equations (15–17) along with the linear approximation, we obtain js = envs −

2eβ2  k(kvs ) mΛ k

1 + cosh (βξk ) . × 3 + 4 cosh (βξk ) + cosh (2βEk )

(18)

As one can see the paramagnetic term obtained above differs significantly from the corresponding BCS one. This is due to the difference between the excitation spectra of both systems. In the investigated model the BogolyubovValatin transformation takes the following form [10] αk1 αk2 αk3 αk4

= uk ak1 − vk a∗k2 a∗k3 a∗k4 , = uk ak2 + vk a∗k3 a∗k4 a∗k1 , = uk ak3 − vk a∗k1 a∗k2 a∗k4 , = uk ak4 + vk a∗k1 a∗k2 a∗k3 .

These new operators generate excited states of the system. However, they possess the different form from those known in the BCS theory. It is clear that this fact must lead to some consequences concerning the spectrum. The most important thing is that, unlike in the case of BCS theory, the excitation energies measured from the ground state are not simply additive. For instance, the excitation energy in a two-quasiparticle state measured from the ground state is not the sum of two one-quasiparticle excitation energies. This can be easily verified by making use of the spectrum from Appendix A. Statistics obeyed by the operators α is not strictly that of fermions because E |{αki , α∗ki }|E  equals either 1 or u2k depending on which vector E | from Table A.1 in Appendix A is used. Properties of the spectrum exclude the simple BCS picture with the Fermi-Dirac distribution function for quasiparticles. In this model there are a few such distribution functions that in fact differ from the Fermi-Dirac function. These functions correspond to one, two, three and four-quasiparticle excited states. The interested reader is referred to [10] for more details. Even now we are able to assess the validity of the expression above. If one performs the limit β → ∞ then js = envs . Of course, we obtained the BCS result but this is not the correct finding because the paramagnetic contribution vanishes at zero temperature. However, we know from [6] and this paper as well that this contribution is nonzero at this temperature. To see this it suffices to put γk = 0 in (14). Unfortunately, the reasons of the failure are not known and one can try to find a way to correct the method. Another possible route to finding the solution is as follows. In [3] Rickayzen added some additional term to the reduced BCS Hamiltonian with broken gauge invariance in order to get a correction to the current density in the self-consisted way. This correction vanishes for the transverse field and one is left with the same formula for the current density as in the case of the theory with broken invariance. However, for the longitudinal field the correction cancels the diamagnetic and paramagnetic terms out yielding the zero current density as it should occur for such a field. In our case we limit ourselves to the zero temperature calculations. At present time we are forced to give only an outline of the method which is still under construction. At the very beginning we face the problem with the choice of terms producing the proper correction to the current density. For the sake of simplicity let us deal here with the four-fermion interaction only. The more general four-particle potential written in terms of the second

P. Tarasewicz et al.: The electromagnetic response of a fermion system at zero temperature

quantization language has the following form  V =− g(k1 , k2 , k3 , k4 , k1 , k2 , k3 , k4 )

Following Rickayzen [3] we define new parameters: ˜ |B ∗ ˜ αik1 k2 :=  G ik1 ,k2 + Bik1 ,k2 | G 

k1 k2 k3 k4 k k k k 1 2 3 4

× a∗k1 + a∗k2 − a∗k3 + a∗k4 − ak4 − ak3 + ak2 − ak1 + , that is invariant under change of gauge. We omitted Λ in the formula above and the subsequent ones. After some momenta reduction procedure one obtains the V4 . Besides, there remain other terms with sums over momenta differing from those in V4 . Some of them can be used to construct the appropriate correction V  . Below we give some examples of possible contributions to V  :  ∗ g(k1 , k2 , k1 , k2 )B1k B1k1 ,k2 , − 1 ,k2 k1 k2 k1 k2

and



−

k1 k2 k1 k2

and ˜ |B ∗ ˜ βik1 k2 :=  G ik1 ,k2 − Bik1 ,k2 | G. Having defined quantities above one can use them to rebuild V  performing the following replacement 

−

g(k1 , k2 , k1 , k2 )

k1 k2 k1 k2

−

1 2

∗ ∗ (Bik Bik1 ,k2 + Bik1 ,k2 Bik )= 1 ,k2 1 ,k2  ∗ g(k1 , k2 , k1 , k2 )[αik1 k2 (Bik + Bik1 ,k2 ) 1 ,k2 k1 k2 k1 k2 ∗ − Bik1 ,k2 )], + βik1 k2 (Bik 1 ,k2

∗ g(k1 , k2 , k1 , k2 )B2k B2k1 ,k2 , 1 ,k2

where ∗ : = a∗k1 + a∗k2 − a∗−k2 + a∗−k2 − , B1k 1 ,k2

B1k1 ,k2 : = a−k2 − a−k2 + ak2 − ak1 +

∗ B2k : = a∗k1 − a∗k2 + a∗−k2 + a∗−k2 − , 1 ,k2

B2k1 ,k2 : = ak2 − a−k2 + ak2 − ak1 + . We do not write all possible terms down because there can be many variants of them. The proposed form of the contributions is justified by the fact that the correction to the current density has to depend on the external field. It can be seen when one calculates the action of the operators B ∗ and B on the ground state vector |G. As a result we obtain the vectors being the products of one and threeparticle states similarly as it is in the case of the operators entering H1 . Of course, we will neglect other terms that are not of the above form. Now let us approximate these terms by  ∗ − g(k1 , k2 , k1 , k2 )(Bik Bik1 ,k2 1 ,k2 k1 k2 k1 k2

∗ + Bik1 ,k2 Bik ), 1 ,k2

with i ∈ {1, ..., L}, where L is the number of terms contributing to the correction to the Hamiltonian. The averages are taken over the perturbed ground state vector ˜  = |G + |CF  + |Cc . |G is the unperturbed ground |G state vector, |CF  is the correction due to the external field and |Cc  denotes the correction due to the new terms involved in the Hamiltonian. Both of the latter terms are calculated from  m|H  |G |CF  = |m EG − Em m=G

and |Cc  =

 n|V  |G |n. EG − En

n=G

31

Of course the quantities αik1 k2 and βik1 k2 have to be determined in the self-consisted way. It is clear that they depend on the external field and enter the correction V  . Finally, one has to calculate the current density ˜  + G ˜ |J ˆ p (r)| G ˜ , ˜ |J(r)| ˆ ˜  = G ˜ |J ˆ d (r)| G J(r) =  G G where indices d and p introduce the distinction between diamagnetic and paramagnetic contributions to the current density. Ultimately, the current density reads ˆ d (r)|G + G|J ˆ p (r)|CF  J(r) = G|J ˆ p (r)|G + corr., + CF |J with the correction term ˆ p (r)|Cc  + Cc |J ˆ p (r)|G. corr. = G|J The correction term depends on the quantities αik1 k2 and βik1 k2 . The first three terms in J(r) lead to the expression (12) for γk = 0. For a transverse field the correction to the current density vanishes because the quantities αik1 k2 and βik1 k2 become zero. It is connected to the fact that the equations for αik1 k2 and βik1 k2 contain terms involving the scalar product a(q) · b, where b is a vector parallel to the vector q. In this case the scalar product must be equal to zero and one is left with a system of homogeneous equations. Obviously, their solutions are zero. In this sense, for a transverse field, we obtain exactly the same result described by (14). We believe that for a longitudinal field the suitable choice of terms in V  will provide a sufficiently good correction to the current density yielding the vanishing of this quantity for nonzero q. The method is still under improvements and seems to be the easiest one. Other methods, e.g. Green’s function approach [32] or starting from a phonon-electron Hamiltonian are much more difficult and so far we have not been able to develop them.

32

The European Physical Journal B

4 Meissner effect

get

To consider whether Schafroth’s criterion is fulfilled let us replace the summation over momenta in (14) by integration with respect ξ. We obtain √ e2 2m K(0, 0) = k(0, 0), (19) cπ 2 3 where 2 3/2 1 k(0, 0) = EF − μ3/2 Δ2 3 3

δ

−δ

dξ . E 2 (4E − 2γ)

(20)

n=

1 G| Λ

kσ

2 2 nkσ |G = vk , Λ

(21)

k

n = 43 ρF EF denotes the average system, ρF the density of states

where density of fermions on the Fermi level, in the EF the Fermi energy. The plots of v(ξk )2 representing the probability that a quadruple state is occupied versus one-electron energy ξk are presented in Figure 5. Let us add that the probability that a quadruple state is empty equals u(ξk )2 = 1 − v(ξk )2 . In Figure 5 one can see some evolution of v(ξk )2 from the standard BCS-like function for γ = 0 to a step function for some nonzero γ. This is easy to understand since for increasing γ the parameter Δ becomes smaller and smaller. Then the function vk2 tends to unity. What is important, the chemical potential differs from the Fermi energy in the case of nonzero γ. The symmetry between states occupied above the Fermi sea and those unoccupied below one becomes broken for nonzero W . It can be easily seen when looking at the equality v(ξ)2 = u(−ξ)2 which holds for γ = 0 only. As we shall see this broken symmetry is the cause of the chemical potential shift downward with respect to the Fermi energy of a free electron gas. For γ = 0 the symmetry is restored and the equality μ = EF is valid. Of course, we are aware of the fact that k(0, 0) depends on γ. One needs to remember that the order parameter and the chemical potential are also functions of γ. Thus, in order to calculate k(0, 0) first a system of equations for these quantities has to be solved. After passing to the thermodynamic limit the system takes the following form: ⎧ ⎪ ⎪ ⎨

gρF

2

μ EF

=

⎪ 3/2 ⎪ ⎩ 34 ρF EF = √ √ξ+μ ρF EF

δ −δ ∞ −μ

2

δ−γ = arc sinh δ+γ Δ + arc sinh Δ ⎩ 4 E 3/2 = 4 μ3/2 + √μ((δ−γ)2 +Δ2 − (δ+γ)2 +Δ2 ). 3 F 3 (22) If γ = 0, then the system of equations above takes the following form  Δ = δ/ sinh gρ1F μ = EF , gρF

μ EF

when |γ|  δ we obtain 

In (20) we used the relation 

⎧ ⎨

dξ E

ρ(ξ)(1 −

ξ+γ E )dξ,

which is the density of states for with ρ(ξ) = the free electron gas. The former equation follows √ from (21). We used here the standard approximation ξ + μ ≈ √ μ for ξ ∈ (−δ, δ). After performing the integration and √ making use of 43 (μ − δ)3/2 ≈ 43 μ3/2 − 2δ μ for μ  δ we

Δ=0 √ 2 3/2 = 23 μ3/2 + δ μ. 3 EF

The second equation has two complex solutions and the real one. The real solution is smaller than the Fermi energy and has appropriate physical meaning. Its value equals 0.995 eV for EF = 1 eV and δ = 0.01 eV and is independent of γ. Results of numerical calculations can be seen in Figures 1–4. Computations were performed for weak, intermediate and strong-coupling regimes, that is gρF ∈ [0.1, 10]. In Appendix C the following examples are selected and listed: a. gρF = 0.1 – weak coupling case (wcc); b. gρF = 1 – intermediate coupling case (icc); c. gρF = 10 – strong coupling case (scc). In all calculations EF and δ were assumed to be 1 eV and 0.01 eV, respectively. The chemical potential decreases linearly with increasing |γ|. The behavior of chemical potential is plotted in Figures 3, 4. Concerning the order parameter, it is seen in Figure 1, 2 that the gap vanishes with increasing |γ|. Thus for sufficiently strong interaction W quadruples decay. This effect can be explained by noticing that both of interactions compete with one another. The competence occurs due to different nature of them. If the four-fermion interaction came from a mechanism standing behind the density-density interaction we could expect some stabilization of the quadruple phase with the increase of γ but this is not the case. In the paper [34] we investigated a system with the conventional BCS interaction and the four-fermion one. It was found there that along with the increase of a BCS gap parameter the quadruple gap increases until the saturation of the latter is achieved. This effect happens due to phononic nature of these interactions. According to Introduction in the present paper we deal with the potential W resulting from the presence of magnetic impurities interacting with electrons via the reduced s–d potential for z-components of spin operators. Having the solutions of the equations for Δ and μ at our disposal we can now return to equation (20) and calculate the integral entering this expression. Finally, we

P. Tarasewicz et al.: The electromagnetic response of a fermion system at zero temperature

33

1 1

0.9

0.8 0.998

0.7

0.996 μ

Δ /Δ (0)

0.6

0.5

0.994

0.4

0.3 0.992

g*ρ =0.2 F g*ρF=0.4 g*ρF=0.6 g*ρF=0.8 g*ρF=1.0

0.2

0.1

0

0

g*ρF=0.2 g*ρF=0.4 g*ρF=0.6 g*ρF=0.8 g*ρ =1.0

0.99

F

0.002

0.004

0.006

0.008

0.01

0.012

−0.008

−0.01

−0.012

0.014

|γ|

Fig. 1. The ratio Δ(γ) plotted for weaker couplings. All curves Δ(0) are plotted according to equation (22).

−0.002

0

1.002

g*ρF=2.0 g*ρF=4.0 g*ρF=6.0 g*ρF=8.0 g*ρ =10.0

1

F

0.8

−0.004

Fig. 3. The chemical potential μ(γ) plotted for weaker couplings. All curves are plotted according to equation (22).

1

0.9

−0.006 γ

g*ρF=2.0 g*ρF=4.0 g*ρF=6.0 g*ρF=8.0 g*ρF=10.0

0.998

0.996

0.6 μ

Δ /Δ (0)

0.7

0.5

0.994

0.4 0.992

0.3 0.99

0.2

0.1

0

0.01

0.02

0.03

0.04

0.05 |γ|

0.06

0.07

0.08

0.09

0.1

Fig. 2. The ratio Δ(γ) plotted for stronger coupling. The Δ(0) curves have been calculated according to equation (22).

obtain, that 

k(0, 0) =

2 3/2 1 3/2 Δ 2Δ  EF − μ 3 3 γ 4Δ2 − γ 2   2δ + γ + 2 (γ + δ)2 + Δ2  × arctan 4Δ2 − γ 2   −2δ + γ + 2 (γ − δ)2 + Δ2  − arctan 4Δ2 − γ 2    (γ − δ)2 + Δ2 γ−δ + + arctan Δ Δ    (γ + δ)2 + Δ2 γ+δ − arctan + . (23) Δ Δ

0.988 −0.1

−0.09

−0.08

−0.07

−0.06

−0.05 γ

−0.04

−0.03

−0.02

−0.01

0

Fig. 4. The chemical potential μ(γ) plotted for stronger couplings. All curves are plotted according to equation (22).

We see from the expression above that the electromagnetic kernel is a sum of two terms: diamagnetic and paramagnetic ones. Unlike BCS theory the paramagnetic term is in general non-vanishing in the ground state and weakens the Meissner effect. This should inevitably lead to the deeper penetration of an external magnetic field. The form of k(0, 0) is rather complex and in order to investigate for which γ Schafroth’s criterion holds, generally there is a need to resort to numerical calculation of k(0, 0). However, let us first see what happens when γ tends to zero. There is no problem with the first term in (23) but the second one requires more effort. After a few standard analytical steps one obtains

lim k(0, 0) =

γ→0

δ 2 3/2 1 3/2 E − EF √ > 0. 2 3 F 6 δ + Δ2

(24)

The European Physical Journal B 1

0.68

0.9

0.67

0.8

0.66

0.7

0.65

0.6

0.64 K(0,0)

V2

34

0.5

0.4

0.61

0.2

0.1

0 −0.05

−0.04

−0.03

−0.02

g*ρ =2.0 F g*ρF=4.0 g*ρF=6.0 g*ρF=8.0 g*ρ =10.0

0.6

γ=−0.0001 Δ=0.00851 γ=−0.0011 Δ=0.00847 γ=−0.0050 Δ=0.00784 γ=−0.0130 Δ=0.00100

0.59

F

−0.01

0 ξ

0.01

0.02

0.03

0.04

0.05

2

Fig. 5. The function v(ξ) plotted for gρF = 1 and four values of γ. An evolution from the conventional BCS-like curve to the step function with the step corresponding to the Fermi level is seen. The chemical potential (ξ = 0) is shifted down with respect to the Fermi level.

0.58

0

0.01

0.02

0.03

0.04

0.05 |γ|

0.06

0.07

0.08

0.09

0.1

Fig. 7. The function k(0, 0) plotted for stronger values of the four-particle coupling constant. It has been calculated according to equation (23). Here one can see the same effect of survival of the Meissner phase as in Figure 6.

5 The analysis of the numerical results

0.7

0.68

0.66

0.64

0.62 K(0,0)

γ<0 0.62

0.3

0.6

0.58

0.56

γ<0 g*ρF=0.2 g*ρF=0.4 g*ρF=0.6 g*ρF=0.8 g*ρF=1.0

0.54

0.52

0.5

0.63

0

0.002

0.004

0.006

0.008

0.01

0.012

|γ|

Fig. 6. The function k(0, 0) plotted for weaker values of the four-particle coupling constant. It has been calculated according to equation (23). In all cases one can see the transition to the London limit (the horizontal lines).

In the vicinity of γ = 0 we have ∂Δ ∂γ ≈ 0. This fact has been used to simplify the calculation of k(0, 0). In the weak coupling limit δ  Δ holds, therefore lim k(0, 0) = γ→0

1 3/2 2 EF .

This result is in agreement with the one obtained in [6]. It is obvious that k(0, 0) differs from the BCS case in 3/2 which k(0, 0) = 23 EF . It will be shown that k(0, 0) varies strongly along with increasing |γ|. The values of the order parameter and the chemical potential for γ = 0 in three selected cases are presented in Appendix C.

In this section we would like to analyze the case more deeply. The movement of pairs with negative γ decreases the total energy of the whole system. The calculations were performed for the set of parameters used in the previous section. The weak, intermediate and strong-coupling regimes were investigated. In Figures 6, 7 the electromagnetic kernel is plotted for different values of coupling constant gρF . For all regimes we find k(0, 0) to be an increasing function of |γ|. In both figures one can observe a crossover from weaker to stronger couplings. In Figure 6 k(0, 0) for gρF = 0.2 is almost a step function of |γ| and for larger values of gρF we can see the change of convexity from up in the weak-coupling regime to down in the strong-coupling regime. Moreover, the value k(0, 0) for γ = 0 is shifted from k(0, 0) = 0.5 in the case of gρF = 0.2 to k(0, 0) = 0.65 in the case of gρF = 10 what can be noticed in Figure 6, 7. This can be understood when looking at equation (24). For larger coupling constant the factor √δ2δ+Δ2 becomes smaller due to the larger Δ. This leads to decreasing the paramagnetic term and as a result increasing the electromagnetic kernel. There exists some γc at which the gap parameter Δ vanishes. Moreover, if γ := γ01 = −2Δ occurs, one can suspect the function k(0, 0) to be divergent. In Appendix B we show a proof of non-existence of the singularity in the function k(0, 0) at γ = −2Δ. The numerical calculations showed that the function behaves in continuous way in the vicinity of that point. In order to get a confirmation for the numerical results we investigated analytically the continuity of the electromagnetic kernel at that point. The values of γ = γ01 , γc and chemical potential are presented in Appendix C where some examples were given. γc always fulfills the inequality γc < γ01 . After simple calculations

P. Tarasewicz et al.: The electromagnetic response of a fermion system at zero temperature

we get k(0, 0) =

2 3/2 E 3 F

for

γ ≤ γ01

which is in accordance with the BCS result. It is strikΔk ), ing that putting Δk = 0 in equation (14), (uk vk = 2E k gives the BCS result. Thus, the investigated system starts to behave similarly to the BCS one when |γ| is increased to |γc | and the gap parameter becomes zero. This result means that the external magnetic field penetrates a specimen with the same value of the penetration length as in the BCS case. The larger |γ| causes that the paramagnetic term decreases. The interaction W destroys gradually the four-fermion correlations and simultaneously binds more strongly fermions with the same momenta and opposite spins into pairs. Finally, the bound pairs replace quadruples in the system. The system tends to the BCS value of the electromagnetic kernel with zero paramagnetic term at γ → γc− . For γ < γc superconductivity seems to decay because the order parameter equals zero and there are no four-particle correlations. Despite the vanishing of the order parameter there is still the Meissner effect due to the fulfillment of the Schafroth’s criterion. For such values 3/2 of γ we obtain k(0, 0) = 23 EF that is the London’s result. Moreover, there is another important point, namely, the penetration depth decreases with increasing |γ| to the BCS value and stay at this level even if |γ| > |γc |. This is a surprise because one expects the Meissner effect to decay completely when the order parameter becomes zero. Here, in this system, the Schafroth’s criterion holds even in absence of four-particle correlations. It is obvious that the potential W is responsible for this strange effect. What is the possible interpretation of that? As is known from [35] the paramagnetic term in the BCS electromagnetic kernel disappears because there is a gap in the excitation spectrum. As one can see, for sufficiently small nonzero γ there are two gaps. One of them is connected with the order parameter Δ. The second one can be found as follows. We know from the Table A.1 in Appendix A that the lowest energy is EkG = 2ξk + 2γ − 2Ek . If we calculate the excitation energies from the EkG then making use of the energies in Table A.1 for Δ = 0 and ξk = 0 we obtain lim (ξk −EkG ) = lim (−ξk−2γ+2Ek) = 4|γ|,

Δ→0 ξk →0

Δ→0 ξk →0

lim (2ξk − EkG ) = lim (−2γ + 2Ek ) = 4|γ|,

Δ→0 ξk →0

Δ→0 ξk →0

lim (2ξk + 2γ − EkG ) = lim 2Ek = 2|γ|,

Δ→0 ξk →0

Δ→0 ξk →0

lim (3ξk + 2γ − EkG ) = lim (ξk + 2Ek ) = 2|γ|,

Δ→0 ξk →0

Δ→0 ξk →0

lim (2ξk + 2γ + 2Ek − EkG ) = lim 4Ek = 4|γ|.

Δ→0 ξk →0

Δ→0 ξk →0

It is seen that there is another gap 2|γ| at the chemical potential. As is known the pseudogap exists above the critical temperature in some HTSC systems. However, this gap occurs at T = 0 and has the different nature from Δ. When this gap increases, the order parameter Δ decays

35

gradually and finally vanishes at γc . If the potential W is stronger and stronger then this process destroys the fourfermion superconductivity but the Meissner phase with the London’s kernel is still present in the system. This behavior of the system points to the existence of the long range order due to the presence of fermion pairs with the same momenta and opposite spins. The question is how to determine an order parameter of such electron pairs. This cannot be done by connecting the parameter with the average in the ground state bk G because this is zero. One should notice that the operator W is diagonal in momentum space and the Hartree-Fock scheme nk+ nk− ≈ nk+ G nk− + nk− G nk+ − nk+ G nk− G seems to be not reasonable. The best choice could be as follows. Let us express the new order parameter in terms of the average b∗k bk G counting the number of pairs of electrons with the same momenta and opposite spins in the ground state, namely    ΔW := γk b∗k bk G = γk G|b∗k bk |G = γk vk2 . k

k

k

After integration one yields ΔW = 2ρF δ

     1 (δ−|γ|)2 +Δ2 − (δ+|γ|)2+Δ2 , |γ| < |γc |; |γ| 1− 2δ 2|γ|,

|γ| ≥ |γc |,

with proviso that gρF ≥ 0.1. This function is depicted for three cases in Figure 8. For |γ| < |γc | one can see a slight deviation from the linear function f (|γ|) = 2|γ|. This is due to the non-zero quadruple gap. The deviation is more visible for stronger coupling gρF . The quadruple gap vanishes at |γc | and there is a transition from the quadruple superconductivity to the new superconducting phase with pairs of electrons having the same momenta and opposite spins. For |γ| ≥ |γc | we obtain the linear ΔW function of |γ|. As is seen 2ρ gives the value of the second Fδ gap 2|γ|. This gap vanishes if γ becomes zero. It is worth noting that the Pauli susceptibility tends to zero at T = 0. Therefore, there are no uncompensated spins in the system at that temperature. We believe that if we take into consideration all the terms in the effective impurity-electron Hamiltonian it could provide a significant modification of k(0, 0). One needs to remember that in the full Hamiltonian from Introduction besides Δ and ΔW there is the parameter x corresponding to magnetization and in fact they are coupled to each other. As was already written the function v(ξ)2 representing the probability that a state with four particles is occupied is shown in Figure 5. With increasing |γ| this function tends to the step function as for the free electron gas with the edge at the point corresponding to EF . The electrons occupy the states a∗k+ a∗k− a∗−k+ a∗−k− |0  with the probability equal 1 as it is in the free electron gas but for ξk ∈ [−δ, δ] their energy equals 4ξk + 4γ instead of 4ξk . As one can see the state |Gk  stops to be the lowest energy state vector. Physically, if |γ| is sufficiently large all

36

The European Physical Journal B 0.12 g*ρF = 0.2 g*ρF = 1.0 g*ρ = 6.0 F

0.1

ΔW /2δρF

0.08

0.06

0.04

0.02

0

0

0.01

0.02

0.03 |γ|

0.04

0.05

0.06

ΔW Fig. 8. The parameter 2ρ plotted versus |γ|. Three cases Fδ are presented, namely, for gρF = 0.2, gρF = 1 and gρF = 6.

electrons in the set S are bound into pairs of electrons with the same momenta and opposite spins and the upper boundary of S corresponds to the Fermi sphere for the free electron gas,it is μ+δ = EF . This stems from the fact that the chemical potential was shifted downward and now is inside the set S. The interaction V4 scatters quadruples from state with k to state with k. Both momenta belong to S. To be scattered pairs must have enough room in momentum space. But all states in the set S are occupied and the Pauli’s exclusion principle makes such scattering impossible. The conclusion is that the stronger W is the fewer states are accessible for particles to be scattered by V4 . This points to the dominant role of the Pauli’s principle in promoting the new order. What intrigues is the question about the physical consequences of all of that. There is the ideal diamagnetism along with some order due to the term W in the system. Pairs with the same momenta and opposite spins replaced quadruples; however is this a real superconducting state? Does Nature realize such a scenario?

6 Conclusions The presence of the Meissner-Ochsenfeld effect has been demonstrated by perturbation theory [30] for a superconductor with a BCS-type four-fermion quadruple binding potential and a two-particle interaction W equivalent to the effective kinetic energy of the movement of fermion pairs with the total momentum 2k and total spin 0. The effect proves to be weaker than in BCS, implying a larger penetration depth λ of the external magnetic field, in the regime of |γ| < |γc |, (γc corresponds to Δ = 0). The important fact is that the Meissner effect does not stop to occur for |γ| ≥ |γc | and moreover its strength gets the BCS one at γc . In other words, for negative γ the Meissner effect tends to the BCS result along with decreasing the gap parameter Δ.

The results seem to be strange. It turns out that not only BCS-type attractions lead to Meissner phase and this can be achieved by another type of attractive interaction coming from influence of magnetic impurities. We guess that the problem presented here is worthwhile to be investigated further, especially in the context of the effective impurity-electron Hamiltonian. Another question is how the Meissner Effect manifests itself at finite temperatures. This problem is now being investigated. The incorporation of the gauge-invariance outlined in Section 2 is being improved as well. Finally, the question is how incorporating the BCS interaction affects the Meissner effect in this system. At present it would be very difficult to show it for finite temperatures. It is now under investigation for the ground state only. As suggested in the Introduction such a scenario could occur where half of flux quanta together with usual ones appear, e.g. in YBa2 Cu3 O7−δ [24]. Finally, the positive γ case is worth investigating because this is a strongly correlated electron system. But this case is much more complicated and puzzling due to a different ground state vector and its enormous degeneracy. In such a case double occupancies (k+, k−) are not preferred and are removed from the reciprocal space. The results obtained here point to the Meissner effect as an effect of a more general nature. Further investigations should shed light on its nature. All numerical calculations were done using the Maple package. Authors would like to thank Prof. J. Ma´ckowiak, Mr D. Baran and Mr D. Borycki for fruitful and helpful discussions and suggestions.

Appendix A Hk acts in the 16-dimensional space of states (a∗k1 )n1 (a∗k2 )n2 (a∗k3 )n3 (a∗k4 )n4 |0 , where ni = 0, 1, i = 1, 2, 3, 4. The eigenstructure of Hk [10]: It is worth noting that the vector number 15 is the momentum k. ground state vector |Gk  of the system for The spin projection operator 2Sk = σnαkσ and α=±1 σ=±1

seniorities Λkσ = nkσ − n−k−σ , σ = ±.

Appendix B There are some controversies regarding the point γ = −2Δ, namely, whether the function k(0, 0) is divergent at this point. To address this question let us define a function F (γ) in the following way: F (γ) :=

L(γ) , f (γ)

(B.1)

P. Tarasewicz et al.: The electromagnetic response of a fermion system at zero temperature

37

Table A.1. eigenvector | 1000  | 0100  | 0010  | 0001  | 1010  | 0101  | 1001  | 0110  | 1100  | 0011  | 1110  | 0111  | 1101  | 1011  uk | 0000  + vk | 1111  uk | 1111  − vk | 0000 

1. 2. 3. 4. 5. 6. 7. 8. 9. 10. 11. 12. 13. 14. 15. 16.

where f (γ) :=

eigenvalue Hk 2S ξk 1 –1 ξk ξk 1 –1 ξk 2ξk 2 2ξk –2 2ξk 0 2ξk 0 2ξk + γk 0 2ξk + γk 0 1 3ξk + γk 3ξk + γk –1 3ξk + γk –1 3ξk + γk 1 2ξk + γk − 2Ek 0 2ξk + γk +2Ek 0

 4Δ2 − γ 2 ,

s1 (γ) s2 (γ) − arctan , f (γ) f (γ)  s1 (γ) := 2δ + γ + 2 (γ + δ)2 + Δ2

L(γ) := arctan

and

s2 (γ) := −2δ + γ + 2

 (γ − δ)2 + Δ2 .

These functions enter the expression (B.1) and are crucial for our considerations. It can be seen at once that there is a difficulty in F (γ) at one of the zeros of 4Δ2 − γ 2 , i.e. γ = −2Δ (let it be denoted by γ01 = −2Δ hereafter). We have to know how this function behaves in the vicinity − of the point γ01 . If γ tends to γ01 we encounter the next problem, namely that the imaginary factor i appears in functions f (γ) and L(γ). We can now write the limit as follows i ˜ L(γ), (B.2) lim− F (γ) = lim− γ→γ01 γ→γ01 f˜(γ) where and

 f˜(γ) := −4Δ2 + γ 2 

s1 (γ) ˜ L(γ) := arctan i f˜(γ)





 s2 (γ) − arctan i . f˜(γ)

Λ+ 1 0 0 –1 1 –1 0 0 1 –1 1 –1 0 0 0 0

Λ− 0 1 –1 0 –1 1 0 0 1 –1 0 0 1 –1 0 0

A proof of the first inequality is very simple. For |γ| ≤ δ it is obvious and for |γ| > δ it suffices to notice that  ∀γ < 0, Δ = 0 δ + γ + Δ2 + (γ + δ)2 > 0  holds. However, 2δ + γ + 2 Δ2 + (γ + δ)2 > δ + γ +  Δ2 + (γ + δ)2 is true and shows that (B.3) is fulfilled for all γ < 0. In order to prove (B.4) let us assume that  γ − 2δ + 2 Δ2 + (γ − δ)2 < 0. This can be transformed into  2 Δ2 + (γ − δ)2 < −γ + 2δ. Both sides of the inequality above are positive and can be squared. After a few simple steps one obtains 4Δ2 + 3γ 2 − 4γδ < 0, but it is in contradiction with positiveness of the l.h.s of above inequality (we should remember that γ < 0) what proves (B.4). Having proved these inequalities one can evaluate the limit (B.2). To this end, let us use the following expressions known from complex analysis: arctan z =

To determine the above limit it is helpful to use the following inequalities, viz., ∀γ < 0, Δ = 0

s1 (γ) > 0

(B.3)

∀γ < 0, Δ = 0

s2 (γ) > 0.

(B.4)

1 + iz 1 ln 2i 1 − iz

(B.5)

and ln z = ln |z| + iϕ,

z = |z| exp iϕ.

(B.6)

These expressions turn out to be very useful for treating ˜ L(γ). In the expression concerning ln z one need not take

38

The European Physical Journal B lim F (γ) =

− γ→γ01

1 2δ + x1 (−2Δ) − x2 (−2Δ) , 2 Δ2 − δ 2 + δ(x2 (−2Δ) − x1 (−2Δ)) − Δ(x1 (−2Δ) + x2 (−2Δ)) + x1 (−2Δ)x2 (−2Δ)

(B.8)

f (γ)(2δ + x1 (γ) − x2 (γ)) , 2(Δ2 − δ 2 + δ(x2 (γ) − x1 (γ)) − Δ(x1 (γ) + x2 (γ)) + x1 (γ)x2 (γ))

(B.10)

L(γ) = arctan

into account the ambiguities connected with i2kπ in this ˜ function because they cancel in L(γ). Making use of (B.5) and (B.6) we obtain   f˜(γ) − s (γ)  1   1 ˜ L(γ) = ln   2i f˜(γ) + s1 (γ)   f˜(γ) − s (γ)    2 − ln  (B.7)  + iϕ1 − iϕ2 , f˜(γ) + s2 (γ) where ϕ1 concerns the first logarithm and ϕ2 the second one. In order to determine ϕ1 and ϕ2 we have to inves˜ ˜ 1 2 and ff˜(γ)−s , respectively. To tigate the signs of ff˜(γ)−s (γ)+s1 (γ)+s2 this end one needs to exploit the inequalities (B.3) and (B.4). It can be seen at once with the help of f˜(γ) > 0 that f˜(γ) + s1 > 0 and f˜(γ) + s2 > 0 hold. The inequalities f˜(γ) − s1 < 0 and f˜(γ) − s2 < 0 prove to hold. The proof of these formulas is as follows. Let us start with the former one. It is obvious that f˜(γ) < |γ| thus it suffices to consider the inequality |γ| < s1 (γ). One finds that it is fulfilled for all γ < 0 and implies that f˜(γ) − s1 < 0 holds. The proof of the latter is analogous to that for the former. Having determined the signs of aforementioned in˜ ˜ 1 2 < 0 and ff˜(γ)−s < 0, thus equalities we see that ff˜(γ)−s (γ)+s1 (γ)+s2 ϕ1 = ϕ2 = π. The limit (B.2) can now be calculated and yields See equation (B.8) above,  Δ2 + (δ + γ)2 and x2 (γ) := where x1 (γ) :=  Δ2 + (γ − δ)2 . + we have F (γ) in the Concerning γ tending to γ01 form (B.1). To find this limit we use the following formula arctan x − arctan y = arctan

x−y , 1 + xy

(xy > −1).

(B.9) With the help of the above expression we are able to write down the function L(γ) as follows See equation (B.10) above, It is allowed to pass to the form (B.9) because the condi2 (γ) tion xy = s1 (γ)s > −1 is satisfied. Thus one obtains f (γ)2 lim F (γ) = lim F (γ).

+ γ→γ01

− γ→γ01

Therefore the function k(0, 0) is not divergent at γ01 .

Appendix C Here we would like to present some important numerical results which could be helpful together with the relevant figures during reading this text. For γ = 0 we obtained: gρF Δ[eV] μ[eV] k(0, 0)[eV3/2 ] −6 1. 0.1 0.908 × 10 1 0.5 2. 1 0.00851 1 0.5397 3. 10 0.09983 1 0.65 and for γ < 0 we have: gρF γc [eV] μ[eV] k(0, 0)[eV3/2 ] 1. 0.1 −0.01 0.99 0.666 2. 1 −0.01309 0.99 0.666 3. 10 −0.09983 0.99 0.666 If γ = γ01 we obtained: gρF γ01 [eV] μ[eV] k(0, 0)[eV3/2 ] −6 1. 0.1 −1.816 · 10 1 0.5715 2. 1 −0.01036 0.992 0.6247 3. 10 −0.08925 0.991 0.664

References 1. J. Bardeen, L.N. Cooper, J.R. Schrieffer, Phys. Rev. 108, 1175 (1957) 2. D. Borycki, Eur. Phys. J. B 65, 2938 (2008) 3. G. Rickayzen, Theory of Superconductivity (Interscience, New York, London, 1965) 4. Ø. Fischer, Magnetic Superconductors in Ferromagnetic Materials, edited by K.H.J. Buschow, E.P. Wohlfarth (Elsevier Science Publishers B.V., 1990), Vol. 5 5. P. Tarasewicz, D. Baran, Phys. Rev. B 73, 094524 (2006) 6. P. Tarasewicz, Eur. Phys. J. B 41, 185 (2004) 7. G.E. Volovik, Exotic Properties of Superfluid 3 He (World Scientific, Singapore, 1992) 8. Yu.M. Bunkov, A.S. Chen, D.J. Cousins, H. Godfrin, Phys. Rev. Lett. 85, 3456 (2000) 9. R. Micnas, J. Ranninger, S. Robaszkiewicz, Rev. Mod. Phys. 62, 113 (1990) 10. J. Ma´ckowiak, P. Tarasewicz, Acta Phys. Polonica A 93, 659 (1998) 11. G.I. Japaridze, A.P. Kampf, M. Sekania, P. Kakashvili, Ph. Brune, Phys. Rev. B 65, 014518 (2002)

P. Tarasewicz et al.: The electromagnetic response of a fermion system at zero temperature 12. 13. 14. 15. 16. 17. 18. 19. 20. 21. 22. 23. 24.

G.I. Japaridze, S. Sarkar, Eur. Phys. J. B 28, 139 (2002) S.P. Strong, Phys. Rev. B 51, 16229 (1995) A.E Sikkema, I. Affleck, Phys. Rev. B 52, 10207 (1995) G. Bouzerar, G.I. Japaridze, [cond-mat/9605161] M. Bossche, M. Caffarel, Phys. Rev. B 54, 17414 (1996) S. Robaszkiewicz, B. Bulka, Phys. Rev. B 59, 6430 (1999) K. Byczuk, J. Spalek, Phys. Rev. B, 50 11403 (1994) Y. Hatsugai, M. Kohmoto, J. Phys. Soc. Jpn 61, 2056 (1992) J. Spalek, W. W´ ojcik, Phys. Rev. B 37 1532 (1988) F. Gebhard, The Mott Metal-Insulator Transition (Springer Verlag, Berlin, 1997) J. Czerwonko, Physica C 235–240, 2337(1994) J. Czerwonko, Mol. Phys. Rep. 12, 79 (1995) C.W. Schneider, G. Hammerl, G. Logvenov, T. Kopp, J.R. Kirtley, P.J. Hirschfeld, J. Mannhart, Europhys. Lett. 68, 86 (2004)

39

25. H. Kleinert, J. Phys. G 8, 239 (1982) 26. G. R¨ opke, A. Schnell, P. Schuck, P. Nozieres, Phys. Rev. Lett. 80, 3177 (1998) 27. T. Schneider, H. Keller, Int. J. Mod. Phys. B 8, 487 (1993) 28. A.A. Aligia, A.P. Kampf, J. Mannhart, Phys. Rev. Lett. 94, 247004 (2005) 29. P. Kleinert, [cond-mat.supr-cond 0802.1459v1] 30. A.A. Abrikosov, Fundamentals of the Theory of Metals (North Holland, Amsterdam, 1978) 31. M. Krzyzosiak, R. Gonczarek, Physica C 426–431 273 (2005) 32. J. Schrieffer, Theory of Superconductivity (W.A. Benjamin, Inc., New York, 1964) 33. G. Rickayzen, Green’s Functions and Condensed Matter (Academic Press, New York, London, 1980) 34. P. Tarasewicz, Eur. Phys. J. B 61, 3345 (2008) 35. J. Bardeen, Phys. Rev. 97, 1724 (1955)