Screening and anti-screening of the pairing interaction in low-density neutron matter
SScreening and anti-screening of the pairing interaction in low-density neutron matter
S. Ramanan ∗ Department of Physics, Indian Institute of Technology Madras, Chennai - 600036, India
M. Urban † Institut de Physique Nucl´eaire, CNRS-IN2P3, Univ. Paris-Sud,Universit´e Paris-Saclay, 91406 Orsay cedex, France
We study pairing in low-density neutron matter including the screening interaction due to theexchange of particle-hole and RPA excitations. As bare force we employ the effective low-momentuminteraction V low k , while the Fermi-liquid parameters are taken from a phenomenological energydensity functional (SLy4) which correctly reproduces the equation of state of neutron matter. Atlow density, we find screening, i.e., pairing is reduced, while at higher densities, we find anti-screening, i.e., pairing is enhanced. This enhancement is mostly due to the strongly attractiveLandau parameter f . We discuss in detail the critical temperature T c in the limit of low densitiesand show that the suppression of T c predicted by Gor’kov and Melik-Barkhudarov can only bereproduced if the cutoff of the V low k interaction is scaled with the Fermi momentum. We alsodiscuss the effect of non-condensed pairs on the density dependence of T c in the framework of theNozi`eres-Schmitt-Rink theory. I. INTRODUCTION
Neutron stars provide a unique laboratory with an in-terplay of a wide range of phenomena. The physics ofthe inner crust of neutron stars, where a dilute gas ofunbound neutrons coexists with nuclear clusters, is par-ticularly interesting [1]. In this work, we focus on theneutron gas, since its superfluid properties are crucialfor the understanding of astrophysical observables suchas pulsar glitches or neutron-star cooling. Glitches arethe observed sudden increase in the rotational frequencyof the pulsars, followed by a long relaxation time andusually they are linked to the neutron superfluidity inthe inner crust [2–4], in particular, to the unpinning ofthe vortices. After the initial rapid cooling via neutrinoemissions, the cooling rate of the neutron star is very de-pendent on the physics of the crust. The superfluidity ofthe neutrons in the crust of the star strongly suppressesthe specific heat and hence influences the cooling rate[5, 6]. In addition, neutron superfluidity allows for novelneutrino emission processes via Cooper pair breaking andformation that affect the cooling rate of the star close tothe transition temperature [7].Even the modelling of uniform matter is theoreticallyvery challenging due to the uncertainties in the nuclearinteractions. In neutron stars, the attractive interactionis provided by the two-body interaction, and the mostimportant channels for neutron pairing turn out to bethe S channel at low densities and therefore occuringin the inner crust, while in the core, the neutrons pairin the triplet P − F channel. Protons can also pair,although a description of proton superfluidity is compli-cated by the asymmetry of matter and the resulting cou- ∗ Electronic address: [email protected] † Electronic address: [email protected] pling of the protons to the denser background [8]. Inaddition to being crucial for the physics of neutron stars,pairing between nucleons plays a very important role inthe spectra of finite nuclei, as well as in description ofneutron rich nuclei close to the drip line.A reliable description of pairing at all densities in in-finite matter is still an open question, although the su-perfluidity in stars has been studied since the early workof Migdal [9] and Ginzburg and Kirzhnits [10, 11] andis needed to explain observations such as the long relax-ation time after a glitch [12]. For a recent review, thereader is referred to [13]. The simplest starting pointfor the study of pairing is the superfluid gap equationwithin the BCS approximation that uses the free-spacetwo-nucleon interaction as input and a free spectrum forthe single-particle energies. However, there is enoughevidence that one needs to go beyond this approxima-tion [8, 14–19]. Medium corrections to the single-particleenergy and to the free-space interaction change the gapdrastically.In this work, we re-visit the issue of building an in-duced interaction that will modify the free space two-body interaction responsible for pairing in the S chan-nel in uniform neutron matter. In the past, several at-tempts have been made to include medium correctionsto the interaction [20–24]. Most of these calculations usemany-body methods analogous to the well-known exam-ple of screening in an electron gas [25], subject to variousapproximations. Because of the exponential dependenceof the gap on the interaction, the final results are alwaysaffected by the details. In view of the persistent un-certainties, of some mistakes in Ref. [23] (see Ref. [24]),and of the simplifying approximation made in Ref. [24]to replace the 3 particle − a r X i v : . [ nu c l - t h ] A ug interaction, V low k , evolved from the AV two-body po-tential. The same interaction is also used for the 3p1hcouplings entering the induced interaction. The main ad-vantage of using V low k is that non-perturbative featurespresent in the bare interaction, such as the short-rangerepulsion that arises from the hard-core and the repul-sive tensor, are softened. The low-momentum effectiveinteraction depends on the renormalization scale (or cut-off) Λ, while the free-space two-body observables suchas scattering phase shifts and energies are independentof Λ. However, in principle the renormalization grouprunning generates also induced three- and higher-bodyforces. In addition, in a many-body calculation, one usu-ally resorts to approximations, which may not hold for allsituations. Therefore, when only the free-space evolvedtwo-body interaction is used as input in a many-body cal-culation, the results may depend on the cutoff and thisdependence gives not only an estimate of the importanceof the missing 3 N force but also indicates the importanceof the missing many-body terms that may become rele-vant [26].For the induced interaction, except at extremely lowdensities, it is necessary to go beyond the exchange ofsimple particle-hole excitations. Following Ref. [24], wesum the particle-hole bubble series (random-phase ap-proximation, RPA) within the Landau approximationand keep only the lowest order ( L = 0) Landau pa-rameters in the particle-hole interaction. In Ref. [24],as in preceding studies [21, 23], the Landau parame-ters were computed microscopically, including the in-duced interaction in a self-consistent manner (so-calledBabu-Brown theory [27]). However, the resulting Lan-dau parameters, in particular F , were much smaller thanwhat one obtains from phenomenological energy-densityfunctionals such as the Skyrme SLy4 or the Gogny D1Nparameterizations, which have both been fitted to theneutron-matter equation of state. Therefore, we follow amore pragmatic but probably more reliable strategy here,namely to determine the Fermi-liquid parameters (Lan-dau parameters and effective mass) directly from thesephenomenological interactions.The medium corrected interaction is then used in theBCS gap equation and the transition temperature is cal-culated. We note that our results show screening at lowdensities and anti-screening at high densities. This is dif-ferent from the results of Cao et al. in Ref. [24], wherethey predict screening for all densities. Our results forscreening, e.g., Fig. 12, are compatible with QuantumMonte-Carlo (QMC) results [28–30] which rule out theextremely strong screening predicted in earlier calcula-tions [20]. Unfortunately, QMC results are not availablein the density range where we find anti-screening.Apart from the induced interaction, there are othereffects that may modify the BCS results for the transi-tion temperature. If the Fermi momentum k F lies ap-proximately between 1 / | a | ∼ .
05 fm − and 1 /r e ∼ . − , where a is the neutron-neutron ( nn ) scatteringlength and r e the effective range, neutron matter is in a strong-coupling situation, in which pair correlations ap-pear already in the normal phase and modify the criticaltemperature T c [31]. This effect is crucial for the un-derstanding of the BCS-BEC crossover as it exists in ul-tracold Fermi gases or in symmetric nuclear matter [32],where one can pass from Cooper pairs to a Bose-Einsteincondensate (BEC) of dimers (deuterons). For a recentreview article, see [33]. The large value of | a | indicatesthat the nn interaction is almost able to produce a boundstate, and in low-density neutron matter the nn Cooper-pair wave function indeed looks almost like a bound-statewave function [34–36]. In fact, one can reach a situationsimilar to the unitary limit, which is the case of a contactinteraction with | a | → ∞ (i.e., 1 / | a | (cid:28) k F (cid:28) /r e ). Therelevance of BEC-BCS cross-over physics for the descrip-tion of dilute neutron matter was pointed out in manyworks, e.g. [29, 30, 34–37].Note that, although the strong-coupling situation isonly reached at densities below ∼ .
01 times nuclear sat-uration density, it is phenomenologically relevant. Neu-tron matter with such low densities is present betweenthe clusters in the inner crust of neutron stars at averagebaryon densities just above the neutron-drip density of ∼ . × − fm − [38, 39]. Since in this region the diluteneutron gas fills almost the entire volume, it represents asizable contribution to the average baryon density evenif its density is a few thousand times smaller than thedensity inside the clusters.In the unitary limit, the Nozi`eres-Schmitt-Rink (NSR)theory of pair correlations in the normal phase [31] pre-dicts a reduction of the transition temperature T c fromthe BCS result ∼ . E F ( E F being the Fermi energy) to ∼ . E F [40]. These numbers do not include screeningeffects, but as shown recently [41], the inclusion of screen-ing on top of the NSR effect leads to good agreement withresults from experiments with ultracold atoms. In a pre-vious work [37], we had studied neutron matter in theframework of the NSR theory. In the present paper, wewill extend that work to see how the NSR correction ischanged by the induced interaction.This paper is organised as follows: in Sec. II, we re-visitthe BCS gap equation and set up the induced interaction.The effect of the induced interaction on the transitiontemperature is discussed in Sec. III. At low densities,one expects a reduction in T c by a factor of (4 e ) − / ,which is the Gor’kov-Melik-Barkhudarov (GMB) result[42], and this region is studied in detail in Sec. IV. Fi-nally, we turn our attention to the correlations withinthe NSR approach in Sec. V. A summary of our resultsis presented in Sec. VI. Some of the details of the calcu-lations have been moved to appendices to facilitate easeof reading. Numerical results for the matrix elementsof the screened pairing interaction are provided in thesupplemental material [43]. II. FORMALISMA. Gap equation and induced interaction
In BCS theory, the S pairing gap ∆ in neutron mat-ter is determined by the gap equation∆( k ) = − π (cid:90) ∞ dq q V ( k, q ) ∆( q ) tanh (cid:0) E ( q )2 T (cid:1) E ( q ) . (1)Here, V ( k, q ) = (cid:104) k | V S | q (cid:105) denotes the nn interaction inthe S partial wave for in- and outgoing momenta q and k , E q = (cid:112) ( (cid:15) ( q ) − µ ) + ∆( q ) is the quasiparticle energywith (cid:15) ( q ) = q / m ∗ , m ∗ is the neutron effective mass, µ isthe effective chemical potential including the mean-fieldenergy shift, and T is the temperature. Except in somerange of low densities, neutron matter is in the weak-coupling limit, in the sense that ∆( k F ) (cid:28) µ , implying µ ≈ k F / m ∗ , with the Fermi momentum k F = (3 π ρ ) / determined by the neutron number density ρ . Equation(1) with V low k as nn interaction has been solved, e.g., inRef. [44].The critical temperature T c is the highest temperaturefor which Eq. (1) has a non-trivial solution. At T = T c ,one can neglect ∆( q ) in E ( q ), so that Eq. (1) reduces toa linear eigenvalue equation φ ( k ) = − π (cid:90) ∞ dq q V ( k, q ) tanh (cid:0) ξ ( q )2 T c (cid:1) ξ ( q ) φ ( q ) , (2)with ξ ( q ) = (cid:15) ( q ) − µ . We will also write this as | φ (cid:105) = K| φ (cid:105) .Hence, in order to find T c , we diagonalize the integraloperator with the kernel K ( k, q ) = − V ( k, q ) tanh (cid:0) ξ ( q )2 T (cid:1) ξ ( q ) , (3)and T c is the temperature where the largest eigenvalue isequal to unity. In weak coupling, T c is directly related tothe gap at T = 0 by T c = 0 .
57 ∆ T =0 ( k F ).It is widely accepted that an important correction toBCS theory consists in adding to the bare interaction inEq. (1) the contribution of the induced interaction V ind due to the exchange of density and spin-density fluctua-tions. In particular, in the weakly interacting limit, thisleads to the famous Gor’kov-Melik-Barkhudarov (GMB)correction, which reduces the gap and the critical tem-perature by approximately a factor of two compared tothe BCS result [42]. In terms of Feynman diagrams, thiscorrection can be represented as in Fig. 1 (a). Note thatthe dotted interaction lines are meant to represent the an-tisymmetrized interaction. This is very important sincethe dominant S interaction acts only between neutronsof opposite spin and therefore cannot contribute to theshown diagram. However, if the outgoing lines are ex-changed in both the interaction vertices, one obtains adiagram to which it contributes.In nuclear physics, except at extremely low density (seeSec. IV), one is never in a weakly interacting regime. − q , − σ q , σ − q , − σ q , σ p − k σ p σ ˜ V j s l l ˜ V j s l l (a) − q , − σ q , σ − q , − σ q , σ p − k σ p σ p − k σ p σ ˜ V j s l l ˜ V j s l l RPA(b)
FIG. 1: Feynman diagrams representing the induced interac-tion. The wiggly line in diagram (b) is meant to include theRPA bubble summation.
Therefore, the simple particle-hole bubble exchanged inFig. 1 (a) is modified by the residual particle-hole interac-tion as shown in Fig. 1 (b). The wiggly line representingthe particle-hole interaction is meant to include the RPAbubble summation to all orders.Throughout this article, “diagram (a)” and “diagram(b)” refer to the diagrams shown in Fig. 1 (a) and (b).When calculating the diagrams, we make the usual ap-proximation to neglect the energy transfer (static approx-imation) which can be justified by the observation thatthe most important contribution to pairing comes fromscattering of particles near the Fermi surface, so that allin- and outgoing particles have energies close to the Fermienergy (cid:15) F = k F / m ∗ . B. Diagram (a): single-bubble exchange
Let us first discuss the vertices coupling the particlesto the particle-hole excitation, represented as dotted linesin Fig. 1. We assume a general (possibly non-local) in-teraction which is expanded in partial waves. Using thenotation of the left part of Fig. 2, the partial-wave ex-pansion of the interaction reads (cid:104) k , σ ; k , σ | V | k (cid:48) , σ (cid:48) ; k (cid:48) , σ (cid:48) (cid:105) = (cid:88) s,m s ,m (cid:48) s (cid:88) l,l (cid:48) ,m l (cid:88) j C sm s σ σ C sm (cid:48) s σ (cid:48) σ (cid:48) C jm j lm l sm s C jm j l (cid:48) m (cid:48) l sm (cid:48) s × (4 π ) i l (cid:48) − l Y ∗ lm l (Ω Q ) Y l (cid:48) m (cid:48) l (Ω Q (cid:48) ) (cid:104) Q | V sll (cid:48) j | Q (cid:48) (cid:105) , (4) k , σ k , σ k , σ k , σ σ σ σ σ p − k p FIG. 2: Elements of Feynman diagrams to clarify the nota-tion. Left: particle-particle interaction. Right: particle-holepropagator. with Q = k − k , Q (cid:48) = k (cid:48) − k (cid:48) , (5) m (cid:48) l = m l + m s − m (cid:48) s , m j = m l + m s . (6)For the Clebsch-Gordan coefficients, we follow the nota-tion of the book by Varshalovich [45].Then it is straight-forward to obtain for the diagram(a) the following expression [the factor ( −
1) comes fromthe closed Fermion loop]: V a ( q, q (cid:48) ) = ( −
1) 14 π (cid:88) σσ (cid:48) C σ − σ C σ (cid:48) − σ (cid:48) (cid:90) d Ω q π (cid:90) d Ω q (cid:48) π (cid:90) d p (2 π ) n ( p − k ) − n ( p ) (cid:15) ( p ) − (cid:15) ( p − k ) × (cid:88) s m s (cid:88) s m s (cid:88) l l (cid:48) m l (cid:88) l ,l (cid:48) ,m l (cid:88) j j C s m s σ − σ C s m (cid:48) s σ (cid:48) − σ (cid:48) C s m s σ σ (cid:48) C s m (cid:48) s σ (cid:48) σ C j m j l m l s m s C j m j l (cid:48) m (cid:48) l s m (cid:48) s C j m j l m l s m s C j m j l (cid:48) m (cid:48) l s m (cid:48) s × (4 π ) i l (cid:48) − l + l (cid:48) − l Y ∗ l m l (Ω Q ) Y l (cid:48) m (cid:48) l (Ω Q (cid:48) ) Y ∗ l m l (Ω Q ) Y l (cid:48) m (cid:48) l (Ω Q (cid:48) ) (cid:104) Q | ˜ V s l l (cid:48) j | Q (cid:48) (cid:105)(cid:104) Q | ˜ V s l l (cid:48) j | Q (cid:48) (cid:105) , (7)with the following abbreviations: k = q − q (cid:48) , (8) Q = q + p , Q (cid:48) = q (cid:48) − k + p , Q = q + k − p , Q (cid:48) = q (cid:48) − p , (9) σ = m s + σ , σ (cid:48) = m s − σ , (10) m (cid:48) s = m s − σ − σ (cid:48) , m (cid:48) s = m s + σ + σ (cid:48) , (11) m (cid:48) l = m l + m s − m (cid:48) s , m (cid:48) l = m l + m s − m (cid:48) s , (12) m j = m l + m s , m j = m l + m s . (13)The tilde in ˜ V indicates that the matrix element isantisymmetrized, i.e., multiplied by a factor of two inthe surviving channels. For the occupation numbers n ( p )and n ( k − p ) entering the integral in Eq. (7), we use thestep function n ( p ) = θ ( k F − p ), which is a very goodapproximation as long as we are in the weak-couplinglimit (∆ , T (cid:28) µ ). Notice that thenlim k → n ( p − k ) − n ( p ) (cid:15) ( p ) − (cid:15) ( p − k ) = m ∗ δ ( p − k F ) , (14)which is useful when evaluating Eq. (7) for q = q (cid:48) , espe-cially in the case q = q (cid:48) = 0. C. Separation of S = 0 and S = 1 contributions It is instructive to split Eq. (7) into contributionsfrom particle-hole excitations having total spin S = 0 (density waves) and S = 1 (spin-density waves). Inorder to do this, consider the particle-hole propaga-tor shown in the right part of Fig. 2, which is givenby G ( p ) G ( p − k ) δ σ σ δ σ (cid:48) σ (cid:48) with G the uncorrelatedsingle-particle Green’s function. This expression appearsalso in diagram (a) if we formally introduce a summationover σ and σ (cid:48) . The spin part can be decomposed usingthe completeness relation of the Pauli matrices σ δ σ σ δ σ (cid:48) σ (cid:48) = ( δ σ σ (cid:48) δ σ σ (cid:48) + σ σ (cid:48) σ · σ σ σ (cid:48) ) , (15)where the two terms correspond, respectively, to S = 0and S = 1. Likewise, this decomposition can also bewritten in terms of Clebsch-Gordan coefficients as δ σ σ δ σ (cid:48) σ (cid:48) = (cid:88) S,m S ( − − σ (cid:48) C Sm S σ − σ (cid:48) ( − − σ (cid:48) C Sm S σ − σ (cid:48) . (16)In the calculation of V a ( q, q (cid:48) ), it is clear that in the S = 1case each of the three spin projections m S of the particle-hole excitation must give the same contribution. We cantherefore compute the S = 1 contribution by restrictingourselves to the m S = 0 term, or, equivalently, by keepingonly the Pauli matrix σ z in the second term of Eq. (15),and multiplying the result by three. This amounts to thereplacement δ σ σ δ σ (cid:48) σ (cid:48) → δ σ σ (cid:48) δ σ σ (cid:48) (cid:2) − − σ − σ (cid:3) . (17)In this way, we arrive at an alternative expression fordiagram (a): V a ( q, q (cid:48) ) = ( −
1) 14 π (cid:88) σσ (cid:48) C σ − σ C σ (cid:48) − σ (cid:48) (cid:90) d Ω q π (cid:90) d Ω q (cid:48) π (cid:90) d p (2 π ) n ( p − k ) − n ( p ) (cid:15) ( p ) − (cid:15) ( p − k ) × (cid:88) s m s (cid:88) s m s (cid:88) l l (cid:48) m l (cid:88) l ,l (cid:48) ,m l (cid:88) j j C s m s σ − σ C s m (cid:48) s σ − σ (cid:48) C s m s σ σ C s m (cid:48) s σ (cid:48) σ C j m j l m l s m s C j m j l (cid:48) m (cid:48) l s m (cid:48) s C j m j l m l s m s C j m j l (cid:48) m (cid:48) l s m (cid:48) s × (4 π ) i l (cid:48) − l + l (cid:48) − l Y ∗ l m l (Ω Q ) Y l (cid:48) m (cid:48) l (Ω Q (cid:48) ) Y ∗ l m l (Ω Q ) Y l (cid:48) m (cid:48) l (Ω Q (cid:48) ) (cid:104) Q | ˜ V s l l (cid:48) j | Q (cid:48) (cid:105)(cid:104) Q | ˜ V s l l (cid:48) j | Q (cid:48) (cid:105)× (cid:2) − − m s − m s (cid:3) , (18)with the same abbreviations k , Q i , Q (cid:48) i , m (cid:48) li , and m ji asbefore [Eqs. (8), (9), (12), and (13)] but: σ = m s + σ , σ = m s − σ , (19) m (cid:48) s = m s + σ − σ (cid:48) , m (cid:48) s = m s − σ + σ (cid:48) , (20) D. Diagram (b): RPA bubble summation
Let us now turn to diagram (b), which includes theRPA bubble summation. In the present work, we willrestrict ourselves to the Landau approximation and keeponly the lowest-order ( L = 0) Landau parameters. Thenthe particle-hole interaction takes the form f + g σ · σ ,which allows one to sum the RPA bubble series separatelyin the S = 0 and S = 1 channels. The resulting particle-hole interactions are then f RPA ( k ) = f − f Π ( k ) , g RPA ( k ) = g − g Π ( k ) , (21)where f and g are the Landau parameters for S = 0and S = 1, respectively, andΠ ( k ) = − (cid:90) d p (2 π ) n ( p − k ) − n ( p ) (cid:15) ( p ) − (cid:15) ( p − k ) (22) is the static ( ω →
0) limit of the usual Lindhard functionΠ ( k, ω ).It is convenient to introduce the dimensionless Landauparameters F = N f and G = N g , where N = m ∗ k F /π is the density of states at the Fermi surface(including the neutron-matter degeneracy factor of two),and the dimensionless Lindhard function ˜Π = Π /N .Then, Eq. (21) can be rewritten as f RPA ( k ) = F /N − F ˜Π ( k ) , g RPA ( k ) = G /N − G ˜Π ( k ) . (23)At zero temperature, the Lindhard function can be givenin closed form [25],˜Π ( k ) = 12 (cid:34) − − ˜ k / k ln (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) − ˜ k/
21 + ˜ k/ (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:35) , (24)with ˜ k = k/k F .When computing diagram (b), we use again the trickexplained in the derivation of Eq. (18) and compute the S = 1 contribution as three times the m S = 0 term, forwhich σ (cid:48) = σ and σ (cid:48) = σ : V b ( q, q (cid:48) ) = 14 π (cid:88) σσ (cid:48) C σ − σ C σ (cid:48) − σ (cid:48) (cid:90) d Ω q π (cid:90) d Ω q (cid:48) π (cid:90) d p (2 π ) n ( p − k ) − n ( p ) (cid:15) ( p ) − (cid:15) ( p − k ) (cid:90) d p (2 π ) n ( p − k ) − n ( p ) (cid:15) ( p ) − (cid:15) ( p − k ) × (cid:88) s m s (cid:88) s m s (cid:88) l l (cid:48) m l (cid:88) l ,l (cid:48) ,m l (cid:88) j j C s m s σ − σ C s m (cid:48) s σ − σ (cid:48) C s m s σ σ C s m (cid:48) s σ (cid:48) σ C j m j l m l s m s C j m j l (cid:48) m (cid:48) l s m (cid:48) s C j m j l m l s m s C j m j l (cid:48) m (cid:48) l s m (cid:48) s × (4 π ) i l (cid:48) − l + l (cid:48) − l Y ∗ l m l (Ω Q ) Y l (cid:48) m (cid:48) l (Ω Q (cid:48) ) Y ∗ l m l (Ω Q ) Y l (cid:48) m (cid:48) l (Ω Q (cid:48) ) (cid:104) Q | ˜ V s l l (cid:48) j | Q (cid:48) (cid:105)(cid:104) Q | ˜ V s l l (cid:48) j | Q (cid:48) (cid:105)× (cid:2) f RPA ( k ) + 3 ( − − m s − m s g RPA ( k ) (cid:3) , (25)with the same abbreviations k , σ i , m (cid:48) si , m (cid:48) li , and m ji as before [Eqs. (8), (19), (20), (12), and (13)] but: Q = q + p , Q (cid:48) = q (cid:48) − k + p , Q = q + k − p , Q (cid:48) = q (cid:48) − p , (26) -0.8-0.6-0.4-0.2 0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1 1.2 1.4 F , G , m * / m k F (fm -1 )m * /mG F SLy4D1ND1
FIG. 3: Fermi-liquid parameters m ∗ /m , G , and F used inthe present work, obtained from different phenomenologicaleffective interactions: Skyrme parameterization SLy4 (solidlines), and Gogny D1N (short dashes) and D1 (long dashes)parameterizations. III. ANTI-SCREENING DUE TO THE RPAA. Parameters
For the nn interaction in the particle-particle chan-nel, we use the low-momentum interaction V low k fromRef. [46], obtained from the AV interaction by a renor-malization group evolution (using a smooth Fermi-Diracregulator with (cid:15) FD = 0 .
5) to a final cutoff of Λ = 2 fm − .For the purpose of comparing with Ref. [23], we alsoperform calculations with the Gogny force, using the D1parameterization [47] and the more recent D1N parame-terization [48]. For a comparison of the matrix elementsof the Gogny force with those of V low k , and the corre-sponding pairing gaps without screening, see Ref. [49].The explicit expressions for the partial-wave expansionof the Gogny force are given in Appendix A.Concerning the Fermi-liquid parameters, we do not at-tempt to compute them from the microscopic theory, butwe take more phenomenological results from the SLy4parameterization of the Skyrme functional [50] or fromthe D1N parameterization of the Gogny force [48]. Theexplicit formulas are given in Appendix B, and the re-sulting Fermi-liquid parameters m ∗ /m , F , and G areshown in Fig. 3. Since both the SLy4 and the D1N effec-tive interactions have been fitted to the neutron-matterequation of state, it is not surprising that they give al-most identical results for the Landau parameter F at lowdensities. But also the G values are quite close to eachother. Above k F ∼ − , however, the Landau param-eters of SLy4 are clearly smaller (in absolute value) thanthose of D1N. Note also that SLy4 systematically yieldsa smaller effective mass m ∗ than D1N. For a comparison -80-60-40-20 0 20 40 60 80 0 0.2 0.4 0.6 0.8 1 1.2 1.4 V i nd ( q , q ) ( M e V f m ) q (fm -1 ) k F = 0.8 fm -1 (a)(b)(a)+(b) j max = 012345 FIG. 4: Convergence of the induced interaction with respectto variation of the maximum angular momentum j max used inthe partial wave expansion of the bare interaction. The figureshows the diagonal matrix elements for k F = 0 . − . Upperthin curves: results for diagram (a), lower thin curves: resultsfor diagram (b), thick curves: sum (a)+(b). The bare interac-tion in this example is V low k , and the Fermi-liquid parametersare those of SLy4. with Ref. [23], we also used the D1 parameterization ofthe Gogny force [47], the resulting Fermi-liquid parame-ters are also shown in Fig. 3. B. Induced interaction
In order to calculate the induced interaction in prac-tice, we restrict the partial-wave expansion in Eq. (4) tosome maximum angular momentum, j ≤ j max . The mul-tidimensional integrals in Eqs. (7) or (18), and (25) arecomputed using Monte-Carlo integration. Data files con-taining tables of the pairing interaction with and withoutthe induced interaction are provided in the supplementalmaterial [43].First, we have to check that convergence w.r.t. j max has been reached. This is indeed the case for j max =3, as can be seen in Fig. 4. As one can see from thisfigure, for the example k F = 0 . − , the net effect ofthe sum of diagrams (a) and (b) is attractive, i.e., thestrong repulsion generated by diagram (a) is more thancompensated for by the attractive diagram (b).This result is in contrast to previous studies [23, 24]where it was found that the contribution of diagram (b)is attractive but not strong enough to compensate forthe repulsion generated by diagram (a). Let us there-fore analyse our result in more detail. It is known thatthe exchange of S = 0 excitations (density fluctuations)is attractive and that of S = 1 excitations (spin-densityfluctuations) is repulsive [24, 51]. This is also the case inour calculation, as shown in Fig. 5, again for the example -80-60-40-20 0 20 40 60 80 0 0.2 0.4 0.6 0.8 1 1.2 1.4 V i nd ( q , q ) ( M e V f m ) q (fm -1 ) k F = 0.8 fm -1 (a)(a)+(b)(a)+(b)(a)(a)+(b) S=0S=1total FIG. 5: Induced pairing interaction due to the exchange of S = 0 (dashed lines) and S = 1 (dotted lines) excitations.The thin lines represent the contributions of diagrams (a)only, and the thick lines are the sums of diagrams (a) and(b). The thick solid line is the sum of S = 0 and S = 1contributions. The parameters are the same as in Fig. 4. k F = 0 . − . If there was only the single bubble ex-change [diagram (a)], the repulsive contribution of S = 1excitations would be three to four times larger than theattractive one of S = 0 excitations. However, the inclu-sion of the RPA [diagram (b)] acts differently in the cases S = 0 and S = 1 because the Landau parameters haveopposite signs. In the S = 0 case, since F <
0, the ef-fect of diagram (a) is enhanced, while in the S = 1 case,since G >
0, the effect of diagram (a) is reduced. There-fore, with the inclusion of the RPA, the attraction dueto the exchange of density waves can finally win againstthe repulsive effect of the spin-density waves.
C. Critical temperature
We can now use the induced interaction V ind = V a + V b and replace the bare interaction V in the gap equation(1) by V + V ind . The resulting critical temperature T c as a function of the Fermi momentum k F is shown inFig. 6. A sample of the results is also listed in Table I.The corresponding pairing gaps ∆( k F ) at T = 0 can beobtained, to a very good approximation, by multiplying T c with 1.76. The dashed line represents the result ob-tained with the bare interaction V low k . The maximumcritical temperature is reached at k F ≈ . − . Whenone includes the induced interaction due to diagram (a)alone, pairing is very strongly suppressed, as shown bythe dotted line. Finally, when including diagrams (a) and(b), one finds that the critical temperature is lowered atlow density, but increased at high density. The changefrom screening to anti-screening is at k F ≈ .
73 fm − ,consistent with our results discussed in Sec. III B where T c ( M e V ) k F (fm -1 )V low k +(a)+(b)V low k +(a)V low k FIG. 6: Critical temperature T c as a function of the Fermimomentum k F , obtained with the V low k interaction (Fermi-liquid parameters from the Skyrme force SLy4). Dashes: re-sult obtained using only the bare interaction; dots: result ob-tained including diagram (a); solid line: full result includingalso diagram (b).TABLE I: Critical temperature as a function of the Fermimomentum k F , obtained with V low k interactions and Fermi-liquid parameters from SLy4. T (bare) c is obtained with thebare interaction, while T (screened) c includes the effect of V ind = V a + V b . The columns marked Λ = 2 fm − correspond to theparameters given in Sec. III A, while for the columns markedΛ = 2 . k F , a V low k interaction with a density dependentcutoff and a different regulator was used, see Sec. IV B.Λ = 2 fm − Λ = 2 . k F k F T (bare) c T (screened) c T (bare) c T (screened) c (fm − ) (MeV) (MeV) (MeV) (MeV)0 .
08 0 . . . . . .
212 0 .
167 0 .
206 0 . . .
752 0 .
523 0 .
743 0 . . .
27 1 .
02 1 .
27 1 . . .
48 1 .
68 1 .
48 1 . . .
18 1 .
94 1 .
18 1 . . .
485 0 .
964 0 .
474 0 . . .
184 0 . . . . we found that at 0 . − the attractive effect of S = 0excitations is stronger than the repulsive effect of S = 1excitations. Whether the net effect of the induced in-teraction is attractive (i.e., anti-screening) or repulsive(i.e. screening), depends of course on the density andon the values of the Landau parameters. With decreas-ing density, the RPA bubbles of diagram (b) become lessimportant and therefore the repulsive effect of diagram(a) wins. This explains why, at very low density, thefull result and the result obtained with only diagram (a)become equal, as one can also see in Fig. 6.To check how sensitive our results are to the details ofthe model, we repeated the calculation with the D1N andD1 Gogny forces. In these cases, the same interaction isused for the bare pairing force, for the vertices enteringthe induced interaction diagrams (a) and (b), and for the T c ( M e V ) k F (fm -1 )D1N+(a)+(b)D1ND1+(a)+(b)D1 FIG. 7: Same as Fig. 6, but here the Gogny D1N and D1interactions are used in the particle-particle channel and forthe Fermi-liquid parameters.
Fermi-liquid parameters. The results are shown in Fig. 7.Of course, since the S matrix elements of the differ-ent interactions are not the same, there is already somedifference at the level of the bare interaction [49]: themaximum is slightly shifted and the gap survives up tohigher density. However, the effect of the induced inter-action is qualitatively the same as in Fig. 6, i.e., the gapis reduced at low density and increased at high density.The change from screening to anti-screening happens atabout the same density as with V low k (with Fermi-liquidparameters from SLy4) in Fig. 6, and compared to the V low k results the anti-screening effect at high density iseven stronger with both the D1N and the D1 Gogny in-teractions. IV. THE LOW-DENSITY LIMIT
As one sees from Fig. 6, with the V low k interactionwith a fixed cutoff of 2 fm − , screening gets weak at lowdensity and finally at k F (cid:46) . − one recovers theBCS result. However, at k F (cid:28) / | a | , the GMB resultshould be valid, predicting a reduction of T c by a factor of(4 e ) − / ≈ .
45. Therefore, let us study the low-densitylimit in more detail.
A. Failure of the weak coupling formula
As we have seen, the contribution of diagram (b) be-comes negligible at low density. Concerning diagram(a), it seems natural to concentrate on matrix elements V a ( q, q (cid:48) ) with q, q (cid:48) (cid:39) k F . If k F becomes small, this meansthat also q and q (cid:48) and hence all the momenta Q etc.that appear in Eq. (7) become small. Therefore, we canreplace (cid:104) Q i | ˜ V s i l i l (cid:48) i j i | Q (cid:48) i (cid:105) q,q (cid:48) ,k F → −−−−−−−→ V (0 , δ s i δ l i δ l (cid:48) i δ j i (27) (the factor of two accounts for the antisymmetrization of˜ V ), and Eq. (7) simplifies tremendously to V a ( q, q (cid:48) ) ≈ − πN | V (0 , | (cid:104) ˜Π (cid:105) . (28)In this expression, we have used the angle-averaged Lind-hard function (cid:104) ˜Π (cid:105) = 12 (cid:90) − d cos θ ˜Π (cid:0)(cid:112) q + q (cid:48) − qq (cid:48) cos θ (cid:1) , (29)see appendix C. In particular, we get V a ( k F , k F ) ≈ πN | V (0 , | ln 4 e . (30)Following well-known weak-coupling arguments [25],the gap and critical temperature should be proportionalto e / [2 πN V ( k F ,k F )] . If we replace V by V + V a in theapproximation given in Eq. (30), we find that the gap andthe critical temperature should indeed be reduced by thefactor (4 e ) − / , in contradiction to our numerical resultswhich show that at low density T c is not modified at allby screening. Obviously the weak-coupling formula doesnot apply in the present case, although we are clearly ina weak coupling situation since T c (cid:28) (cid:15) F . Note that thereare a couple of cases in nuclear physics where the weakcoupling formula is known to fail [52].When using the weak coupling formula, one assumesthat the kernel K ( k, q ) given in Eq. (3) is sharply peakedat q = k F and that this peak gives the dominant contri-bution to the integral in the gap equation. However, wewill show that the contribution of the peak is not domi-nant at low density, and this is the reason why the weakcoupling formula fails in this case.Remember that the critical temperature is given bythe temperature where the largest eigenvalue η of thekernel K ( k, q ) given in Eq. (3) is equal to unity. Thecorresponding eigenvector | φ (cid:105) can be found by numericaldiagonalization, its representation in momentum space, φ ( q ) = (cid:104) q | φ (cid:105) , is a smooth function of q which has ap-proximately the shape of V ( q, k F ). If we normalize theeigenvector to (cid:104) φ | φ (cid:105) = (2 /π ) (cid:82) dq q | φ ( q ) | = 1, we canwrite the eigenvalue η as η = (cid:104) φ |K| φ (cid:105) = 4 π (cid:90) dq q (cid:90) dk k φ ( k ) K ( k, q ) φ ( q ) . (31)To measure the importance of the peak of the kernel at q = k F , we can look at this integral as a function of itsupper limit q max , I η ( q max ) = − π (cid:90) q max dq q φ ( q ) tanh (cid:0) ξ ( q )2 T (cid:1) ξ ( q ) × (cid:90) ∞ dk k φ ( k ) V ( k, q ) . (32)At T = T c , we know that I η → q max → ∞ since η = 1. For the weak coupling formula to be valid,the main contribution to the integral should come from I η ( q m a x ) q max (fm -1 )k F (fm -1 )0.0121.45 0 0.1 0 0.01 0.02 FIG. 8: Measure of the contribution of different momenta tothe gap equation as defined in Eq. (32), for two different den-sities ( k F = 0 .
012 fm − (solid line) and 1 .
45 fm − (dashes),indicated by the thin vertical lines). The integrals were cal-culated with the V low k interaction m ∗ from the Skyrme forceSLy4) at the respective critical temperatures. q ≈ k F , i.e., I η should be close to the step function θ ( q max − k F ). In Fig. 8 we show the behavior of I η fortwo cases, k F = 1 .
45 fm − (dashes) and k F = 0 .
012 fm − (solid line). In both cases, we are in the weak-couplinglimit, in the sense that T c /E F is very small (of the or-der of 10 − ). In the case k F = 1 .
45 fm − , we see thatabout 80% of the integral come from momenta close to k F , so that in this case T c is indeed determined to alarge extent by V ( k F , k F ). But in the low-density case, k F = 0 .
012 fm − , the situation is completely different.Although there is again a sharp rise of I η at q ≈ k F (vis-ible in the zoom), its contribution to the total integral isless than 10%. The largest contribution to the integralcomes from momenta that are considerably larger than k F .Let us now look at the matrix elements V ( q, k F ) for k F = 0 .
012 fm − with and without screening, which aredisplayed in Fig. 9. The screening correction is limitedto the tiny region q (cid:46) .
05 fm − ∼ k F , because ofthe strong momentum dependence of the angle-averagedLindhard function. But as we have seen before, this smallregion contributes only about 10% to the integral in thegap equation, and therefore the screening correction haspractically no effect on the gap or T c .The observation that the screening effect disappearsat low density is not a singular feature of our calcula-tion, but it can also be found in the existing literature[24]. However, as we will discuss below, there are otherproblems with the low-density limit. Taking these intoaccount, we will eventually retrieve the GMB result. B. Failure of perturbation theory anddensity-dependent cutoff
When calculating diagrams (a) and (b), we use the bareinteraction V perturbatively to describe the vertex cou- -80-60-40-20 0 0 0.5 1 1.5 2 V ( q , k F ) ( M e V f m ) q (fm -1 ) k F = 0.012 fm -1 screenedapprox.V low k -80-79 0 0.05 FIG. 9: Matrix elements V ( q, k F ) of the bare ( V low k ,dashes) and of the screened (solid line) interaction for k F =0 .
012 fm − . For comparison, we display also the screened in-teraction obtained with the analytical approximation Eq. (28)(dots, almost indistinguishable from the solid line). Thescreening correction is so tiny that it is almost invisible onthe big graph, see the inset for a zoom. The thin vertical lineindicates q = k F .FIG. 10: Higher-order ladder diagrams in the 3p1h verticeswhich are not included in the present work. pling the particles to the particle-hole excitations. Sincewe are using renormalized interactions whose matrix ele-ments decrease rapidly with increasing relative momenta Q i and Q (cid:48) i , which are typically of the order of k F , thismay be a good approximation at higher densities. How-ever, for small Q i and Q (cid:48) i , as they appear at low densities,we know from the large value of the nn scattering length a that the perturbative treatment must fail [53, 54].When looking at the historical work by GMB [42], oneobserves that they compute the correction in a differentway. Namely, instead of using the potential V in thedashed interaction vertices of diagram (a), they use a/m .This amounts to including, at least approximately, theresummation of ladder diagrams as shown in Fig. 10.In contrast to the Gogny interaction, therenormalization-group evolved V low k interaction givesus the additional freedom to change the cutoff Λ. Onthe one hand, by lowering the cutoff, the interactiongets obviously “more perturbative”. In this sense, it istempting to lower the cutoff as much as possible. Infact, for q, q (cid:48) < Λ and Λ →
0, the matrix elements getmore and more attractive and flow towards the constant0 -600-400-200 0 0 0.01 0.02 0.03 0.04 0.05 V ( q , k F ) ( M e V f m ) q (fm -1 ) k F = 0.012 fm -1 Λ = 0.03 fm -1 screenedV low k FIG. 11: Same as Fig. 9 but now calculated with an inter-action V low k evolved to a much lower cutoff Λ = 0 .
03 fm − =2 . k F (and with an exponential instead of Fermi-Dirac regu-lator, see text). a/m as V ( q, q (cid:48) ) ≈ (cid:16) ma − m Λ π (cid:17) − . (33)This means that the contribution of higher-order lad-der diagrams gets progressively included, via the renor-malization group flow, in the two-body matrix elements,while the loop integrals become suppressed, and as a re-sult, it should be possible to work with a Born approx-imation to the T matrix at low cutoffs. On the otherhand, one of course must not lower the cutoff below therelevant momentum scale of the order of k F .The cutoff dependence of the gap (without screen-ing corrections) was investigated in Ref. [55]. Numer-ically, we obtain cutoff independent results for T c atthe BCS level in the whole range of densities for Λ (cid:38) . k F , if we use an exponential regulator of the formexp( − ( k / Λ ) n exp ) with n exp = 5. (With the Fermi-Dirac regulator and with (cid:15) FD = 0 . − that we usedbefore we would need somewhat larger cutoffs.)So, let us see what we find when we choose instead ofa constant cutoff Λ = 2 fm − the lowest possible cutofffor each value of k F , i.e., Λ = 2 . k F .As an example, let us consider as in Fig. 9 the case k F = 0 .
012 fm − . If we evolve the cutoff to the lowestpossible value for this k F , i.e., to Λ = 2 . k F = 0 .
03 fm − ,we obtain the matrix elements V ( q, k F ) shown in Fig. 11.As in Fig. 9, the dashed line represents V low k withoutscreening and the solid line has screening included. Themost obvious difference between Figs. 9 and 11 is that,when the cutoff is lowered, the V low k matrix elements(dashed lines) get more attractive, cf. Eq. (33). How-ever, the renormalization group flow does not only en-sure that the low-energy scattering in free space remainsunchanged, but also the gap and T c at the BCS level(i.e., without screening) remain the same, as mentionedabove. But the results with screening change. Now, themodification of the interaction due to screening (differ-ence between the solid and the dashed lines in Fig. 11) T c ( s c r ee n e d ) / T c ( b a r e ) k F (fm -1 ) Λ = 2 fm -1 Λ = 2.5 k F FIG. 12: Reduction of the critical temperature due to thescreening correction V ind = V a + V b as a function of k F , ob-tained with the constant cutoff Λ = 2 fm − (dashes) and withthe density-dependent cutoff Λ = 2 . k F (solid line), respec-tively. extends over the whole momentum range up to ∼ Λ, andtherefore the screening will reduce T c , contrary to whathappened in the case Λ = 2 fm − .Since the results for T c obtained without the screeningcorrection is the same as the one we obtained before forΛ = 2 fm − , we can concentrate on the correction of T c due to screening. In Fig. 12, we therefore display theratio of T c with screening to T c without screening as afunction of k F . The red dashes correspond to the resultsshown already in Fig. 6, obtained with a constant cutoffΛ = 2 fm − , and we clearly see that the effect of thescreening correction vanishes at low density, as explainedin Sec. IV A. The new results obtained with the variablecutoff 2 . k F are shown as the blue solid line. We seethat now the reduction of T c due to screening survives atlow densities, and in the limit k F → e ) − / ≈ .
45 predicted by GMB.Note that the original GMB paper [42] considers k F | a | (cid:28)
1, i.e., in the case of neutron matter, k F (cid:28) .
05 fm − . V. EFFECT OF THENOZI`ERES-SCHMITT-RINK CORRECTIONA. Brief summary of the formalism
In our previous work [37], we had studied neutron mat-ter within the NSR approach using only the free-spacerenormalized effective interaction V . In the presentwork, we will revisit the inclusion of preformed pairsabove T c , including the induced interaction V ind shownin Fig. 1. For the sake of completeness, we summarizebriefly the key ideas and formulas of the NSR approach.For more details, we refer the reader to Ref. [37].Within the NSR approach, for a given chemical poten-tial µ , the density of the interacting neutrons is enhancedby the pair correlations that build up as a precursor ef-fect to the superfluid phase transition already above T c .1Therefore, the total density of neutrons, ρ tot , can be writ-ten as ρ tot = ρ + ρ corr . (34)The uncorrelated neutron density ρ is given by ρ = 2 (cid:90) d k (2 π ) f ( ξ ( k )) , (35)where f ( ξ ) = 1 / ( e βξ + 1) is the Fermi-Dirac distributionfunction (with β = 1 /T ) and the factor of 2 arises dueto the spin degeneracy. The correlated density, ρ corr , inthe imaginary-time formalism [25], is calculated to firstorder in the single-particle self-energy Σ as ρ corr = 2 (cid:90) d k (2 π ) β (cid:88) ω n (cid:0) G ( k , ω n ) (cid:1) × [Σ( k , iω n ) − Re Σ( k , ξ ( k ))] , (36)where ω n are the fermionic Matsubara frequencies and G = 1 / ( iω n − ξ ( k )) is the uncorrelated single-particleGreen’s function. The subtraction of the on-shell self-energy in the square bracket of Eq. (36) is absent in theoriginal NSR approach. It takes into account the factthat G includes already the in-medium quasiparticle en-ergy ξ ( k ) which therefore must not be shifted by theself-energy [32, 56].Let us consider the first term without the subtraction.Σ( k , iω n ) is calculated within the ladder approximation,i.e.,Σ( k , iω n ) = (cid:90) d K (2 π ) β (cid:88) ω N G ( K − k , ω N − ω n ) × (cid:10) K − k (cid:12)(cid:12) T ( K , iω N ) (cid:12)(cid:12) K − k (cid:11) , (37)where T ( K , iω N ) is the in-medium T matrix for thebosonic Matsubara frequency ω N and total momentum K . The T -matrix is subsequently expanded in a partialwave basis and we pick out only the s -wave contribution.Following the steps outlined in [37] and analytically con-tinuing to real ω , one obtains for the correlated densitywithin the NSR approach: ρ corr , = − ∂∂µ (cid:90) K dK π (cid:90) dωπ g ( ω ) Im Tr log (cid:0) − V G (2)0 (cid:1) . (38)Here, g ( ω ) = 1 / ( e βω −
1) is the Bose function, thetrace is taken w.r.t. the relative momentum q , G (2)0 = Q ( K, q ) / ( ω − K / m ∗ − q /m ∗ +2 µ ) is the angle-averaged(since we consider only the s wave) retarded two-particleGreen’s function, with Q ( K, q ) the Pauli-blocking factor1 − f ( ξ ( K / − q )) − f ( ξ ( K / q )) averaged over the an-gle between K and q . Working in the basis where V G (2)0 is diagonal, one can write Eq. (38) as ρ corr , = − ∂∂µ (cid:90) K dK π (cid:90) dωπ g ( ω ) × (cid:88) ν Im log(1 − η ν ( K, ω )) , (39) V i nd ( q , q ) ( M e V f m ) k F = 0.2 fm -1 Λ = 0.5 fm -1 K = 0K = k F K = 2 k F K = 4 k F -1 )k F = 0.2 fm -1 Λ = 2 fm -1 F = 0.8 fm -1 Λ = 2 fm -1 FIG. 13: Dependence of V ind on the momentum K of thecenter of mass. For a density k F = 0 . − (left panels),the K -dependence is extremely weak even for momenta K exceeding 2 k F . For k F = 0 . − (upper right panel), the K dependence is somewhat stronger but still too weak to makea significant contribution. where η ν are the (complex) eigenvalues of V G (2)0 .However, as mentioned below Eq. (36), one needs tocorrect for the shift of the quasiparticle energies thatcomes from the real part of the single-particle self-energy.Following [37], we approximate Σ( k , ξ ( k )) by the first-order (Hartree-Fock) self-energy and finally arrive at thefollowing correction: ρ corr , = ∂∂µ (cid:90) K dK π π (cid:90) q dq g (cid:0) K m ∗ + q m ∗ − µ (cid:1) × V ( q, q ) Q ( K, q ) , (40)which is added to Eq. (39).In Ref. [37], the interaction V that was used inEqs. (39) and (40) was the V low k interaction obtainedfrom AV via the free-space renormalization group evo-lution. But it seems straight-forward to include in addi-tion the medium corrections from diagrams (a) and (b),i.e., to use V = V + V ind . The only complication is thatso far we calculated V ind only for a pair at rest, while weshould now take into account the finite center of massmomentum K of the pair.To obtain the screening correction V ind for finite K ,some minor modifications of Eqs. (7) and (25) are neces-sary. Details are given in Appendix D. We have checkedthat, at least for T = T c , the contributions to the inte-grals in Eqs. (39) and (40) come only from K (cid:46) k F . Asseen in Fig. 13, numerically it turns out that the K de-pendence of V ind is very weak for K < k F in the rangeof k F where the NSR correction can be expected to beimportant. We will therefore neglect this K dependenceand use in Eqs. (39) and (40), the screening correctioncalculated for K = 0.There are a couple more points that need to be dis-cussed. For instance, now one has two different den-sities, the uncorrelated one, ρ , and the corrected one,2 F (fm -1 )00.10.20.30.40.5 ρ c o rr ( - f m - ) Λ = 2.0 fm -1 Λ = 2.5 k F0 -4 -3 -2 -1 FIG. 14: The un-subtracted correlated density ρ corr , as afunction of the Fermi momentum k F with and without thescreening correction, calculated at the respective critical tem-peratures T c . Here, the black solid lines and the red dashedlines show the results for the two cutoffs of Λ = 2 . − andΛ = 2 . k F . The thin lines contain only V , while the thicklines include the induced interactions. The inset in the figuremagnifies the cutoff dependence in ρ corr , at low densities.The Fermi-liquid parameters are calculated using the SLy4interaction. ρ tot . The question arises which density one should use inthe calculation of the induced interaction V ind . Since V ind is computed with uncorrelated propagators and occupa-tion numbers, it seems more appropriate to take only theuncorrelated density ρ into account in the calculationof V ind . From the derivation of Eqs. (39) and (40) itis also clear that the derivatives ∂/∂µ should be takenwith the interaction V ind kept constant (and the effectivemass m ∗ , too). This points to fundamental problems ofthe present approach, which is clearly not a fully consis-tent treatment of both particle-particle and particle-holefluctuations. Nevertheless, we expect to get at least arough idea about the change of the NSR effect when thepair correlations are modified by screening. B. Results
Before discussing the critical temperature as a func-tion of density, let us look at the density correction. Theun-subtracted correlated density, ρ corr , as a function ofthe Fermi-momentum corresponding to the uncorrelateddensity ρ , denoted here as k F = (3 π ρ ) / , is shown inFig. 14. The black solid lines and the red dashed linesrepresent two different cutoff choices, a constant cutoffΛ = 2 . − and a density dependent cutoff Λ = 2 . k F .The thin lines show the correlated density ρ corr , withonly the free-space interaction V . Analogous to Fig. 5 ofRef. [37], we see that ρ corr , with only V is independent ofthe cutoff. With the inclusion of the induced interaction F0 (fm -1 )00.20.40.60.811.21.4 ρ c o rr ( - f m - ) Λ = 2.0 fm -1 Λ = 2.5 k F0 Thin lines: V Thick lines: V + V ind FIG. 15: Subtracted correlated density ρ corr as a function of k F with and without screening, calculated at the respectivecritical temperatures. See Fig. 14 for details. (thick lines) we note that the cutoff dependence of ρ corr , is again negligible, except at very low densities (see in-set), where we found stronger screening with the variablecutoff compared to the fixed cutoff (cf. Fig. 12). In addi-tion, up to k F ∼ . − , the correlated density ρ corr , with the induced interaction is smaller than the corre-lated density without the induced interaction, consistentwith the earlier observation that the induced interactionscreens V . However, in the range of Fermi-momentawhere the induced interaction anti-screens V , the cor-related density ρ corr , is larger than the correspondingquantity without the induced interaction.Let us now turn our attention to the correlated den-sity with the first-order (Hartree-Fock) subtraction, ρ corr .The dependence of ρ corr on k F is shown in Fig. 15. Asin Fig. 14, the black solid lines and the red dashed linesshow results for the two different cutoffs: the constantcutoff Λ = 2 . − and the density dependent cutoffΛ = 2 . k F , respectively. For low k F , we see that thecorrelated density with the inclusion of the induced in-teraction (thick lines) is smaller than in the V -only case(thin lines) which is consistent with the screening of V by V ind and similar to the trend seen in Fig. 14. However,what is surprising is that even in the region where V ind anti-screens V , the correlated density gets smaller withthe inclusion of V ind compared to the V -only case. Fur-ther, one notices strong cutoff dependence in the low k F region if one compares the solid black line with the reddashed line, both with and without the inclusion of theinduced interaction. Both these observations are com-pletely different from Fig. 14 and are clearly the effect ofthe Hartree-Fock subtraction. For the density dependentcutoff, at low-densities, this subtraction should get betteras the interaction gets more perturbative at smaller cut-offs. However, at high densities, where the subtraction ρ corr , is almost of the same magnitude as ρ corr , itself,3 F (fm -1 )00.511.52 T c ( M e V ) F (fm -1 ) V + V ind V + V ind (NSR)BCS Λ = 2.0 fm -1 Λ = 2.5 k F0 FIG. 16: T c versus k F : (Left panel) Results with fixed cutoffΛ = 2 fm − ; (Right panel) density dependent cutoff 2 . k F .The green dashed-dotted lines are the full results includingthe induced interaction V ind and the correlated density ρ corr in the NSR framework. For comparison, we also show theBCS result (only V and ρ , black solid lines) and the resultsobtained with the induced interaction V ind but without theNSR correction (red dashed lines). the Hartree-Fock approximation is not precise enough togive a reliable result for the subtracted ρ corr . Hence,the suppression of the correlated density for higher k F in Fig. 15, once the induced interaction is included, isprobably unphysical. Fortunately, in this region, ρ corr iscompletely negligible compared to ρ .Now we are in the position to discuss the final resultsfor the critical temperature T c as a function of k F , dis-played in Fig. 16. Note that in the NSR framework, T c asa function of µ is computed as usual, and only the relationbetween µ and k F (and ρ ) is changed. Here, k F denotesthe Fermi momentum corresponding to the total densityincluding ρ corr , i.e., k F = (3 π ρ tot ) / (green dashed-dotted lines). As a consequence, the presence of the cor-related density ρ corr shifts the curve slightly to the right.In order to make easy comparisons, we also show the BCSresult (solid line) and the results obtained with V ind butwithout the NSR correction (red dashed lines). In bothpanels, we note that the pair correlations lower the tran-sition temperature compared to the one with screeningalone at the same k F . However, the trends already ob-served with the medium corrections (Figs. 6 and 12), i.e.,screening at low densities and anti-screening at high den-sities, remain unchanged, since the NSR effect is muchweaker than the screening or anti-screening effect of V ind .Please notice that the relation ∆ T =0 ( k F ) = 1 . T c fora given k F , mentioned in Sec. III C, is not valid for theNSR results. VI. CONCLUSIONS
It has been known for a long time that screening cor-rections have a very strong effect on the superfluid tran-sition temperature of neutron matter. Also the fact thatthe RPA, diagram (b), reduces the effect of diagram (a),has been known before [24]. However, in Ref. [24] the ef-fect of diagram (b) was too weak to overcome the strongscreening generated by diagram (a), while we find that,around n (cid:38) . − .
02 fm − , the net effect of V ind is at-tractive and screening turns into anti-screening. A simi-lar effect was found in Ref. [21], but only at much higherdensities ( n (cid:38) .
07 fm − ). There are three main differ-ences between our calculation and that of Ref. [24]. First,we are using V low k while in [24] the Br¨uckner G matrixwas used in the vertices. Second, while we keep the fullmomentum dependence of the non-local interaction, thevertices in [24] were replaced by an average matrix ele-ment. Probably the most important difference, however,is the choice of the Landau parameters. Here, we takethem from a phenomenological energy density functional(SLy4). Since this functional was fitted to QMC resultsfor the neutron matter equation of state, we assume thatthe Landau parameters are rather well determined. Theanti-screening effect arises primarily from the enhance-ment of the attractive density ( S = 0) fluctuations dueto the strongly negative f parameter. In [24], on thecontrary, the Landau parameters were obtained follow-ing the so-called Babu-Brown theory as explained in [21].This results in particular in a much smaller (less nega-tive) value of the f parameter, and as a consequence, thedensity fluctuations are not strong enough to compensatefor the repulsive effect of the spin-density ( S = 1) fluc-tuations.We addressed in some detail the problem of the lowdensity limit. When a constant (density-independent)potential V is used in the vertices of diagram (a), thescreening effect disappears at low density, although fromthe weak-coupling formula one would conclude that thegap should be reduced by the factor (4 e ) − / predicted byGMB [42]. We explained why the weak-coupling formulafails in this particular case. We then observed that GMBused the full T matrix instead of the potential V in thevertices of diagram (a). This allowed us to finally recoverthe GMB result, namely by using for each density a V low k interaction evolved to a cutoff Λ that scales with k F . Inthis way, one ensures that, on the one hand, one doesnot cut the relevant degrees of freedom ( q (cid:46) k F ), andon the other hand, the Born term is already a reasonableapproximation to the full T matrix at momenta of theorder of q ∼ k F .In the last part of the paper we discussed the effectof preformed pairs on the critical temperature T c in theNSR framework. In spite of some cutoff and regulator de-pendence in the detailed study of the correlated density ρ corr , one can clearly see that due to ρ corr the criticaltemperature T c for a given density is slightly reduced.But this effect is much less important than the induced4interaction. Compared to ultracold atoms in the uni-tary limit or even on the BEC ( a >
0) side of the BCS-BEC crossover, neutron matter remains more or less ina weakly coupled regime at all densities.There remain obviously many open questions. For in-stance, as discussed in [24], the reduction of the quasi-particle residue
Z < T c , andthis effect has not been included in the present study.Another point that clearly needs to be improved is theLandau approximation in the RPA. In principle, it is onlyvalid for momentum transfer k (cid:28) k F , but in the inducedinteraction, the relevant range of momentum transfers is0 ≤ k ≤ k F . In the framework of Skyrme interactions itis actually straight-forward to solve the RPA beyond theLandau approximation, and this issue will be addressedin a future study.Concerning the meaning of the density dependent cut-off introduced in Sec. IV B, one might wonder how thisis related to the so-called functional renormalization-group approach in which one solves flow equations inthe medium, integrating out all momenta except theFermi surface. Such approaches have been used to in-clude screening corrections in a non-perturbative way forneutron matter [57] and ultracold atoms [58, 59]. In thecontext of the small cutoff, one should also mention thatlowering the cutoff induces three- and higher-body inter-actions. These are neglected in V low k since it is obtainedfor two particles in free space. A better approach in thisrespect would be the in-medium similarity renormaliza-tion group [60], which allows one to include many-bodyeffects at least approximately into the effective two-bodyinteraction.Because of the extreme sensitivity of the gap and thecritical temperature to the details of the effective inter-action, it seems likely that large theoretical uncertaintieswill remain. Maybe astrophysical observations of neutronstars can help to decide which theory is correct. Appendix A: Partial wave expansion of the Gognyforce
We expand the Gogny force as given in Ref. [47] intopartial waves, neglecting the spin-orbit term as in [23].The resulting matrix elements in the nn channel read: (cid:104) Q | V ls | Q (cid:48) (cid:105) = 14 π (cid:88) i =1 , [ W i − H i + ( − S ( M i − B i )] × ( √ πµ i ) e − ( Q + Q (cid:48) ) µ i / i l ( QQ (cid:48) µ i /
2) (A1)where i l ( z ) = (cid:112) π/ zI l +1 / ( z ) is a modified sphericalBessel function of the first kind [61]: i ( z ) = sinh( z ) /z ,etc. The antisymmetrized matrix elements are then ob-tained by (cid:104) Q | ˜ V ls | Q (cid:48) (cid:105) = [1 + ( − l + s ] (cid:104) Q | V ls | Q (cid:48) (cid:105) . The den-sity dependent contact term of the Gogny force does notcontribute since it acts only in the neutron-proton chan-nel. Concerning the values of µ i , W i , H i , B i , and M i , weuse either the parameterization D1 [47] to compare withRef. [23] or the more recent parameterization D1N [48]. Appendix B: Fermi-liquid parameters
In this work, we use the Fermi-liquid parameters fromthe SLy4 parameterization of the Skyrme functional [50]or from the D1N parametrization of the Gogny force [48].The explicit expressions in terms of the Skyrme-force pa-rameters t i , x i ( i = 0 . . . σ read [62]1 m ∗ = 1 m + [ t (1 − x ) + 3 t (1 + x )] ρ , (B1) f = t (1 − x ) + [ t (1 − x ) + 3 t (1 + x )] k F + t (1 − x )(1 + σ )(2 + σ ) ρ σ , (B2) g = t ( x −
1) + [ t ( x −
1) + t (1 + x )] k F + t ( x − ρ σ . (B3)In the case of the Gogny force, one obtains the followingexpressions for the Fermi-liquid parameters [23]:1 m ∗ = 1 m + m √ πk F (cid:88) i =1 , µ i ( W i + 2 B i − H i − M i ) × z i e − z i i ( z i ) , (B4) f = (cid:88) i =1 , ( √ πµ i ) W i + B i − H i − M i ) − ( W i + 2 B i − H i − M i ) e − z i i ( z i )] , (B5) g = (cid:88) i =1 , ( √ πµ i ) B i − M i ) − ( W i − H i ) e − z i i ( z i )] , (B6)where z i = k F µ i / Appendix C: Angle-averaged Lindhard function
For q, q (cid:48) (cid:54) = 0, the general explicit expression for theangle-averaged Lindhard function defined in Eq. (29)reads (cid:104) ˜Π (cid:105) = −
13 + k F qq (cid:48) [ F (2 − x − ) + F (2 + x − ) − F (2 − x + ) − F (2 + x + )] , (C1)with F ( x ) = x (6 − x ) ln | x | and x ± = | q ± q (cid:48) | /k F . Inthe special case of interest q = q (cid:48) = k F mentioned inthe main text this gives (cid:104) ˜Π (cid:105) = − ln 4 e ≈ − . q (cid:54) = q (cid:48) = 0 or q (cid:48) (cid:54) = q = 0 reads (cid:104) ˜Π (cid:105) = x − x artanh (cid:16) x (cid:17) − , (C2)with x = q/k F or q (cid:48) /k F , respectively. In the special case q = q (cid:48) = 0, one obtains (cid:104) ˜Π (cid:105) = −
1. For q (cid:29) k F or q (cid:48) (cid:29) k F , (cid:104) ˜Π (cid:105) tends to zero.5 Appendix D: Computation of the screeningcorrections for pairs with finite total momentum
In Fig. 1 and the corresponding Eqs. (7) and (25), wehave considered from the beginning a pair at rest (withrespect to the medium). However, for the NSR correc-tion, one needs pairs with finite total momentum K . Inorder to compute the screening corrections V ind = V a + V b for K (cid:54) = 0, one has to change the definitions of the vec-tors Q , Q (cid:48) , Q , and Q (cid:48) that appear in Eqs. (7) and(25). For diagram (a), one has to replace Eq. (9) by Q = q + p − K , Q (cid:48) = q (cid:48) − k + p − K , Q = q + k − p K , Q (cid:48) = q (cid:48) − p K . (D1)For diagram (b), the definition (26) has to be replacedby Q = q + p − K , Q (cid:48) = q (cid:48) − k + p − K , Q = q + k − p K , Q (cid:48) = q (cid:48) − p K . (D2) However, for diagram (b), this is not sufficient, becausewe used the isotropy to replace the sum over the threespin projections m S = − , , S = 1 particle-holeexcitation by the contribution of m S = 0, multiplied bythree. But for K (cid:54) = 0, the isotropy is lost and thereforethe contributions of the three spin projections will not beequal any more. Nevertheless, after summation over m S ,the final result for V b can only depend on K = | K | andnot on the direction of K . Hence, we can average over theangle of K . By doing so, we have restored the isotropyand it is therefore again sufficient to compute only thecontribution of m S = 0 and to multiply the result bythree. [1] N. Chamel and P. Haensel, Living Rev. Relativity, 11, 10(2008).[2] P. W. Anderson and N. Itoh, Nature , 25 (1975).[3] D. Pines and M. A. Alpar, Nature , 27 (1985).[4] B. Haskell and A. Sedrakian, arXiv:1709.10340 [astro-ph.HE].[5] D. G. Yakovlev and C. J. Pethick, Ann. Rev. Astron.Astrophys. , 169 (2004).[6] M. Fortin, F. Grill, J. Margueron, D. Page, and N. San-dulescu, Phys. Rev. C , 065804 (2010).[7] D. Page, J. M. Lattimer, M. Prakash and A. W. Steiner,Astrophys. J. , 1131 (2009).[8] A. Gezerlis, C. J. Pethick and A. Schwenk, in: K. H.Bennemann and J. B. Ketterson (eds.), Novel Superflu-ids, Volume 2 (Oxford University Press, 2014)[9] A. B. Migdal, Zh. Eksp. Teor. Fiz.
249 (1959) [Sov.Phys. JETP , 176 (1960)]; Nucl. Phys. , 655 (1959).[10] V. L. Ginzburg and D. A. Kirzhnits, Zh. Eksp. Teor. Fiz. , 2006 (1964) [Sov. Phys. JETP , 1346 (1965)].[11] V. L. Ginzburg, Usp. Fiz. Nauk , 601 (1969) [Sov.Phys.-Uspekhi , 241 (1969)]; J. Stat. Phys. , 3 (1969).[12] G. Baym, C. J. Pethick, D. Pines, and Malvin Ruderman,Nature , 872 (1969).[13] A. Sedrakian and J. W. Clark, arXiv:1802.00017 [nucl-th].[14] D. J. Dean and M. Hjorth-Jensen, Rev. Mod. Phys. ,607 (2003).[15] S. Srinivas and S. Ramanan, Phys. Rev. C , 064303(2016).[16] C. Drischler, T. Kr¨uger, K. Hebeler and A. Schwenk,Phys. Rev. C , 024302 (2017).[17] P. Papakonstantinou and J. W. Clark, J. Low. Temp.Phys. , 361 (2017).[18] A. Rios, A. Polls and W. H. Dickhoff, J. Low. Temp.Phys. , 234 (2017). [19] A. Rios, D. Ding, H. Dussan, W. H. Dickhoff, S. J. Witteand A. Polls, J. Phys. Conf. Ser. , 012014 (2018).[20] J. Wambach, T. L. Ainsworth and D. Pines, Nucl. Phys.A , 128 (1993).[21] H. J. Schulze, J. Cugnon, A. Lejeune, M. Baldo andU. Lombardo, Phys. Lett. B , 1 (1996).[22] C. Shen, U. Lombardo and P. Schuck, Phys. Rev. C ,061302 (2003).[23] Caiwan Shen, U. Lombardo, and P. Schuck, Phys. Rev.C , 054301 (2005).[24] L. G. Cao, U. Lombardo, and P. Schuck, Phys. Rev. C , 064301 (2006).[25] A. L. Fetter and J. D. Walecka, Quantum Theory ofMany-Particle Systems (McGraw-Hill, New York, 1971).[26] K. Hebeler and A. Schwenk, Phys. Rev. C , 014314(2010).[27] S. Babu and G. E. Brown, Ann. Phys. (NY) , 1 (1973).[28] S. Gandolfi, A. Yu. Illarionov, S. Fantoni, F. Pederiva,and K. E. Schmidt, Phys. Rev. Lett. 101, 132501, (2008).[29] T. Abe and R. Seki Phys. Rev. C , 054003 (2009).[30] A. Gezerlis and J. Carlson, Phys. Rev. C , 025803(2010).[31] P. Nozi`eres and S. Schmitt-Rink, J. Low. Temp. Phys. , 195 (1985).[32] M. Jin, M. Urban, and P. Schuck, Phys. Rev. C ,024911 (2010).[33] G. C. Strinati, P. Pieri, G. R¨opke, P. Schuck, and M.Urban, Phys. Rep. , 1 (2018).[34] M. Matsuo, Phys. Rev. C , 044309 (2006).[35] J. Margueron, H. Sagawa, and K. Hagino, Phys. Rev. C , 064316 (2007).[36] B. Y. Sun, H. Toki, and J. Meng, Phys. Lett. B , 134(2010).[37] S. Ramanan and M. Urban, Phys. Rev. C , 054315(2013) [38] J. W. Negele and D. Vautherin, Nucl. Phys. A , 298(1973).[39] M. Baldo, E. E. Saperstein, and S. V. Tolokonnikov,Phys. Rev. C , 025803 (2007).[40] C. A. R. S´a de Melo, M. Randeria, and J. R. Engelbrecht,Phys. Rev. Lett. , 3202 (1993).[41] L. Pisani, A. Perali, P. Pieri, and G. C. Strinati, Phys.Rev. B , 014528 (2018).[42] L. P. Gor’kov and T.K. Melik-Barkhudarov, Zh. Eksp.Teor. Fiz. , 1452 (1961) [Sov. Phys. JETP , 1018(1961)].[43] See ancillary files of the present arxiv submission fordata files containing matrix elements with and withoutinduced interaction.[44] J. Kuckei, F. Montani, H. M¨uther, and A. Sedrakian,Nucl. Phys. A , 32 (2003).[45] D. A. Varshalovich, A. N. Moskalev, V. K. Khersonskii, Quantum Theory of Angular Momentum (World Scien-tific, Singapore, 1988).[46] S. K. Bogner, R. J. Furnstahl, S. Ramanan, and A.Schwenk, Nucl. Phys. A , 79 (2007).[47] J. Decharg´e and D. Gogny, Phys. Rev. C , 1568 (1980).[48] F. Chappert, M. Girod, and S. Hilaire, Phys. Lett. B , 420 (2008).[49] A. Sedrakian, T. T. S. Kuo, H. M¨uther, and P. Schuck,Phys. Lett. B , 68 (2003).[50] E. Chabanat, P. Bonche, P. Haensel, J. Meyer, and R.Schaeffer, Nucl. Phys. A , 710 (1997). [51] H. Heiselberg, C. J. Pethick, H. Smith, and L. Viverit,Phys. Rev. Lett. , 2418 (2000).[52] J. W. Clark, in R. A. Broglia and V. Zelevinsky (eds.), Fifty Years of Nuclear BCS (World Scientific, Singapore2013), 360.[53] S. K. Bogner, A. Schwenk, R. J. Furnstahl and A. Nogga,Nucl. Phys. A , 59 (2005).[54] S. K. Bogner, R. J. Furnstahl, S. Ramanan andA. Schwenk, Nucl. Phys. A , 203 (2006).[55] K. Hebeler, A. Schwenk and B. Friman, Phys. Lett. B , 176 (2007).[56] R. Zimmermann and H. Stolz, Phys. Status Solidi B ,151 (1985).[57] A. Schwenk, B. Friman and G. E. Brown, Nucl. Phys. A , 191 (2003).[58] K. B. Gubbels and H. T. C. Stoof, Phys. Rev. Lett. ,140407 (2008).[59] S. Floerchinger, M. Scherer, S. Diehl, and C. Wetterich,Phys. Rev. B , 174528 (2008).[60] H. Hergert, S. K. Bogner, J. G. Lietz, T. D. Morris, S. J.Novario, N. M. Parzuchowski, and F. Yuan, Lect. NotesPhys. , 477 (2017).[61] M. Abramowitz and I. A. Stegun (eds.), Handbook ofMathematical Functions (Dover, New York 1965).[62] J. Margueron, J. Navarro, and Nguyen Van Giai, Phys.Rev. C66