# Temperature dependent BCS equations with continuum coupling

## I Introduction

The effect of temperature upon pairing correlations in nuclei was studied long ago , after the theory of superfluidity was introduced to describe the appearance of a gap in the low-energy spectrum of finite nuclei . The extension of the theory of superfluidity to the finite temperature case was prompted by the search of evidence about the nuclear phase transition from a superfluid to a normal state, in analogy with the one found in condensed matter systems . The vanishing of pairing correlations with temperature is related to the fact that the excitation energy breaks pairs of particles which block the single particle levels close to the Fermi surface, where the pairing correlations are initiated.

The change of the pairing properties with the temperature was also studied in competition with the angular momentum. In this case a new effect was observed, which is due to the fact that at zero temperature the single-particle states close to the Fermi level are already partially blocked by the un-paired nucleons which form the finite angular momentum of the nucleus. Thus, when the temperature is switched on the first effect is a depletion of the partially blocked single-particle states which, in turn, induces an enhancement of the pairing correlations. This behaviour is different for the case of thermal excitations of nuclei with zero angular momentum, where the pairing correlations decrease monotonically with the increase of the temperature.

The temperature effects on pairing correlations were studied both in the BCS  and HFB  approximations, usually with a constant single-particle level density approximation. The available theoretical evidence shows that for medium and heavy-mass nuclei the pairing correlations disappear for a critical temperature of the order of 0.5-1.0 MeV. The excitation energies corresponding to this critical temperature are quite small for nuclei close to the beta stability line and therefore in all such calculations the coupling with the continuum spectrum was neglected. The situation is different for nuclei which are far from the stability line. In this case the Fermi level lies close to the continuum threshold and the coupling with the continuum becomes important. We shall focus here on the study of this coupling since, to our knowledge, the effect of the continuum coupling on thermodynamical properties of superfluid nuclei has not been investigated so far.

The main problem in dealing with the continuum coupling in statistical calculations of excited nuclei is the fact that the particles moving in the continuum have a finite probability to be emitted from the nucleus. In other words, such processes are time dependent and they are difficult to accomodate in stationary models as BCS or HFB. But to extent stationary many body theories to time-dependent formalisms is not an easy undertaking. In fact it may even be a not well defined task since the initial conditions may induce chaotic solutions. These features were already recognized in the beginning of quantum mechanics. Thus, Gamow and Wigner tried to reconciliate the outgoing character of the decaying process with the conveniences of stationarity by solving the Schrödinger equation with outgoing boundary conditions (for references see ). The corresponding solutions are related to the complex poles of the S-matrix, which define the so-called Gamow resonances. If the resonances are narrow the real parts of the complex poles of the S-matrix give the positions of the resonances while the corresponding imaginary parts give the decay widths. The narrow resonances are very important to describe nuclear correlations, especially in unbound or excited nuclei, because in these states the nucleons could move within the nuclear volume during a certain minimum time, so that they can interact with each other. One thus expects that a basis formed by bound and narrow Gamow resonances would provide a convenient framework to describe decaying processes . However, the drawback of the calculations based on this representation is that one gets complex probabilities which are not always easy to interpret .

Another alternative to describe unstable nuclei is to use a basis consisting of scattering states, instead of Gamow states, in the vecinity of the resonant poles . In this case all quantities are real and one does not have problems of interpreting complex probabilities. Within this representation the widths of the resonances are obtained by evaluating the derivative of the corresponding phase shifts. This also defines the continuum level density  commonly used to estimate the contribution of the continuum to nuclear partition functions . The escape of the particles which move in resonant states is thus treated in these calculations as a stationary process, reflected by a constant particle density at large distances from the nucleus . The underling picture  is an excited nucleus in dynamical equilibrium with an external nucleonic gas, whose contribution should be eventually extracted. Recently a similar framework was used to include continuum coupling in BCS equations at zero temperature . In this paper we will extend this method to study temperature dependent BCS equations (TBCS).

## II Formalism

The standard procedure to derive the temperature dependent BCS equations is to minimize the grand potential corresponding to a pairing Hamiltonian. For a bound single particle spectrum with energies $`ϵ_j`$ and a constant pairing interaction of strength G, the gap and particle number equations at a finite temperature $`T`$ are given by ,

$$\frac{2}{G}=\underset{j}{}\frac{12f_j}{2E_j},$$
(1)   
$$N=\underset{j}{}\frac{1}{2}[1\frac{ϵ_j\lambda }{E_j}(12f_j)]$$
(2)   
where $`E_j=[(ϵ_j\lambda )^2+\mathrm{\Delta }^2]^{1/2}`$ is the quasiparticle energy, $`f_j=(1+exp(\beta E_j))^1`$ is the Fermi distribution function and $`\beta =1/kT`$.

In principle the single particle energies $`ϵ_j`$ depend also upon temperature because the average mean field is a function of the nuclear excitations. However, in self-consistent temperature dependent Hartree-Fock calculations one finds  that for temperatures below $`T=1`$ MeV, which is the range explored in temperature dependent BCS calculations, the single-particle spectrum is virtually the same as the one at zero temperature. We will also assume here that the single particle energies are those at zero temperature.

The contribution of the continuum on thermodynamical properties of finite nuclei was studied mainly in connection with the problem of the liquid-gas phase transition  as well as in temperature dependent shell corrections , but without including pairing correlations. In these calculations the effect of the continuum was introduced into the thermodynamical quantities through the level density. Thus the grand potential for a non-interacting system was taken as 

$$\mathrm{\Omega }=T(g_b(ϵ)+\stackrel{~}{g}(ϵ))ln(1+exp[\beta (ϵ\lambda )])𝑑ϵ$$
(3)   
where $`g_b(ϵ)=_j\delta (ϵϵ_j)`$ is the level density of the bound spectrum and $`\stackrel{~}{g}(ϵ)`$ is the level density associated with the positive energy spectrum. In Ref.  it is shown that the grand potential (3) describes a nucleus in dynamical equilibrium with a nucleonic gas. As discussed above, this is due to the fact that in a stationary treatment the nucleons scattered in the continuum are permanently emitted from the nucleus. To obtain the proper grand potential, i. e. the one corresponding to the nucleus itself, one should take away from Eq. (3) the contribution of the nucleonic gas. This can be done by subtracting from the grand potential (3) the grand potential of the free nucleonic gas , or by replacing in Eq. (3) the level density $`\stackrel{~}{g}(ϵ)`$ by the quantity 

$$g(ϵ)=\stackrel{~}{g}g_{free}=\frac{1}{\pi }\underset{j}{}\frac{d\delta _j}{dϵ}$$
(4)   
where $`g_{free}`$ is the level density in the absence of the mean field and $`\delta _j`$ is the phase shift. The quantity $`g(ϵ)`$ is the continuum level density . It takes into account the contribution of the resonant part of the continuum spectrum. The continuum coupling can thus be included by replacing in the grand potential (3) the density $`\stackrel{~}{g}`$ by the continuum level density $`g(ϵ)`$. This is a general recipe which can be applied to all quantities derived from the grand potential, as e. g. the energy and entropy of the system.

The incorporation of the continuum in interacting systems can readly be performed following a similar prescription. Thus the contribution of the continuum to the grand potential of an excited superfluid nucleus can be expressed in term of the continuum level density as in Eq. (3), with the difference that now instead of the single particle energies one should use the quasiparticle energies. The corresponding TBCS equations can be obtained directly from Eqs. (1) and (2) by replacing the level density of the bound states with the total level density, i. e. by $`g_b(ϵ)+g(ϵ)`$. Thus the TBCS equations with continuum coupling become,

$$\frac{2}{G}=\underset{j}{}\frac{12f_j}{2E_j}+g(ϵ)\frac{12f(ϵ)}{2E(ϵ)}𝑑ϵ,$$
(5)   
$$N=\underset{j}{}\frac{1}{2}[1\frac{ϵ_j\lambda }{E_j}(12f_j)]+g(ϵ)\frac{1}{2}[1\frac{ϵ\lambda }{E(ϵ)}(12f(ϵ))]𝑑ϵ$$
(6)   
where the second term gives the contribution of the continuum to the pairing correlations. In the limit T=0 one gets the same equations as in Refs . For temperatures higher than the critical temperature the gap vanishes and the particle number equation is similar to the one used in thermodynamical calculations of non-interacting systems .

The contribution of the continuum to the energy and to the entropy can be introduced in a similar fashion, as mentioned above. One gets,

$$E=\underset{j}{}n_jϵ_j+g(ϵ)n(ϵ)ϵ𝑑ϵ\frac{\mathrm{\Delta }^2}{G},$$
(7)   
$$S=\underset{j}{}\{ln(\stackrel{~}{f}_j)+\beta f_jE_j\}+g(ϵ)\{ln(\stackrel{~}{f}_\nu (ϵ))+\beta f_\nu (ϵ)E_\nu (ϵ)\}𝑑ϵ$$
(8)   
where $`n_j`$ is the occupancy of the state of energy $`ϵ_j`$, given by

$$n_j=\frac{1}{2}[1\frac{ϵ_j\lambda }{E_j}(12f_j)]$$
(9)   
and $`\stackrel{~}{f}_j=1+exp(\beta E_j)`$. A similar notation is used for the corresponding quantities in the continuum.

In the limit of vanishing widths the continuum level density becomes a sum of Dirac delta functions and the resonances act as bound states (”quasibound states”). This is the case for, e. g., protons trapped in a high Coulomb barrier. If a resonance is not narrow then a pair scattered to that resonance have a large probability to escape from the system. Therefore, its contribution to the pairing correlation is small as compared with the corresponding contribution from a quasibound state with similar energy and angular momentum .

Although narrow resonances play a fundamental role in the enhancement of pair correlations (and in most other measurable physical processes as well) their proper inclusion in the applications may be difficult. The reason for this is that in the region of a narrow resonance the level density increases abruptly and the numerical evaluation of the integrals in the equations above may require extremely small mesh intervals. One can circumvent this problem by noticing that the resonance is narrow because there is a pole of the $`S`$ matrix which is very close to the real energy axis. Therefore, one can evaluate the integrals by changing the integration path, that is by choosing a contour $`C`$ in the complex energy plane which embodies the real axis around the narrow resonance, and by thereafter applying the Cauchy theorem . Eq. (7) thus becomes

$$E=\underset{i}{}n_iϵ_i+\underset{\nu }{}n(𝐄_\nu )𝐄_\nu +_Cg(ϵ)n(ϵ)ϵ𝑑ϵ\frac{\mathrm{\Delta }^2}{G}$$
(10)   
where $`n(𝐄_\nu )`$ are the occupation probabilities calculated in the complex poles $`𝐄_\nu `$ enclosed by the path $`C`$. If one neglects the contribution of the integral over the contour $`C`$ in Eq. (10) then the energy $`E`$ would become complex. This is the case in representations based on Gamow resonances only . As already mentioned, within such representations one obtains complex physical quantities which have to be interpreted. For the energy such a task is rather easy, since already by looking at the temporal evolution of the wave function one realizes that the real part corresponds to the actual energy of the system while the imaginary part is related to the corresponding decay probability. This interpretation is valid if the resonance is narrow, i. e. if the ratio between the width and the energy of the resonance is small . But the interpretation of complex probabilities in general is not so straightforward . These problems do not appear here because in Eq. (10) the contribution of the integral over the contour $`C`$ is included and therefore the energy $`E`$ is real.

## III Numerical application

In order to illustrate how the continuum affects the properties of superfluid nuclei close to the drip line, in what follows we present the results given by the TBCS equations for the isotope $`{}_{}{}^{84}Ni`$.

We calculated the single particle spectrum by using the HF approximation and the Skyrme III interaction . The resonant energies are defined as the energies where the phase shift passes through $`\pi /2`$ with a positive slope . The width is extracted from the value of the energy where the derivative of the phase shift is half of its maximum value. The bound states outside the closed shell N=50 and the resonant states considered in the TBCS calculations are listed in Table 1. It can be seen that these states form the equivalent of the major shell N=50-82. One notices that the relative positions of the single-particle states differ substantially with the ones corresponding to beta stable nuclei. Thus, the states with low angular momenta are shifted down as compared with the ones with high angular momenta. This is a general feature related to the diffusivity of the mean field in nuclei close to the drip line , which is larger than the corresponding one in stable nuclei .

The TBCS equations (5) and (6) are solved starting with the HF spectrum calculated at zero temperature. In the absence of experimental information on heavy $`Ni`$ isotopes in the open shell $`N=5082`$ (the heaviest known $`Ni`$ isotope is the double magic nucleus <sup>78</sup>Ni) we use for the strength of the pairing force the standard value $`G=25/A`$ MeV, where $`A`$ is the mass number.

The variation of the gap with temperature is shown in Fig. 1. In the limit of $`T=0`$ the gap has the value $`\mathrm{\Delta }(0)=0.955`$ MeV and decreases to zero at the critical temperature $`T_c=0.524`$ MeV. This critical temperature is smaller by about $`3.8\%`$ than the one obtained from the relation $`T_c=0.57\mathrm{\Delta }(0)`$, which is the value predicted by a constant level density approximation .

In Fig. 1 it is also shown the gap one would get if the width of the resonances are neglected, i.e. if in Eqs. (5-6) one would take instead of the continuum level density a sum of Dirac delta functions. The gap at zero temperature is in this case bigger than in the previous case . As discussed above, the wider is a resonance the smaller is its contribution to the pairing correlation. That is, a pair scattered to a wide resonance spends less time in that state as compared to the time that the same pair would spend in a quasibound state and, therefore, the wider resonance contributes less to the pairing correlations. This has the important consequence that neglecting the widths of the resonances one increases the correlations of the system and therefore increases the critical temperature. This is what happens, e. g., with calculations that quantize the continuum by using an impenetrable box if the dimensions of the box are not extremely large (so that one gets a dense enough spectrum in the energy regions of resonant states). However, one has to mention that this is not a serious problem if the resonances are very narrow.

Neglecting the widths one gets for the critical temperature the value $`T_c=0.722`$ MeV, which is only $`1\%`$ larger than the value one would obtain in a constant level density approximation. This indicates that this approximation would work quite well for nuclei close to the proton drip line, where the width of the proton single-particle resonances may be so narrow that they can be neglected.

In Fig. 2 we show the dependence of the excitation energy upon temperature up to the point where the superfluid phase vanishes. It can be seen that if the effect of the width of resonant states is neglected, then the slope of the excitation energy becomes much smaller. This is also a manifestation of the fact that by neglecting the widths, the ground state becomes more correlated and therefore stiffer to thermal excitations. The same effect is observed for the entropy, as seen in Figure 3.

At critical temperature the excitation energy has the value 2.185 MeV and the entropy is 5.708. For a constant level density $`g`$ the excitation energy at critical temperature is given by  $`E_caT_c^2+0.5g\mathrm{\Delta }(0)^2`$, where $`a=\frac{\pi ^2}{3}g`$ is the level density parameter, which in this case acquires the value $`a=3.958`$ MeV<sup>-1</sup>. We found that neglecting the widths the level density parameter increases to the value $`a=4.143`$ MeV<sup>-1</sup>. This implies that the difference in the excitation energies seen in Fig. 2 is esentially due to the mean occupancy of the resonant states (which is larger if the widths are neglected) and not due to the effective level density parameter.

In conclusion, in this paper we have extended the temperature dependent BCS equations by introducing the coupling with the single-particle continuum. The contribution of the continuum is given by the resonant states and their effect is taken into account through the continuum level density.

We found that the widths of the resonances affect significantly all physical quantities. In particular, the pairing correlations are diminished and this modifies significantly the value of the critical temperature at which the supefluid phase disappears. Also the dependence upon temperature of the excitation energy and of the entropy is considerably affected by the widths of the resonances, i.e. by their lifetime.

However, in the case of proton superfluidity the single-particle resonant states close to the continuum threshold may have a very small width and therefore they can be treated in the TBCS calculations as quasibound states, neglecting their widths altogether.

(O.C.) is partially supported by the grant (PICT0079) and by the CONICET, Argentina. (N.S.) is supported by the Wenner-Gren Foundation.

## Figure Captions

Figure 1: Dependence of the gap parameter $`\mathrm{\Delta }`$ upon the temperature $`T`$. The dashed line corresponds to the case when the resonant states are considered as quasibound states.

Figure 2: Excitation energy plotted as a function of temperature. The dashed line corresponds to the case when the resonant states are considered as quasibound states. The energy is plotted up to the critical temperature.

Figure 3: Entropy plotted as a function of temperature. The dashed line corresponds to the case when the resonant states are considered as quasibound states. The entropy is plotted up to the critical temperature.