Energy density functional on a microscopic basis
In recent years impressive progress has been made in the development of highly accurate energy density functionals, which allow to treat medium-heavy nuclei. In this approach one tries to describe not only the ground state but also the first relevant excited states. In general, higher accuracy requires a larger set of parameters, which must be carefully chosen to avoid redundancy. Following this line of development, it is unavoidable that the connection of the functional with the bare nucleon-nucleon interaction becomes more and more elusive. In principle, the construction of a density functional from a density matrix expansion based on the effective nucleon-nucleon interaction is possible, and indeed the approach has been followed by few authors. However, to what extent a density functional based on such a microscopic approach can reach the accuracy of the fully phenomenological ones remains an open question. A related question is to establish which part of a functional can be actually derived by a microscopic approach and which part, on the contrary, must be left as purely phenomenological. In this paper we discuss the main problems that are encountered when the microscopic approach is followed. To this purpose we will use the method we have recently introduced to illustrate the different aspects of these problems. In particular we will discuss the possible connection of the density functional with the nuclear matter Equation of State and the distinct features of finite size effects proper of nuclei.
The microscopic calculation of ground state energy and particle density of medium and heavy nuclei based on realistic nucleon-nucleon interaction requires the solution of a formidable many-body problem. For this reason effective nucleon-nucleon interactions have been introduced, like the Skyrme  and Gogny  forces. They simplify enormously the problem since, by construction, they must be used at the mean field level, and the calculation of the mean single particle potential and of the ground state energy becomes easily manageable. The number of parameters which enter these effective forces is typically around ten and they are adjusted to reproduce finite nuclei and some equilibrium nuclear matter properties. However, recently also data from a theoretically determined neutron matter Equation of State (EOS)  have been used as input (see also an older attempt in this direction in ). Generally these forces give rise to an effective nucleon mass with typically in the non-relativistic framework. With these ingredients nuclear mean field theories are very successful to describe nuclear properties as, e.g. binding energies, radii. Their use can be extended to the evaluation of more realistic single level schemes, of the nuclear excitation spectra, in particular giant resonances, of the fission barriers, of the nucleon-nucleus optical potential, and so on. The price to be payed is that the connection with the bare nucleon-nucleon interaction is not apparent. The first attempt to connect the effective NN interaction and the underlying bare interaction was the density matrix expansion of Negele and Vautherin , based on the assumption that the effective NN interaction can be identified essentially with the Brueckner G-matrix calculated in nuclear matter and on the gradient expansion of the non-local one-body density matrix. See also another attempts in the same direction in Ref. [6, 7, 8]. This approach was recently generalized to include also the vector (spin) part of the density matrix at the same level of accuracy as the scalar part . The challenging program of deriving the effective interaction to be used in nuclear structure studies from the bare NN interaction has still to be completed .
An approach similar to the Skyrme-like one is the energy density functional (EDF) approach, where the basic quantity to start with is directly the functional that expresses the energy in terms of the matter density and its gradients. The functional must be minimized to obtain the actual ground state energy and matter density profile. Following the method of Kohn and Sham [10, 11, 12, 13, 14, 15, 16], developed in atomic, molecular and solid state physics, the minimization procedure can be performed by introducing a set of auxiliary single particle states and taking for the kinetic energy the Slater form in this basis. In this way the minimization procedure gives Hartree-like equations for the single particle states, where the interaction part includes in an effective way the overall exchange and correlation contributions. The latter, in condensed matter physics, is taken from accurate calculations for the homogeneous electron system or in a purely phenomenological fashion. This is not the usual way of proceeding in nuclear physics, although it exists, to our knowledge, an earlier attempt in this direction, see ref. .
In this paper we will discuss to what extent the EDF approach in nuclear structure can be based on microscopic many-body results on nuclear matter Equation of State (EOS) and which accuracy can be reached. To this purpose we will rely on recent achievements along this line of research and we will analyze the main open problems that hinder the development of the microscopic many-body theory of the nuclear EDF.
2 Basic formalism
Despite the fact that the applicability of the KS-DFT approach to self-bound systems, as nuclei, is not obvious [18, 19, 20, 21, 22], it is commonly believed that the basis of KS-DFT lies in the Hohenberg-Kohn (HK) theorem , which states that for a Fermi system, with a non-degenerate ground state, the total energy can be expressed as a functional of the density only. Such a functional reaches its variational minimum when evaluated with the exact ground state density. Furthermore, in the standard KS-DFT method one introduces an auxiliary set of orthonormal single particle wave functions , where is the number of particles, and the density is assumed to be given by
where and stand for spin and iso-spin indices. The variational procedure to minimize the functional is performed in terms of the orbitals instead of the density. As in condensed matter and atomic physics the HK functional ] is split into two parts: . The first piece corresponds to the uncorrelated part of the kinetic energy and within the KS method it is written as
The other piece contains the potential energy as well as the correlated part of the kinetic energy.
Then, upon variation, one gets a closed set of Hartree-like equations with an effective potential, the functional derivative of with respect to the local density . Since the latter depends on the density, and therefore on the ’s, a self-consistent procedure is necessary. The equations are exact but they only can be of some use if a reliable approximation is found for the otherwise unknown density functional . It has to be stressed that in the KS-DFT formalism the exact ground state wave function is actually not known, the density being the basic quantity. In nuclear physics, contrary to the situation in condensed matter and atomic physics, the contribution of the spin-orbit interaction to the energy functional is very important. Non-local contributions have been included in DFT in several ways already long ago (see  for a recent review of this topic). Consequently, the spin-orbit part also can be split in an uncorrelated part plus a remainder. The form of the uncorrelated spin-orbit part is taken exactly as in the Skyrme  or Gogny forces . We thus write for the functional in the nuclear case , where we explicitly split off the Coulomb energy because it is a quite distinct part in the Hamiltonian. It shall be treated, as usual, at lowest order, i.e. the direct term plus the exchange contribution in the Slater approximation, that is , and with and the proton/neutron density.
Let us now discuss the nuclear energy functional contribution which contains the nuclear potential energy as well as additional correlations. We shall split it in a finite range term to account for correct surface properties and a bulk correlation part that we take from a microscopic infinite nuclear matter calculation . Thus our final KS-DFT -like functional reads:
For the finite range term we make the simplest phenomenological ansatz possible
with proton/neutron and the volume integral of . The substraction in (4) is made in order not to contaminate the bulk part, determined from the microscopic infinite matter calculation. Finite range terms have already been used earlier, generalizing usual Skyrme functionals (see e.g [28, 29, 30]). In this study, for the finite range form factor we make a simple Gaussian ansatz: . We choose a minimum of three open parameters: , and . The only undetermined and most important piece in (3) is then the bulk contribution . As already mentioned, we obtain from microscopic infinite matter calculations, using a realistic bare force, together with a converged hole-line expansion . We first reproduce by interpolating functions the correlation part of the ground state energy per particle of symmetric and pure neutron matters, and then make a quadratic interpolation for asymmetric matter. Finally the total correlation contribution to the energy functional in local density approximation reads:
where and are two interpolating polynomials for symmetric and pure neutron matter, respectively, at the density , and is the asymmetry parameter. The interpolating polynomial for symmetric matter has been constrained to allow a minimum exactly at the energy = - 16.00 MeV and Fermi momentum = 1.36 fm, i.e. = 0.16 fm. This is within the uncertainty of the numerical microscopic calculations of the EOS. It has to be stressed that the use of a polynomial in density is just for practical reasons.
The constrained fit was performed by keeping the EOS as smooth as possible, thus allowing for some very small deviations from the microscopic calculations below saturation density. An interpolating fit which goes exactly through the calculated EOS, as performed in , gives a not good enough saturation point (typically E/A = -15.6 MeV, =1.38 fm).
As discussed in , the low density behavior of the nuclear matter EOS is quite intricate and usually not reproduced by Skyrme and Gogny functionals ( see also ref. ), missing quite a substantial part of binding. We show our EOS for nuclear and neutron matter in Fig. 1. Since we want to construct the EDF, as much as possible, on the basis of the microscopic calculations, the bulk part of the functional, directly related to bare NN and TBF (three-body forces), is determined once and for all and we will use it in (3) together with LDA.
The only open parameters are, therefore, the ones contained in the finite range surface part, Eq.(4), and the strength of the spin-orbit contribution, that we fit to reproduce finite nuclei properties. We, thus, follow exactly the strategy employed in condensed matter. On the contrary, in nuclear physics, almost exclusively a different strategy is usually adopted (see, however, Ref. with some ingredients similar as in the present approach ): functionals like the one of (3), i.e. bulk, surface, etc., were globally parametrized with typically of the order of ten parameters which, then, were determined fitting simultaneously some equilibrium nuclear matter (binding energy per particle, saturation density, incompressibiliy, etc.) and finite nuclei properties. However, in this way, bulk and surface are not properly separated and early attempts used to miss important infinite nuclear matter properties, as, e.g. stability of neutron matter at high density  and other stability criteria. Modern Skyrme forces , like the Saclay-Lyon (SLy) ones, explicitly use the high density part () of microscopic neutron matter calculations for the EOS in the fitting procedure  and thus avoid collapse. Therefore, modern functionals usually reproduce reasonably well microscopically determined EOS for neutron and nuclear matter . Examples are, among others, the SLy-forces  (see Fig.1 in ) and the Fayans functional DF3 (see e.g. Fig. 3 in ). In this presentation we follow the just discussed alternative approach, that is different on a qualitative level from the usual and allows, via the fit, to reproduce very accurately the microscopic infinite matter results in the whole range of densities considered. This may be important for surface properties and neutron skins in exotic nuclei, what shall be investigated in the future.
For open shell nuclei, we still have to add pairing. The formal generalization of the rigorous HK theorem to paired system has been given in Ref. . In the present work our main objective is to discuss the KS-DFT scheme for the non pairing part, thus we add pairing in a very simple way within the BCS approach. For this we simply take the density dependent delta force defined in Ref  for with the same parameters and in particular with the same cutoff. As far as this amounts to a cutoff of 10 MeV into the continuum for finite nuclei, we have to deal with single-particle energy levels lying in the continuum. We have simulated it by taking in the pairing window all the quasi-bound levels, i.e. the levels retained by the centrifugal (neutrons) and centrifugal plus Coulomb (protons) barriers. This treatment of the continuum works properly, at least for nuclei not far from the stability valley as it has been extensively shown in . In this way we obtain two-neutron () and two-proton () energy separations for magic proton and neutron numbers in quite good agreement with the experiment (see also below).
In our calculations the two-body center of mass correction has been included in the self-consistent calculation using the pocket formula, based on the harmonic oscillator, derived in Ref. which nicely reproduces the exact correction as it has been shown in . Our functional is now fully defined and, henceforth, we call it BCP-functional. Preliminary results based on this functional have been presented in ref. [36, 37].
3 Fitting procedure and general results
In this section we discuss the fitting protocol to determine the open parameters of our functional and the accuracy that can be reached following the scheme described in the previous section. To this purpose we make a comparison between our recent results and the ones that can be obtained with the Gogny functional, which is one of the most accurate phenomenological functional. In this paper we follow a novel strategy, recently proposed , for the fitting procedure to obtain the set of parameters for the surface and spin-orbit part of the functional, that we call the BCP09 functional. Contrary to our previous works we no longer use the binding energies of spherical nuclei for the fit but, following the suggestion of Ref. , we use a set of deformed nuclei carefully chosen in the rare earth, actinide and super-heavy regions instead. The underlying philosophy behind this choice  is that for spherical nuclei the fluctuations in the two most relevant collective variables, namely, the quadrupole deformation and the pairing correlations are large (i.e. the minimum in those variables is broad in some sense) and therefore correlations beyond mean field are important and difficult to evaluate. On the other hand, for deformed nuclei the minimum as a function of the quadrupole degree of freedom is stiffer than its counterpart in spherical nuclei and the additional correlations (rotational energy correction mainly) are not so difficult to compute in a mean field model.
In order to fit the three free parameters of the surface part (, and ) of BCP and the spin-orbit strength we have taken as in  84 well deformed nuclei in the rare earth, actinide and super-heavy regions where the experimental binding energies are known. The four parameters are determined by minimizing the mean square root deviation of the binding energy . The theoretical binding energy has been computed as the Hartree- Fock- Bogoliubov (HFB) mean field energy plus a rotational energy correction, which differs from the standard one  by a phenomenological factor to account for the approximate calculation of the Yoccoz moment of inertia. We also consider an additional correction factor (not relevant in the present calculation) to deal with the weak deformation regime (see, for instance, ). Also a correction to the binding energy ensuing from the finite size of the harmonic oscillator basis used, as estimated in , has been included. Axial symmetry has been preserved in the HFB calculation as the nuclei chosen are not expected to develop triaxiality. On the other hand reflection symmetry is allowed to be broken in the solution of the HFB equation, to allow for octupole deformation in the ground state which is relevant for a few of the actinides considered. The calculation of the 84 ground states for a set of parameters can be carried out in a powerful personal computer in around half an hour. As a consequence of this figure, an unconstrained and blind search of the four parameters leading to the absolute minimum of is a task out of the scope of the present exploratory considerations. Fortunately, the expression for the binding energy depends linearly on three out of the four free parameters (, and ) and the method suggested by Bertsch et al.  in those situation should work well. We have implemented Bertsch’s procedure for the three parameters , and and performed systematic calculations for as a function of the other parameter . Usually, for a given value and with a reasonable choice of starting parameters, the linearized Bertsch’s method leads to a local minimum of in a couple of iterations. The reason is that the correlation matrix of Bertsch’s method has in our case two eigenvalues which are two to three orders of magnitude smaller than the remaining one. This fact points to the existence of only one free parameter (a linear combination of , and with weights corresponding to the components of the eigenvector of the largest eigenvalue).
At this point it has to be mentioned that the fact that the expression of the binding energy is linear in , and by no means imply the value of binding energy to be globally a linear function of those quantities. The reason is that the binding energy is obtained after a self-consistent HFB calculation and the self-consistency breaks the linear dependence of the wave function on the parameters mentioned above. As a consequence, Bertsch method is only valid locally what implies that the procedure will end up in different minimum depending upon the values of the initial parameters. A careful analysis of the results obtained so far indicates that the spin-orbit strength is the most sensitive starting parameter and therefore a search on this degree of freedom should be carried out in addition to the search on .
Using this new protocol we have performed a fit using for the bulk density dependent part of the functional new interpolating polynomials for the nuclear matter equation of state in a wider range of density with respect to the previous ones (BCP1 and BCP2 in refs. [36, 37]). The new fit has still a binding energy per nucleon MeV at saturation density fm (for a compressibility of 220 MeV). For the spin-orbit strength we choose to stay around the MeV value and perform a minimization of for different values of from 0.85 fm up to 1.00 fm. The minimum value of was obtained for fm. The values of the other parameters are MeV, MeV and MeV fm for a value of of 0.545 MeV. In order to check for the suitability of the spin-orbit strength value used we performed another two minimizations of fixing at 0.9 fm but taking MeV fm and 102.51 MeV fm. The values of obtained are 0.569 MeV and 0.566 MeV respectively, supporting the assignment of MeV fm as the value leading to the lowest (this check is incomplete as a more extensive search in the two parameters should be performed; work in this direction is in progress and will be reported in a near future).
In Fig. 2 we show the individual differences in binding energy () between the theoretical results and the experimental ones for the 84 nuclei considered. The quantity has been plotted as a function of the number of neutrons and the values corresponding to the same isotope are joined with lines. Two sets of parameters have been considered, namely one with fm producing the lowest value of (full dots) and other with a slightly higher value MeV and obtained fixing fm (circles). It is worth pointing out that the absolute value of never exceeds 1.2 MeV. We also observe that the change in parameters only produce minor changes in the ’s. Finally, we notice that the actinides and superheavies look a little bit too underbound as compared to the rare earth isotopes. Finding the origin of such a relative underbinding will help to improve the quality of the fit reducing the value of . Although the value of was obtained by considering an optimal set of deformed nuclei ( those with a well established deformation ) it is very likely that the value for a complete set of nuclei from proton drip line to neutron drip line will be promising.
In order to compare the results obtained with BCP 09 with others of other interactions we have performed calculations with the Gogny force and the most recent parametrization D1M that was made up with the idea of producing a mass table. Unfortunately, the details of how the theoretical binding energy was obtained are not detailed enough as to allow a reproduction of them (the main uncertainty is in the quadrupole motion zero point energy correction). We have decided to use the same protocol as with BCP 09, that is, the rotational energy correction and the finite size of the basis effect are included, and we have left out the zero point energy correction mentioned above. As a consequence, our values for of D1M are too high and in order to make a more fair comparison we just shifted all the binding energies computed by a constant quantity (2.8 MeV) obtained as to minimize . With this readjustment we obtain for D1M and the 84 nuclei considered the value MeV which is slightly better than our value but, on the other hand, of similar quality to the one of the BCP 09 functional. In Fig. 3 we have plotted the individual values of both for BCP 09 and Gogny D1M. We notice how the behavior of for the Gogny D1M as a function of for each isotopic chain is different from the one of BCP 09; decreasing with increasing in the case of Gogny D1M whereas it is increasing with increasing for BCP 09. It is also worth to point out that Gogny D1M force produces results in the actinide and super-heavy nuclei which are more spread out than the ones of the BCP 09 functional. By comparing both results we can conclude that, at least for the 84 deformed nuclei considered, the performance of BCP 09 is as good as the one of Gogny D1M and at a fraction of the computational cost (the absence of exchange terms and non local pairing in BCP makes the numerical evaluation of the quantities entering the HFB equation much faster than for the Gogny D1M case). This is an encouraging starting point for a more detailed study of the viability of BCP 09 in the description of low energy nuclear structure.
Next we want to explore the ability of our proposed BCP09 energy density functional to describe nuclear ground-state properties of other nuclei not included in the fitting protocol used to obtain the free parameters that describe the surface and spin-orbit contributions to the functional.
To this end we have computed the ground-state energies of 161 even-even and 306 odd spherical nuclei. These nuclei are chosen to be spherical according to the deformation properties of the compilation of Möller and Nix . To deal with odd nuclei, we have used the blocking approach on top of the BCS calculation. The differences , where are the theoretical predictions of the BCP09 functional and the experimental values taken from , are displayed in Fig. 4. The agreement found between the theoretical prediction and the experimental values is fairly good finding an energy rms 1.3 MeV. It should be mention that both, even-even and odd nuclei, give basically the same rms when they are considered separately. From the figure we can see that the theoretical ground-state energies are scattered, rather symmetrically, around within a window of about 3 MeV with the only exception of the nucleus Er. Combining these values for 467 spherical nuclei with the results obtained for the 84 deformed nuclei used in the fitting protocol, we find a global rms 1.2 MeV, a value that can be improved if the analysis of the ground-state energies include more deformed nuclei.
In Fig. 5 we show the differences between the calculated and experimental charge radii . In this Figure we have displayed the differences for 88 even-even and 111 odd spherical nuclei for which the experimental charge radii are known . From this Figure we see that the BCP09 energy density functional predicts, on the average, smaller charge radii than the experimental values. With few exceptions, the theoretical charge radii are scattered within a window of 0.04 fm around an average values = -0.02 fm. The global quality of this fit of the charge radii is given by rms =0.030 fm and again the rms of the even-even and odd nuclei are basically the same. This global rms can be compared with one provided by the HFB-8 model  which is =0.0275 fm, however computed using 782 experimental data. As compared with the results provided by the earlier BCP2 functional, that included some experimental charge radii in the fit, the values are roughly the same, but in the case of the BCP2 functional the theoretical charge radii are rather overestimated, at least for light nuclei.
After having defined a performing functional in the realm of nuclear masses one has to check its suitability in describing other nuclear properties like deformations. A detailed study of deformations with BCP 09 is underway and will be reported in the future, but in order to give a taste of the results already obtained we show in Fig. 6 the typical outcome of a fission barrier calculation for the nucleus Pu that has also been discussed in Ref.  for early versions of the BCP functional. In fission studies the energy of the system obtained in the HFB framework is obtained as a function of the quadrupole moment by constraining the HFB solution. In this way a potential energy curve (PEC) is obtained that in the present case is displayed in the lower panel of Fig. 6 for the BCP1 and BCP09 functionals and the Gogny D1S force. The BCP1 and BCP09 curves have been shifted by 0.65 MeV and 4.5 MeV respectively as to make the first minimum at b ( 0.26 ) coincide in energy. After the first (ground state) minimum there is a barrier (first fission barrier) that is higher for D1S than for the BCP functionals. The height of this barrier is known to be very sensitive to triaxial effects which are still not accounted for in the present calculation. After the first fission barrier there is an additional excited local minimum at b (0.73) that corresponds to the fission isomer. The excitation energy of the fission isomer obtained with the two BCP functionals is more or less the same with a value of 3.5 MeV. For the Gogny D1S the result is slightly higher (around 4.5 MeV). The following maximum corresponds to the second fission barrier which is definitely much higher for Gogny D1S than for the BCP functional. After the second fission barrier the typical decrease of the energy consequence of the Coulomb repulsion between the nascent fragments (the real shape of the nucleus is depicted at the corresponding value as a contour plot of the real density at a value of the contour of 0.08 fm). As the two fragments are still connected by the neck the nuclear part of the interaction still has an effect, explaining the differences observed between the three calculations. Finally, in the upper panel of Fig. 6 the octupole () and hexadecapole () moments of the mass distribution obtained with the three functionals/interactions are shown. The values lie on top of each others indicating (as the multipole moments considered are the most relevant) that the shapes of the nucleus along the whole fission path are essentially the same for the three calculations.
4 Open problems
The problems to be clarified in the EDF approach can be in general classified in two categories. On one hand, the use of the nuclear matter EOS must be taken with care since the microscopic calculations cannot reach the accuracy usually required in nuclear structure calculations and contains intrinsic and unavoidable uncertainties discussed below. On the other hand, finite size effects must be included in the functional with a minimal set of parameters. The procedure of doing that is not straightforward and requires extensive analysis.
4.1 Nuclear matter EOS
Following the KS scheme, the nuclear EOS calculated microscopically should be included in the EDF mainly by an analytic fit that reproduces accurately both the nuclear matter and neutron matter results. In this way the saturation point, as well the whole density dependence of the binding energy, can be incorporated in the EDF. However, the microscopic calculations cannot be considered without any uncertainty. Besides the point of numerical accuracy, it is well known that it is necessary to introduce three-body forces in the microscopic calculations in order to reproduce the correct saturation point, which however contains necessarily a phenomenological uncertainty. Therefore the TBF can be adjusted to reproduce different saturation points. Furthermore the microscopic data do not fix uniquely the fitted EOS. We found that the quality of the nuclear data fitting is very sensitive to the saturation point of the interpolated EOS. The optimal value of the energy per particle at saturation seems to be close to -16 MeV, but even a shift as small as 0.1 MeV of this quantity can deteriorate appreciably the quality of the nuclear data fit. This should be a general feature, independent of the particular EDF that is used. Less sensitivity seems to be present for the saturation density. For the EOS obtained from the BHF theory, the TBF contribution below saturation is quite small and the BHF is expected to be accurate at low density, so that the main problem can be considered as how to fix the saturation point.
Of course the just discussed sensitivity on the saturation point could be considered a not serious problem, in view of the phenomenological uncertainty that affects its position.
All these problems appear absent for the neutron matter EOS, for which the TBF contribution is smaller and correlation energy is substantially reduced.
A separate discussion is needed for the type and strength of the TBF that are used in microscopic nuclear matter calculations. The numerical calculations of the binding energy of nuclear few-body systems, like triton and alpha particle, can be performed numerically very accurately. If realistic two-body forces, like the Argonne v potential, are used the binding energies turn out systematically underestimated. To remedy to this drawback usually one adds a semi-phenomenological TBF, like the Urbana model , that is adjusted to reproduce the experimental binding energies. The unpleasant discovery is that the use of the same TBF in microscopic calculations for nuclear matter does not produce a good saturation point. The TBF appear too repulsive. Indeed, in general semi-phenomenological TBF’s contain a repulsive and attractive part, and the standard procedure is to reduce the repulsive component in order to get a good saturation point. At phenomenological level this means that one needs at least a four-body force in a non-relativistic many-body scheme.
It has to be mentioned in this respect that non-local two-body forces  which are able to reproduce the binding energy of few-nucleons systems fail to reproduce the correct saturation point, that then turns out to be at too high density and too low energy .
Recently  it has been shown that one can get a reasonable saturation point by keeping the TBF obtained by fitting few-body systems if for the effective two-body force the so-called V interaction is used. The latter is obtained by projecting out the high momenta of the bare NN interaction and renormalizing accordingly the interaction at low momenta. The resulting interaction is then phase-equivalent to the original bare NN interaction up to the momentum cutoff, but it is much softer, and therefore it can be treated perturbatively. Nuclear matter calculations with V are not saturating, most probably because they misses the effects of the Pauli principle and of the self-consistent single particle potential, which are distint features of the G-matrix. The final saturation point is a consequence of the compensation between the too large attraction of V and the too large repulsion of the TBF. All these uncertainties are of course embodied in the final EOS and affect the detailed properties of the resulting EDF.
Unfortunately there is a more basic problem for the low density EOS. At low enough density symmetric matter is expected to be unstable towards cluster formation. The picture of an homogeneous matter cannot be kept at very low density. The real EOS is therefore altered by the appearance of clusters of different sizes. However the EOS which includes cluster formation cannot be used in the EDF, since such a proliferation of light nuclei cannot be present in the low density region of finite nuclei. This is a general problem for any EDF. The absence of clusters at the nuclear surface is due to the small value of the diffusion of the density profile, which prevents any long range correlation to dominate. This fact offers a partial justification of taking the homogeneous matter EOS at the BHF level. Indeed, short range correlations are well taken into account by the BHF procedure, and therefore the corresponding EOS should be able to describe the properties of the bulk contribution to the nuclear EDF even in the surface region.
This problem still needs further analysis before it can be considered satisfactorily solved.
4.2 Incorporating the nuclear matter EOS in the EDF
In the standard KS procedure the kinetic energy contribution is taken at the independent particle level with an effective mass equal to the bare one. This is in agreement with the BHF scheme, where the kinetic energy part is also kept at the free value and the whole correlation contribution is included in the interaction energy part coming from the G-matrix. This means that also the correlated part of the kinetic energy is included in the G-matrix contribution. Indeed the momentum dependence of the G-matrix is clearly the origin of the deviation of the effective mass value from the bare one. It would not be difficult to separate this effect from the G-matrix contribution, so that the kinetic energy would include a (density dependent) effective mass, while leaving the rest of correlations as a genuine interaction energy. This leaves some freedom to the way the nuclear matter EOS is actually included in the EDF, which allows to go beyond the usual KS scheme. It can be expected that different procedures are not necessarily equivalent. This problem has still to be analyzed in detail.
4.3 Finite size effects
The density matrix expansion of Negele and Vautherin  suggests that the surface contribution to the binding energy of nuclei may also be extracted from nuclear matter calculations. This would be very welcome, since this would further reduce the number of adjustable parameters and link the functional even more to the microscopic approach. Essentially only the spin-orbit and pairing contributions would then need some phenomenological adjustments. However, the G-matrix also has spin-spin interaction terms and, in principle, the spin orbit also could be extracted if the vector part of the one-body density matrix is properly taken into account . To what extent such a goal can be achieved remains an open problem. Let us, however, point to some possible improvements over the past procedures. In first place we want to point out that it is our believe that neither the DME of Negele- Vautherin, nor the one of Campi-Bouissy  are the optimal procedures. We have long standing experience with the semiclassical DME of Wigner and Kirkwood which is based on a systematic expansion in powers of . Most accurate results are obtained, at least for the scalar part of the one- body density matrix, with this method . A second problem stems from the way how the G-matrix of an infinite matter calculation is used for a finite nucleus. The standard procedure is the LDA, that is one replaces the dependence by a density dependence via the standard relation . The dependence enters mostly into the Pauli operator . However, dependence also can be in the self consistent single particle energies. We want to point out that just replacing by is not necessarily the best procedure. Our argument goes as follows. The two particle propagator entering the G-matrix equation
has an obvious semiclassical, i.e. limit, namely we just replace the single particle Hamiltonians by their classical counterparts with the phase space variables of position and momenta. With this approximation the two particles move according to their respective classical Hamiltonians, each. Supposing for simplicity a local mean field (the argument can be generalized to nonlocal mean field potentials), we see that in the semiclassical two body propagator two local chemical potentials appear: and () and, therefore, two Fermi momenta and, thus, two densities at positions appear in the semiclassical propagator and correspondingly in the G-matrix. Since the effective range of the interaction is over 2 fm wide, this non locality in the densities may be of quite some importance. One may call this the Non Local Density Approximation (NLDA). The standard LDA is recovered in putting . This NLDA effect may be more important, or at least of equal importance, than to keep gradient expansions of the density itself. It could be interesting to test this. Since finite size effects are surface effects, one may argue that only relatively low densities will be involved. It has recently been shown, that at least in neutron matter for low densities a separable force works very well for the G-matrix , reducing very much the numerical effort. So eventually one could combine the NLDA approach with a separable force to get to a surface term. Of course, one will not imagine that with such a procedure one will hit the good result right on the spot. But even if the result is only semi quantitatively correct, it still may lead to a reduction of the open parameters and to more insight into the underlying physics. This may be then a first step in the direction to get everything from the microscopic input, i.e. from the underlying bare forces, also for finite nuclei.
5 Conclusions and prospects
We have presented a formulation of the energy functional method that enables to keep, to a certain extent, the connection with the nuclear matter Equation of State as calculated microscopically from many-body theory and realistic nucleon-nucleon interactions. The functional is constructed following closely the Kohn and Sham method, where the bulk part of the functional is taken directly from the microscopic calculations and kept fixed in the fitting procedure of the adjustable parameters of the functional. The latter include the surface part of the functional and the spin-orbit strength, ending up with a total of four free parameters. The pairing interaction and the corresponding strength has been taken from the simplest standard choices. Despite the analysis was not carried out up to the optimal level, leaving space for refinements, the accuracy of the results can compete with the one obtained with the best purely phenomenological functional like Gogny. These promising results open the possibility of building an energy density functional closely connected with the bare nucleon-nucleon forces. Developments in this direction could be obtained along the density matrix expansion method of ref. [5, 9] or its refinements, which could further reduce the number of adjustable parameters. However, we pointed out the various open problems that must be solved before this project can be completed. The incorporation of the nuclear EOS in the density functional can be done in different ways, and this leaves some ambiguity in the direction of a universal nuclear functional. Furthermore the very low density part of the EOS cannot be taken literarily from microscopic calculations of homogeneous nuclear matter. In fact cluster formation is expected to dominate this region of the EOS. This appears a serious problem for any density functional theory which is constructed in such a way to be compatible with the nuclear matter case. Further studies and analysis are needed to proceed further in this project on the nuclear density functional.
Work supported in part by MICINN (Spain) grants FPA2007-66069 and FPA2008-03865-E/IN2P3, and by the Consolider-Ingenio 2010 program (Spain) CPAN, CSD2007-00042. X. V. also acknowledges the support from FIS2008-01661 (Spain and FEDER) and 2009SGR-1289 (Spain). This work was supported by CompStar, a Research Networking Programme of the European Science Foundation.
- D. Vautherin, D.M. Brink, Phys. Rev. C5 (1972) 626.
- J. Dechargé, D. Gogny, Phys. Rev. C21 (1980) 1568.
- E. Chabanat, P. Bonche, P. Haensel, J. Mayer and R. Schaeffer, Nucl. Phys. A627 (1997) 710; A635 (1998) 231.
- D. Q. Lamb, J. M. Lattimer, C. J. Pethick, D. G. Ravenhall, Nucl. Phys. A360, 459 (1984).
- J. W. Negele and D. Vautherin, Phys. Rev. C5, 1472 (1972).
- X. Campi and A. Bouyssy, Phys.Lett. B73, 263 (1978).
- J. Meyer, J.Bartel, M. Brack, P. Quentin, S. Aicher, Phys.Lett. B172, 122 (1986).
- V.B. Soubbotin, X. Viñas, Nucl.Phys. A665, 291 (2000)
- T. Duguet and T. Lesinski, arXiv:0907.1043 (2009).
- W. Kohn and L. J. Sham, Phys. Rev. 140, A1133 (1965).
- R. O. Jones and O. Gunnarsson, Rev. Mod. Phys. 61, 689 (1989).
- H. Eschrig, The Fundamentals of Density Functional Theory (B. G. Teubner, Stuttgart, 1996).
- A. E. Mattsson, Science 298, 759 (2002).
- J. P. Perdew and K.Schmidt, Density Functional Theory and Its Applications to Materials, V. Van Doren et al. Eds. CP577 (American Institute of Physics, Melville NY, 2001)
- Jianmin Tao, J. P. Perdew, V. N. Staroverov and G. F. Scuseria, Phys. Rev. Lett. 91, 146401 (2003).
- Jianmin Tao and J. P. Perdew, Phys. Rev. Lett.95, 196403 (2005).
- S.A. Fayans, JETP Lett. 68, 109 (1998); S.A.Fayans,D.Zawischa, Int.J.Mod.Phys. B15, 1684 (2001).
- T. Kreibich and E. K. U. Gross, Phys. Rev. Lett. 86, 2984 (2001).
- J. Engel, Phys. Rev. C75, 014306 (2007).
- J. Dobackzewski, M. V. Stoitsov, W. Nazarewiczand P.-G. Reinhard, Phys. Rev. C76, 054315 (2007).
- B. G. Giraud, B. K. Jenning ands B. R. Barret, arXiv:0707.3099 (2007); B. G. Giraud, arXiv:0707.3901 (2007)
- First ”FIDIPRO-JSPS Workshop on Energy Density Functionals in Nuclei”, Jyvaskyla, October 25-27, 2007; ”Formal aspects of Nuclear Energy Density Functional Methods”, Saclay, November 5-7, 2007.
- P. Hohenberg and W. Kohn, Phys. Rev. 136, B864 (1964).
- J. Engel, Lecture Notes in Physics 620, 56 (2003), C. Fiolhais, F. Nogueira and M. Marques Eds. (Springer-Verlag, Berlin, 2003).
- M. Baldo, C. Maieron, P. Schuck and X. Viñas, Nucl. Phys. A736, 241 (2004), and references therein.
- J. Margueron, E. van Dalen and C. Fuchs, Phys. Rev. C76, 034309 (2007).
- A. Akmal, V. R. Pandharipande and D. G. Ravenhall, Phys. Rev. C58, 1804 (1998).
- P. Bonche, S. Koonin and J. W. Negele, Phys. Rev. C13, 1226 (1976).
- A. S. Umar, M. R. Strayer, P.-G. Reinhard, K. T. R. Davies and J.-S. Lee, Phys. Rev. C40, 706 (1989).
- S. A. Fayans et al Nucl. Phys. 676, 49 (2000).
- L. N. Oliveira, E. K. U. Gross and W. Kohn, Phys. Rev. Lett. 60, 2430 (1988).
- E. Garrido, P. Sarriguren, E. Moya de Guerra, and P. Schuck, Phys. Rev. C 60, 064312 (1999)
- M. Del Estal, M. Centelles, X. Viñas, Phys. Rev. C63, 044321 (2001).
- M.N. Butler, D.W.L. Sprung and J. Martorell, Nucl. Phys. A422, 157 (1984).
- V.B. Soubbotin, V.I. Tselyaev and X. Viñas, Phys. Rev. C67 (2003) 014324.
- M. Baldo, P. Schuck and X. Viñas, Phys. Lett. B663, 390 (2008)
- L.M. Robledo, M. Baldo, P. Schuck and X. Viñas, Phys. Rev. C77, 051301(R) (2008)
- T. Niksic, D. Vretenar and P. Ring, Phys. Rev. C78, 034318 (2008)
- P. Ring and P. Schuck, ”The nuclear many body problem”, (Springer Berlag, Berlin, 1980)
- J.L. Egido and L.M. Robledo, Lecture Notes in Physics, 641, 269 (2004)
- S. Hilaire and M. Girod, Eur. Phys. J. A33, 237 (2007)
- G. F. Bertsch, B. Sabbey, and M. Uusnäkki, Phys. Rev. C71, 054311 (2005)
- S. Goriely, S. Hilaire, M. Girod, and S. Péru, Phys. Rev. Lett. 102, 242501 (2009)
- P. Möller, J.R. Nix and K.L. Kratz, At. Data Nucl. Data Tables 66,131 (1997)
- G. Audi and A.H. Waspra, Nucl. Phys. A729, 337 (2003)
- I. Angeli, At. Data Nucl. Data Tables 87, 185 (2004)
- F. Buchinger and J.M. Pearson, Phys. Rev. C72, 057305 (2005)
- J. Carlson, V.R. Pandharipande and R.B. Wiringa, Nucl. Phys. A401, 59 (1983)
- P. Doleschall and I. Borbely, Phys. Rev. C62, 054004 (2000); P. Doleschall, I. Borbely, Z. Papp and W. Plessas, Phys. Rev. C67, 064005 (2003); P. Doleschall, Phys. Rev. C69, 054001 (2004)
- M. Baldo and C. Maieron, Phys. Rev. C72, 034005 (2005)
- S.K Bogner, R.J. Furnstahl, A. Nogga and A. Schwenk, arXiv:0903.3366 (2009).
- M. Baldo and C. Maieron, Phys. Rev. C77, 015801 (2008)