Nuclear density functional theory

The goal of nuclear structure physics is to provide a complete understanding of the static properties of atomic nuclei, their excitation spectra, their response to external fields and their decays. While it is hard to achieve these goals within a single framework, so that there is no nuclear ‘standard model’, it is clear that nuclear Density Functional Theory (DFT) has probably the widest range of applicability so far. In this paper, we try to put DFT in a broader context, with frequent comparisons to electronic DFT. We also include a discussion of the relationships with ab initio methods and Effective Field Theories (EFTs) in general, as well as a short survey of the quite large number of applications. Although written with a personal and possibly biased perspective, the paper aims at fostering cross-fertilizations with other domains of science. ARTICLE HISTORY Received 29 December 2019 Accepted 3 March 2020


Introduction
Nuclear physics has the well-deserved reputation of being an intricate, demanding, and sometimes painful subject in physics.
At the phenomenological level, one can start by considering the huge variety of properties that nuclear systems display. Their existence or nonexistence is already a non-trivial property. Today we know % 250 stable nuclei, that is, given combinations of Z protons and N neutrons (that can be arranged in a plane like in Figure 1) that have infinite lifetime. Other combinations can be bound and yet decay into other forms on a very long or very short timescale; very short-lived nuclei are hard to detect experimentally but we are experiencing a continuous progress in this respect, as testified by the fact that 3302 stable and unstable nuclei have been reported to exist at the end of 2018 [1], and that 13 nuclei/year have been discovered on average in the years 2103-2016 and 34 in 2017 [2]. Thus, we can affirm that the nuclear landscape is changing continuously.
Nuclei are described, still, in many introductory courses [3,4] as liquid drops. Their binding energy BE is defined as BEðN; ZÞ ¼ MðN; ZÞc 2 À Zm p c 2 À Nm n c 2 ; (1) where the nuclear mass M and the proton and neutron masses m p and m n are introduced. For typical medium-heavy nuclei of mass number A ¼ N þ Z, the binding energy per nucleon BE=A is % 8 MeV. This value is small with respect to m p or m n : in other words, the mass of the nucleus is mainly due to the mass of the constituents. 1 The fact that the binding energy per nucleon in standard nuclei is pretty constant tells us that each nucleon interacts with nearest neighbours and the number of interactions is not proportional to the number of nucleon pairs AðA À 1Þ, but simply to A. This is called 'saturation' of the nuclear force. 2 Saturation produces stable nuclei that have an inner density ρ (sum of proton and neutron densities, ρ p and ρ n ) which is also pretty constant. This so-called saturation density is ρ 0 % 0.16 fm À3 . In such conditions, the average distance between nucleons is somewhat larger than the range of the nuclear ZÞ plane] is shown together with some predictions of the proton and neutron drip lines. The drip lines correspond to the limits of nuclear existence, and their precise definition can be found around Equation (3) and (4) below. These predictions have been drawn by using the tools available at the URL http://massexplorer.frib.msu.edu/ content/Plotting_Tools.html and including results from all the models. On the right side, the extrapolation from finite nuclei to neutron stars is shown. As briefly mentioned in the main text, neutron stars are essentially systems of % 10 55À56 neutrons that are bound by gravity.
force. This is mainly an effect of Pauli principle correlations. The nucleon mean free path in the nucleus has been found (see, e.g. Figure 7 of [5]) to be larger than the nuclear radius, that is, with r 0 % 1:2 fm. It must be noted that this large value of the nucleon mean free path is consistent with the well-known fact that, to a first approximation, nuclei can be described in terms of nucleons that move independently in an average potential. In fact, nuclei display a clear shell structure characterised by neutron and proton orbitals; closed-shell nuclei are known as magic nuclei. We may be tempted, then, to estimate the typical kinetic energy of nucleons at saturation density by using the Fermi gas formula, and we obtain % 35 MeV. From this value and the value of the mean free path, we can also infer that the time between two successive nuclear interactions is at least of the order of 10 À22 s. We shall use this estimate below. The picture of independent nucleons in a spherical average potential is a first and rough one. Quadrupole and pairing correlations play important roles in atomic nuclei. Away from closed neutron and proton shells, it is known since the early days of nuclear physics [6] that many nuclei are conveniently described by assuming a quadrupole-deformed intrinsic shape in which nucleons are confined, and that undergoes rotational motion with respect to the laboratory frame. 3 Evidences of nuclei with octupole deformation (i.e. pear-like shapes) are more recent [7]. The quest for more exotic shapes is still ongoing [8]. Non-spherical shapes arise because many nuclei gain energy when nucleons are placed in the orbitals associated with such shapes and are also favoured by nucleon-nucleon (in particular, protonneutron [6,9]) correlations. An example is shown in Figure 2. Pairing  24 Mg, calculated using one of the models already employed to draw the left panel of Figure 1, as a function of the quadrupole deformation parameter β 2 ; 3 4 ffiffi π 5 p δ, where δ is the relative difference between the major and minor semi-axis of the ellipsoid. (Right) Total ground-state density (in fm À3 ) as a function of cylindrical coordinates r; z (z is along the symmetry axis of the ellipsoid).
correlations are responsible for the fact that nucleons outside closed shells exhibit, as a rule, a superfluid behaviour; at the same time, we have to add that several (and by far not all explored yet) superfluid phases exist in nuclear and neutron matter [10][11][12], in analogy with the case of electronic superconductors. We shall relate, in Sec. 3, deformation and pairing to spontaneous symmetry breaking.
Although 8 MeV of binding energy per nucleon is a typical value, nuclei can become weakly bound. As we mentioned at the start of the manuscript, we still do not know all nuclei that exist in nature. The limit of nuclear existence is reached when hitting the so-called drip lines. The proton and neutron drip lines correspond, more precisely, to the loci where the separation energy S for that kind of particle, defined by S p ; BEðN; ZÞ À BEðN; Z À 1Þ; (3) and by S n ; BEðN; ZÞ À BEðN À 1; ZÞ; changes its sign [the binding energies BE are defined in Equation (1)]. The lines depicted in Figure 1 are obtained by using some of the theoretical models that will be discussed below. Our current uncertainties on the drip lines are nicely discussed in Refs [13,14]. The neutron (proton) drip lines are reached by increasing (decreasing) the number of neutrons with respect to the typical values for stable nuclei, that are those which maximise the binding energy (1). Protons and neutrons can be seen as two states of the same particle with different z-projections of the so-called isospin. The strong interaction commutes with isospin, that is, preserves isospin symmetry; nonetheless, this interaction is isospin-dependent and stronger among protons and neutrons, than among particles of the same kind. This difference is also larger than the effect of the Coulomb interaction, that breaks isospin symmetry. We now briefly discuss some consequences of these facts.
The Coulomb repulsion prevents medium-heavy nuclei from having substantial proton excess; nuclei close to the proton drip line are 'neutron deficient', in the sense that they have less neutrons than stable nuclei in the same mass region. On the other hand, nuclei with substantial neutron excess exist, and the neutron drip line extends much further apart from the stable nuclei than the proton drip line. Due to the strong proton-neutron component of the nuclear force, nuclei close to the proton drip line display strongly bound neutron orbitals and proton orbitals close to the continuum; for the very same reason, close to the neutron drip line one expects strongly bound proton orbitals and neutron orbitals close to the continuum. Moreover, in this latter case, when neutron binding is weak and the angular momenta of the last neutrons are small, the wave functions are not limited by the centrifugal barrier and extend far from the nucleus, giving rise to regions where the density is smaller than the typical saturation density ρ 0 . These regions form a so-called neutron 'skin' or 'halo' of dilute neutron matter. 'Exotic' nuclei with skin or halo are of great interest for a number of reasons. Their diffuse surface may change the sequence of the nucleon orbitals and even the magic numbers [15], and the excitation modes of exotic nuclei may be different from the well-known modes in standard nuclei [16].
Low-density neutron matter can become quite different from ordinary nuclear matter. In fact, very dilute neutron matter is expected to be the prototype of a universal system, the so-called unitary gas. To understand this point, we remind that, in a gas of particles, the low-energy elastic cross section can be written as [17] σ ¼ 4π where a is the s-wave scattering length and k is the relative momentum. One usually identifies a 'natural' case in which the relative momentum k can be close to zero and k ( jaj À1 , so that the cross section is σ ¼ 4πa 2 . The neutron-neutron scattering length is 'unnaturally' large, a ¼ À18:7 fm [18]. Since in a uniform system k F ,ρ 1=3 , there exist finite, observable values of the density for which jaj À1 ( k and one reaches the so-called unitary limit, 4 in which the cross section attains the value σ ¼ 4π k 2 , which is independent of the system under study. Other many-body systems may be close to the unitary limit. In fact, it has been suggested recently that cold atoms may explore the whole regime from unitarity to naturalness, and that close to unitarity measurements in cold atoms may actually also provide the neutron matter equation of state (EoS) [19,20].
A recent review about neutron matter from low to high density can be found in Ref [21]. In fact, the properties of high-density nuclear matter are even far less known. The interest in these extreme conditions stems from the fact that, as the right part of Figure 1 hints, neutron matter does exist in nature: neutron stars [22,23] contain mainly neutrons, with just enough protons and electrons so that the system is globally neutral and in thermodynamical equilibrium with respect to weak processes. Neutron stars can be bound due to gravity (cf. the pedagogical exercise at p. 226 of [3]): the socalled Tolman-Oppenheimer-Volkov (TOV) equation assumes that neutron stars are spherical, and that each spherical shell is in hydrostatic equilibrium as a result from the balance between the gravitational attraction and the nuclear repulsion (that arises since the inner shell has larger density ρ than the outer shell and, thus, the pressure P from inside is larger than the pressure from outside). The TOV equation can be solved if a nuclear model provides the EoS in terms of the function PðρÞ.
Terrestrial nuclei only inform us about the properties of nuclear matter close to ρ % ρ 0 , and with a number of neutrons that does not exceed too much the number of protons. Then, extrapolations to neutron star matter are plagued by huge uncertainties. At the surface of the star (the outer crust, at very low density) bound nuclei immersed in an electron gas are the obvious stable form that minimises the energy, but the transitions that take place as the density increases are not well known. There is consensus that, qualitatively, one should move to bound nuclei surrounded by free neutrons (in the so-called inner crust) and then to various forms of clusterised neutron matter, before reaching the uniform neutron matter phase around ρ 0 in the outer core. In the so-called inner core, namely for densities larger than ρ 0 , the formation of new particles, like hyperons, 5 can be considered natural in keeping with the fact that the energy density allows such creation. Nobody knows, at present, how large the density can be at the centre of the star, and if a gas of deconfined quarks can be formed or not. See Figures 1 and 3 in keeping with this description of the layers of a neutron star.
In summary, the constituents of atomic nuclei, protons and neutrons, can live in very different conditions of total density, ρ, and asymmetry, β ; ρ n À ρ p =ρ. This gives rise to an extremely rich host of phenomena, in terms of overall shapes of the bound systems that are formed, shell  [77], where a description of the precise models behind this representation is provided. From the scale on the right side, it is clear that phases similar to the halo of exotic nuclei or to dilute nuclear matter are encountered. structure, energetics, etc. The exploration, and possibly the design, of properties of matter at this scale (fm and MeV) is not less rich and is definitely more challenging than doing the same at the electronic, nanometric, or mesoscopic scale. Some of the differences and analogies with molecules or cold atoms have been briefly mentioned above. Actually, we have only so far addressed the phenomenology of ground states. As energy and angular momentum are given to nuclei, one can explore an even richer world of complex spectra, nuclear reactions and transitions. Nuclei can respond to three of the four elementary forces (strong, electromagnetic and weak) and be laboratories to test their properties, namely, they can absorb and emit electromagnetic radiation, scatter or absorb leptons as well as undergo βdecay, scatter and emit nucleons or other hadrons.
The lowest excited state can be as low as 7.6 eV (in 229 Th) or as large as 6 MeV (in 16 O). Low-lying collective rotations, that is, states in which a deformed nucleus rotates as a whole, are characterised by slow times (10 À20 -10 À19 s) as compared with those of successive nuclear collisions (10 À22 s as we discussed above). As a consequence, a kind of Born-Oppenheimer factorization can be applied to separate the collective motion from that of the single nucleons. As the angular momentum and rotational energy increase, nuclei undergo several transitions to exotic states before, eventually, fissioning [3,6]. Spherical nuclei cannot rotate but have a rich vibrational spectrum. The vibrational states have energies of the order of % MeV or larger; at % tens of MeV, highly collective motion takes place under the form of the so-called 'giant resonances' [24]. Here, the separation of scales does not take place.
All these features should make nuclear structure of interest for condensed matter physicists, chemists, biologists. Fostering the understanding of different types of many-body systems through mutual exchanges between different communities has been always a valuable strategy and should not be abandoned. In this spirit, this paper is devoted to nuclear applications of Density Functional Theory (DFT). One reason is certainly the fact that DFT has a broader range of applicability, as we discuss in the next Section, with respect to other theoretical approaches to nuclei. The other reason is that within DFT an unified language can be found to exchange ideas with scientists dealing with electronic systems. The reader who is interested in going deeper into nuclear DFT can exploit the excellent and recent book edited by N. Schunck [25].

Observables, energy scales and theories
Thinking a theorist, one of the first obvious statements is that nuclei are hard to study because they are many-body systems and the inter-particle interaction is very involved. There is, anyway, a steady progress in the socalled ab initio methods [26][27][28]. This wording refers to theoretical approaches where controlled, in principle exact, many-body approximation schemes are implemented in conjunction with an interaction between nucleons that produces accurate results when used to solve the two-body problem [properties of the deuteron and nucleon-nucleon (NN) scattering], the three-body and possibly other few-body problems. It is well known that one of the main features of the nuclear many-body problem is that threebody forces are very relevant and many-body forces cannot be ruled out.
While electrons are, to the best of our knowledge, structureless and pointlike, and we know well the electromagnetic interaction among them, this is not the case for nucleons. Everybody is inclined to think that the forces among nucleons should be viewed as residual forces between systems of bound quarks. The analogy with the Van der Waals residual forces between complex molecules could be invoked here, but it is of no help. Nonpractitioners may ask straightforwardly 'Why don't we calculate nuclear properties starting from the fundamental theory of strong interactions?'. Quantum Chromo Dynamics (QCD) is the strong interaction theory in the same manner as Quantum Electrodynamics (QED) is the fundamental theory for electromagnetism. But at variance with QED, whose coupling constant α ¼ 1 137 is small at the energies of interest, QCD is not perturbative at low energies as the corresponding coupling constant becomes large [29].
A non-perturbative solution for the case of baryon-baryon interactions can in principle be obtained by means of lattice QCD calculations; however, in the case of the NN interaction, these calculations have not yet been able to provide reliable results (see [30] and references therein for a discussion of the difficulties). Most of the ab initio practitioners have turned their attention to effective chiral Lagrangians. These include the quantum fields associated with the nucleons, and this feature is shared by other QFT Lagrangians that also include effective mesons like pions or heavy ones, and their coupling with the nucleons. However, unique aspects of chiral Lagrangians are the systematic application of the separation of energy scales as well as, obviously, chiral symmetry.
These concepts have been introduced by S. Weinberg [31], who proposed the equivalence between chiral Lagrangians and the QCD Lagrangian at low energy. Low energy means here 'low with respect to the energy scale associated with the mass of the nucleon and of the heavy mesons' (Λ χ % GeV). If we restrict our considerations to the u and d quarks, their masses are very light (of the order of MeV), and if they are set to zero the QCD Lagrangian becomes invariant under chiral transformations. The fact that this symmetry is partially kept, and partially spontaneously broken in the QCD ground state is confirmed by the hadron spectrum (cf., e.g. the discussion in Ref [32]). Moreover, it is known that the presence of a continuous symmetry which is spontaneously broken dictates the existence of a massless boson due to the Goldstone theorem [33]. The pion is not massless, due to the approximate nature of this symmetry, and yet its mass ( % 100 MeV) is much smaller than those of the other mesons. In such a situation, Weinberg's idea is that the most general Lagrangian that includes nucleons and pions (the light, almost massless Goldstone bosons), and is consistent with chiral symmetry, should provide the same results as the QCD Lagrangian, at least in the low-energy regime up to Λ χ .
The strategy of an effective field theory (EFT) based on these ideas (in short, chiral EFT or χEFT) consists in writing the chiral Lagrangian including nucleons and pions. The contributions to the NN (NNN, etc. . . . ) force emerge from the set of Feynman diagrams with two (three, etc. . . . ) nucleon asymptotic states. For the theory to be tractable, these diagrams must be ordered in a hierarchy according to some perturbative parameter. The widely used choice for this parameter is Q Λ χ , where Q is the typical nucleon momentum. This choice is not completely obvious and is debated in the literature [34]. The power ν at is stopped defines lowest-order (LO), next-tolowest-order (NLO), etc. . . . potentials. Details are discussed in [27,28,32]. The existence of a NN bound state (the deuteron) shows that this system is not perturbative. Weinberg (see [35] and references therein) suggested that one can obtain a sensible NN potential by applying perturbation theory with (only) nucleon intermediate states. Therefore, chiral NN potentials are constructed solving the Lippmann-Schwinger equation on top of a LO, NLO, etc. . . . formulation. It has to be noted that this requires introducing a regulator (or ultraviolet cutoff) Λ that avoids divergences. This cutoff is a parameter that should not be confused with the energy scale Λ χ . Chiral NN potentials can at present be tested against the nuclear phenomenology by using many-body methods that, albeit approximate, have control on the nature and the quantitative impact of the approximations, and can in principle be improved up to exact results, as we mentioned at the start of this Section. These methods include Quantum Monte Carlo (QMC) approaches [36,37], the In-Medium Similarity Renormalization Group (IM-SRG) method [38], the Coupled Cluster (CC) approach [39], the Self-Consistent Green's Function Theory (SCGFT) [40], nonrelativistic and covariant Brueckner-Hartree-Fock (BHF) theory [41], and the No-Core Shell Model [42].
So far, the results of these ab initio approaches are sensitive to the many-body scheme, to the order at which the chiral potential is determined (LO, NLO, etc. . . . ), but also to the cutoff Λ. This is the open issue for researchers in this field.
We shall quote a few paradigmatic examples along this line. As for finite nuclei, Ref [43] highlights that chiral potentials within QMC calculations provide results for binding energies and radii of light systems (up to 16 O) that are characterised by uncertainties related to the order truncation and to the cutoff which is employed. The same kind of uncertainties also shows up in QMC calculations for neutron matter [44]. IM-SRG calculations have also clearly shown, so far, a dependence of the results in medium-heavy nuclei on the resolution scale Λ used to build the Hamiltonian [45]. This is a consequence of the difficulty in reproducing the empirical saturation point of nuclear matter. This reproduction remains a challenge for chiral potentials, and may require an expansion to higher order in the perturbative parameter as it has been done so far, as the very recent many-body perturbation theory results based on N 3 LO potentials seem to hint [46]. All these issues can be looked at from a complementary point of view: the chiral potential NLO sat gives successful results, and at the same time it betrays the original spirit of ab initio calculations by including data from medium-heavy nuclei in the parameter fit. The fit has something in common with that of nuclear Energy Density Functionals (EDFs), that we will now discuss.

Nuclear DFT: differences with the case of electronic DFT
DFT is based on the fundamental theorem that has been proposed by P. Hohenberg and W. Kohn back in 1964 [47]. The theorem asserts that, for a system of fermions that are subject to an external potential v ext , the total energy can be written as a functional of the particle density ρðrÞ. More precisely, In fact, the first equality is the definition of the total energy, assuming thatĤ is the Hamiltonian that governs the many-fermion system, to which the external potential is added. The theorem is related to the second equality, where the contribution of the external potential is singled out: the functional F should be universal, namely, it is related to the kind of fermions and not to v ext . Equation (6) also implies that the total functional E v displays a minimum at the exact ground-state density, and its value is the exact ground-state energy.
It is hard to fail appreciating the asset provided by Equation (6) and the theorem behind it. For a N-body system, one should in principle deal with a complex function of 3N variables. 6 This is unmanageable except for small N, but we have just discovered that dealing with a real function of 3 variables is enough! In fact, as discussed in textbooks, every observable related to the ground-state of the system can be deduced from the density alone. The drawback comes from the fact that the proof of the HK theorem is not a constructive proof; in other words, it does not give any clue on how to build the universal functional F and/or how to start from a simple ansatz and improve it systematically.
A commonly used scheme is the one introduced by W. Kohn and L.J. Sham [48], in which one exploits the fact that the density ρðrÞ can be represented in terms of auxiliary single-particle wave functions (orbitals) ϕ j ðrÞ [49], namely The kinetic energy is then written as that is, it is a functional of τðrÞ ¼ P j jÑϕ j ðrÞj 2 as one easily can realise if an integration by parts is performed 7 starting from [8]. In the Kohn-Sham scheme, not only the kinetic energy but also the direct Hartree energy (the classical electrostatic interaction energy between distributed charges) is singled out and one is left with an unknown part of the energy functional which is dubbed exchange-correlation. In the case of Coulomb systems, the success of this scheme is also due to the compensation of errors in the exchange and correlation contributions.
In nuclear systems, there are a number of differences that we wish to stress in what follows.
[1] There a historical difference. Nuclear EDFs were born as Hartree-Fock expectation values of Hamiltonians, and gradually developed into something more similar to a real EDF [50].
In the 1970s and 1980s, various groups have started to introduce effective HamiltoniansĤ eff with the idea of using them in Hartree-Fock (HF), i.e. mean-field calculations. In other words, the idea was to calculate the nuclear ground-state by minimising the energy where jΦi is the most general Slater determinant consistent with the symmetries displayed by the system and H eff ¼ T þ V eff is a properly designed effective Hamiltonian. In principle, an effective interaction V eff that accounts for the short-range, mainly attractive character of the NN interaction in the medium could work as the saturation point can result, in this case, from the balance with the repulsive kinetic energy. In practice, the situation is more involved. In particular, it has been soon realised that only density-dependent interactions can provide realistic results for the properties of atomic nuclei. This emerged from the experience of several decades, as reviewed in [50]. The most widely used effective interactions Vðr;r 0 Þ have been those of the Skyrme type [51][52][53], and those of the Gogny type [54,55]. These interactions are strongly spin-and isospin-dependent. In terms of their radial dependence, the terms of the Gogny forces are Gaussians, while Skyrme forces are proportional to Dirac δ-functions and mimic the finite range with a dependence on the relative momentum between the nucleons. Both kinds of forces include % 10-15 terms, so this is the numbers of parameters to be fit against a set of data. Both include at least one densitydependent term. Is this term avoidable or not? So far, no density-independent Hamiltonian has the same performances of a density-dependent one. Pedagogical examples along this line can be found in [56], and a simple yet powerful argument why naïve mean-field with a density-independent interaction cannot work for nuclei is provided in [57]. A thorough analysis has been carried out in Ref [58].
A density-dependent Hamiltonian H eff ½ρ can be thought as an object that poses conceptual problems, unless one merely considers it as a tool to generate an EDF through Equation (9) [59]. Starting from the 1990s, a philosophical and practical shift has taken place: existing EDFs based on effective Hamiltonians have been conceived to be per se a realisation of nuclear DFT, and practitioners have started to generalise their form without resorting to underlying Hamiltonians. In the pioneering work of Ref [60], a spin-orbit term has been introduced in the EDF without any reference to a spin-orbit potential. This idea has been since then become the standard, and we presently should talk of a nuclear EDF without any reference to a generating Hamiltonian. [2] In electronic systems, the fixed ion positions constrain the position of centre-ofmass of the whole system in the laboratory frame. In nuclei, that are self-bound systems, this is not the case. It has been argued by several authors [61] that the HK theorem in its standard formulation is irrelevant to the nuclear case because it concerns the laboratory density: as the centre-of-mass can translate, this density is almost zero everywhere. Experiments probe actually the intrinsic density (relative to the nuclear centre-of-mass). The solution of this puzzle can be found, based on the fact that it has been proven [62] that, given a Hermitian operatorQ, one can build an energy functional depending on QðrÞ ; hQðrÞi that is universal in the HK sense and has its minimum at the correct value of Q with the correct energy. Thus, one can replace the laboratory density with the intrinsic density in the HK theorem [61,63,64]. [3] Symmetry breaking is one of the essential elements of nuclear DFT.
The intrinsic nuclear density can break some symmetries that characterise the nuclear Hamiltonian. The typical example is that of deformed nuclei, that are so defined since the intrinsic density is non-spherical and a clear signature is the existence of rotational bands. However, also nuclear superfluidity can be cast into the framework of symmetry breaking [10]. These aspects have been already mentioned in Sec. 1.
It is a common practice to refer to these cases as examples of spontaneous symmetry breaking (SSB). We shall not dwell on the differences between SSB breaking in either infinite or finite systems, as this is beyond the scope of the paper. We just remind that in finite systems, quantum fluctuations allow coupling symmetry-breaking states with one another and restore the symmetry [25,65]. Indeed, in the case of the rotational symmetry, the nucleus displays a zero-point rotational motion so that the shape in the laboratory system is spherical (i.e. consistent with the symmetry).
The formal way to start from a symmetry-breaking solution and restore the symmetry is to use projection methods. The operators which are introduced project on the relevant quantum numbers, like angular momentum in the case of the rotational symmetry. Superfluid nuclei are better described by approaches, like HF-Bardeen-Cooper-Schrieffer (HF-BCS) or HF-Bogoliubov (HFB), in which the particle number symmetry is broken and can eventually be restored through projection. An excellent review on projection methods can be found in Ref [66]. The aspect of symmetry breaking and restoration as a unique feature of nuclear DFT is emphasised in the lecture notes of Ref [67] (see, nonetheless, Ref [68]).
The issue of symmetry breaking and restoration in finite systems cannot yet be said to be fully understood. In Ref [69], the formalism of decoherence histories is applied to study it. Within such framework, correlations between the particles at different times become a key issue, and one can appreciate the fact that symmetry restoration is more important in nuclei than in molecules because in nuclei the adiabatic approximation is less good. In simpler terms, the rotational motion of a deformed system (that is responsible for the symmetry restoration) is slower than the intrinsic motion 8 but not extremely slower. The difference is, as a rule, larger in electronic systems. [4] The structure of the nuclear EDF is richer, in terms of relevant operators, compared to the electronic case.
Let us consider the Skyrme Hamiltonians that we have previously introduced. Because of their zero-range character, Skyrme interactions when inserted in Equation (9) will always produce a local functional, that is, an energy density that depends on densities at the very same point. If we now wish to use it as a starting point and build a more general and flexible ansatz, we ought to remember that the most general local functional will depend on all possible local densities.
All local densities can be built in a systematic way, and the interested reader can follow Refs [70,71]. and/or Sec. 7.3.4 of Ref [67]. Local densities are obtained by letting operators acting on the non-local density ρðr;r 0 Þ, and equatingr ¼r 0 afterwards. Due to the strong spin-dependence of nuclear properties, we have to distinguish the non-local density from the non-local spin density. The orbitals ϕ j ðrÞ that we have introduced around Equation (7) are actually spinors ϕ j ðr; sÞ. Then, we are forced to introduce ρðr;r 0 ; s; s 0 Þ ; X j ϕ Ã j ðr; sÞϕ j ðr 0 ; s 0 Þ; and this latter quantity can be written as ρðr;r 0 ; s; s 0 Þ ¼ ρðr;r 0 Þ1 þsðr;r 0 Þ Áσ; where the Pauli matricesσ form a complete basis, together with the identity matrix 1, in order to express the 2 Â 2 matrix defined by the spin indices s and s 0 . Derivative operators can be applied in a systematic way and, up to second order, the procedure is described in Ref [71]. Since the functional must be invariant under space translations, only derivatives with respect to the relative coordinate,Ñ ÀÑ 0 , must be considered. Due to the identitỹ the action of the second derivative boils down to that of the operatorÑ ÁÑ 0 . Consequently, starting from ρðr;r 0 Þ one obtains the local current and the kinetic energy density τðrÞ ¼Ñ ÁÑ 0 ρðr;r 0 Þ h ĩ r¼r 0 : There are the corresponding derivatives of the spin densities. Moreover, one must consider both neutron and proton densities, 9 so that the number of densities is doubled when the complete set obtained at a given order 10 is written.
Once all possible local densities are defined, the most general EDF will be a combination of scalar terms that depend on these densities, and are invariant under parity and time-reversal. Even if one limits oneself to terms that are quadratic in the densities, the number of these terms can grow up dramatically if one allows higher-order derivatives (cf. Ref [72], and in particular Table XXVII therein). In front of each term one may write either a number, or a scalar function of ρ. One notices already here the lack of a perturbative parameter, or of a strategy for systematic improvement. On top of it, if the structure of the EDF becomes involved, fitting all parameters becomes practically impossible. Last but not least, the question arises whether correlations exist among different terms of the EDF. Note that we have left aside the question about the non-locality, by assuming that gradient terms take care of it. Ref [73] is useful to compare with the electronic case. [5] Is a covariant EDF needed or is a nonrelativistic formulation appropriate?
Relativistic effects could be considered in the kinematics (as we stated in Sec. 1, typical kinetic energies of nucleons in the nucleus are % 35 MeV which amounts to v=c % 0.28); however, the success of nonrelativistic theory means that these effects can be reabsorbed when defining the parameters of the nuclear EDF. A relativistic, or covariant, theory has quite different advantages, namely it permits addressing the question about the spin-orbit terms as well as other questions on the structure of functionals. Excellent review papers on covariant nuclear EDFs can be found in Refs [74,75].
Historically, these EDFs have also evolved in a similar manner as we discussed above. They started as energy functionals generated through effective Lagrangians that include nucleons, effective mesons, and their mutual couplings. Recently, they are also built without any reference to an underlying Lagrangian. One builds instead all possible densities in terms of the Dirac spinors ψ i associated with the single-nucleon orbitals (as we did with the Schrödinger orbitals ϕ i above). The results are the following 16 covariant local quantities: Γ ¼ 1 and Γ ¼ γ μ provide the familiar scalar and vector densities. The most general covariant EDF will be a combination of all possible terms that are built with the covariant local densities and respect symmetries. Although a number of efforts have been devoted to studying the nonrelativistic reduction of covariant EDFs, 11 there is not a systematic way to connect existing covariant and nonrelativistic EDFs. Despite this, if one looks at results obtained along the two or three last decades, it is hard to detect any significant overall difference between the performances of nonrelativistic and covariant functionals. One can only stress again that the latter contain in general less parameters, and that the spin-orbit phenomenology is fixed by the Dirac structure. However, a lot of physics has been learnt through the systematic comparison of the different types of EDFs. [5] There are intrinsic limitations in EDFs, related to the choice of nucleon densities as the only degrees-of-freedom even at low-or high-density. What about formation of other particles at high densities, and clusterisation at low densities?
We have already mentioned, along Sec. 1, that we expect hyperons to form in neutron stars simply because the energy density above 3 or 4 times saturation density will allow so. In general, it is unclear which is the upper limit in density for a purely nuclear EDF, and how to match with a similarly successful model that includes other particles that are created above a given threshold.
In a different sense, and yet with a kind of analogy, there are limitations for EDFs at low densities. Recently, clusterisation phenomena have captured interest in nuclear physics (see, e.g. the review paper [76]). In light nuclei, there are several evidences of phenomena in which nuclei clusterise by forming (mainly although not exclusively) α-particle clusters. Intuitively, one can grasp the fact that by forming α-clusters nuclear matter can gain energy at some low density, because of the large binding energy of the αparticle. Deuteron clustering may also take place in keeping with the relatively strong proton-neutron attraction.
Clustering, once again, lie between nuclear physics and astrophysics as neutron stars also explore low densities. For astrophysical simulations, some multi-purpose EDFs have been designed that include some aspects related to clustering (see [77,78] and references therein). Despite this, further investigations are needed to clarify which clusterisation phenomena are already implicitly included in existing EDFs [79], and which phenomena call for an extension of the degrees of freedom that one usually puts in EDFs so far. This is another difference with electronic DFT, where no purely electronic bound states show up.

Applications
The format of this paper prevents displaying an exhaustive list of DFT applications; in this Section, we simply focus on few illustrative examples.

Ground-state properties
One of the main observables that one aims at calculating within DFT is the total energy, that is, the binding energy (1). Typical errors that affect the existing functionals, in comparison with the experimental data, are of the order of % 1-2 MeV. To achieve a better accuracy, it is mandatory to add extra terms that are outside the pure DFT philosophy and have a mere empirical justification (cf. Table I of Ref [80] and the related discussion). These accuracies may be considered as satisfactory, with respect to the total binding energy of % 10 2 À 10 3 MeV. However, one should bear in mind that predicting masses is instrumental to obtain the Q-values of many reaction processes that depend on mass differences, 12 and these can either benefit from error cancellations or be harmed by their amplification. The transition rates can depend on powers of the Q-value, like in case of β-decay.
Therefore, the effort of pushing down the error of mass models is still continuing.
In Figure 4 we show two examples of comparison between the result of EDF calculations and the experimental data. The Skyrme functional SLy4 [53] is the latest of the standard Skyrme functionals built from a Hamiltonian. Basically, only data from magic nuclei are considered in the parameter fit. A more sophisticated fitting protocol has been used to fit UNEDF0 [81], and both spherical and deformed nuclei are included in that protocol. This helps in improving the overall trend of the residuals. However, these are not really normally distributed as one notices immediately the arch-like behaviour in both panels. This is common to many EDFs. The arches span shells: in other words, it seems that, still, reproducing with similar accuracy both closed-shell and open-shell nuclei represents a non-trivial challenge.
The prediction of masses is closely related to the prediction of the proton and neutron drip lines. In Figure 1 we have shown the drip lines obtained with some EDFs. A thorough analysis of the uncertainties that plague currently our knowledge of the drip lines can be found in Refs [13,14]. Experimentally, the neutron drip line is known only up to Z ¼ 8. Going towards heavy nuclei, the theoretical predictions for the neutron drip line become increasingly model-dependent: the systematic uncertainties can reach % ten mass units. Proton drip lines suffer from less uncertainties because they lie closer to the stability valley.
Other quantities that testify to the quality of a given EDF are the neutron and proton density distributions (cf. Sec. 5.2 of [80]). The most relevant moments of these distributions are the neutron and proton r.m.s. radii, namely the square roots of Unfortunately, we do not know experimentally very well the neutron densities so far, due to the lack of a model-independent probe. Charge density distributions have been in many cases obtained thanks to elastic electron scattering experiments. Laser spectroscopy also allows extracting from the atomic transitions some moments of the nuclear charge distributions. There exist precise relationships between the charge density r.m.s. radii and the proton density r.m.s. radii, that involve the nucleon electromagnetic form factors; nevertheless, the approximate relationship hr 2 charge i ¼ hr 2 p i þ 0:8 2 ½fm 2 has been shown, as a rule, to be accurate enough.
Some examples of comparison of the results of DFT calculations with experiment, in the case of charge radii, are shown in Figure 5. Starting our short discussion from the heavier systems, the two different EDFs have a very small spread, among themselves and with experiment, in the case of the Sn isotopes. In the case of the Pb isotopes, if we restrict ourselves to neutron numbers N larger than 116, HFB24 is closer to the experimental data in terms of absolute values but DD-PC1 predicts better the kink after N = 128, which has been known for some time to be attributable to the different spin-orbit terms associated with the functionals. 13 In the case of Ni isotopes, DD-PC1 is closer to the data whereas, finally, if we come to the case of Ca, the bell-shaped trend of the radii between 40 Ca and 48 Ca is not reproduced by essentially any of the available EDFs. Albeit short, this discussion highlights how to charge radii can really inform us about several features of EDFs. Figure 5. Comparison of charge radii from experimental measurements and from DFT calculations. The experimental data from Ref [82] are represented by symbols. The lines refer to different calculations, performed either with the nonrelativistic HFB24 mass model [83] or with the covariant functional DD-PC1 introduced in Ref [84]. See the text for a short discussion. Figure taken from Ref [87].
There are many other open questions in nuclear science, to which the attention of DFT practitioners has turned. In addition to neutron and proton drip lines, one of the hot topics is related to superheavy nuclei. Intensive research has been devoted to identifying these elements, and the periodic table now extends far beyond Z ¼ 92 [85]. The existence of these high-Z elements cannot be understood by means of liquid drop-like formulas, as the Coulomb repulsion would be too strong. No other microscopic approach apart from DFT can be employed because of the size of the systems that makes computation quite demanding. The role of EDFs for present and future challenges in this domain is discussed in Ref [86].

Collective excitations
Nuclear DFT can be extended, beyond the static regime, to the study of nuclear excitations. It has been known for some time that nuclei respond in quite different ways to external fields, the response at high frequencies (like in the giant resonance region) being more elastic and that at low frequencies (associated to energies of the order of % MeV) being more plastic. Pedagogical introductions to the DFT extension to the time-dependent frame, in the elastic regime (i.e. the small-amplitude limit) can be found in [56,87]. A broad survey of all relevant, recent results are provided by Ref [80]. In general, many works have highlighted that the excitation modes of nuclei are sensitive to features of the EDFs that do not manifest themselves in the nuclear ground states.
To give a simple example, the nuclear so-called 'breathing' mode, 14 in which the nucleus expands and contracts, gives information on the nuclear compression modulus. Although not in a straightforward manner, the compression moduli of finite nuclei can be extrapolated to the compression modulus of infinite, uniform nuclear matter. The compressibility is usually defined as where P and V are, respectively, pressure and volume, and the second equality holds for a fixed number of particles A, when the density ρ ¼ A=V is the variable. From this equation, it becomes apparent that if we know the compressibility, we have a further constraint related to the dependence of energy on the density ρ. In nuclear physics, it is customary to define at the saturation point of nuclear matter the so-called incompressibility given by K 1 ; 9 ρ 0 χ À1 . Recent findings confirm that nuclear matter is more incompressible than steel by almost 22 orders of magnitude [88]. In general, other kinds of responses provide valuable information on the Equation of State (EoS) of infinite matter. Equation (15) sets a constraint on the dependence of the energy per particle on the total density ρ. However, nuclear matter is composed of neutrons and protons. One can use as variables either ρ n and ρ p , or the total density ρ and the relative neutron-proton asymmetry β ; ρ n À ρ p =ρ. Starting from the energy per particle in symmetric nuclear matter, that is, E A ρ; β ¼ 0 ð Þ, one can perform a Taylor expansion in β. Odd terms are prohibited by isospin symmetry, 15 and one can write Many calculations have actually shown that the quartic term in β 4 does not play a significant role around saturation density and, to some extent, above it. The quantity SðρÞ is called symmetry energy. Neglecting the Oðβ 4 Þ term in Equation (16), the symmetry energy is also the difference between the energy per particle of neutron matter and symmetric matter. Symmetric nuclear matter is the most stable, as neutrons and protons can occupy the lowest possible quantum levels, while increasing the number of neutrons implies promoting particles to highest levels due to the Pauli principle; in addition, one misses attractive contributions from the strong proton-neutron interaction. This is why we expect the symmetry energy to be positive. A schematic picture is shown in Figure 6.  (16) and in the text below. The curves labelled by SNM (PNM) correspond to the energy per particle in symmetric nuclear matter (pure neutron matter), namely to β ¼ 0 (β ¼ 1), The symmetry energy SðρÞ is the difference between the two curves.
The collective modes in which neutrons and protons oscillate out of phase, so that the density stays approximately constant but β changes, inform us about δ 2 E δβ 2 which is the symmetry energy [80]. One example is the so-called Giant Dipole Resonance (GDR), in which all protons vibrate against all neutrons under the effect of a uniform electric field, like in the plasmon excitations of electronic systems. This kind of excitations is called isovector as the local change of protons into neutron implies a local change of a unit of isospin. Charge-exchange excitations in which a proton turns into a neutron, or vice-versa, are also very interesting isovector excitation modes of nuclei and they can also spontaneously occur under the form of β-decay [89,90].
While most of these modes represent examples of small-amplitude motion, applications of nuclear DFT to the large-amplitude case are considered in [57]. Reactions between two colliding nuclei [91] or fission processes [92] can be and have been calculated using EDFs.

Neutron star calculations
Neutron stars, that we have introduced in Sec. 1, have a density profile that spans 2-3 orders of magnitude. We have also reminded in Sec. 3 that EDFs need still to be improved, in order to deal with the changes that nuclear matter undergo at low and high density. Nevertheless, a large number of DFT-based calculations have been devoted to neutron star properties.
The symmetry energy, that we have defined in Equation (16) as the energy per particle to change protons into neutrons starting from symmetric matter and moving to neutron matter, affects almost every property of neutron stars either directly or indirectly. The very existence of a neutron star, and its bulk properties like the mass and the radius, depend on it: if the symmetry energy grows enough as a function of the density, more and more gravitational energy can be counterbalanced and the neutron stars become more massive. For quite some time, it was believed that neutron stars should have masses around 1.4 M (M being the symbol for the solar mass). Along the last decade, new observational techniques have become possible, and stars with mass larger than 2 M have been identified; moreover, the determination of radii has been constantly improved and is still improving (see the references mentioned in [23]). Masses and radii put firm constraints on the EoS of neutron matter. If hyperons, let alone deconfined quarks, were present in the inner part, the conclusions concerning the EoS would change; but currently, there are neither positive nor negative evidences related to such presence.
The symmetry energy also affects the composition of the outer crust (where less neutron-rich nuclei are present if the symmetry energy increases) [93,94], the transition density to the inner crust (the density at which neutrons spill out of nuclei) [95], and the proton fraction in the inner crust and outer core [94]. Last but not least, in the inner crust, neutrons are superfluid and the calculations must account for this exotic form of superfluidity [96]. The quantities that we have mentioned in this paragraph are not directly observable, at variance with masses and radii, and yet affect in a significant manner the neutron star behaviour.
Neutron stars are formed when a star with M > 8M ends its life, at a moment when only nuclei in the Fe region remain: these have maximal values of BEðN; ZÞ [Equation (1)], and no nuclear reaction can release energy and compensate the gravitational pressure. Then, collapse takes place and is followed by a supernova explosion and the formation of a proto-neutron star. The simulation of this process requires the nuclear EoS [97], including possibly some clustering aspects at low density as we mentioned briefly in Sec. 3. The processes explore high temperatures and this must be accounted for. While a number of nuclear EoS are already available in tabulated form to be used for such simulations, the analysis of their accuracy, of the sensitivity of the results to specific inputs, and of the possible missing physics, is still ongoing.
Many types of electromagnetic signals can be detected from neutron stars. It is useful to remind that neutron stars were identified under the form of pulsars, namely stars that rotate along their magnetic axis by emitting dipole radiation. There are also signals that testify to the existence of cataclysmic events, like the so-called 'giant flares'. Oscillations in the X-ray flux that follows these flares are believed to be associated to the fact that neutron stars can undergo starquakes: these vibrational modes have been also used as a signature for the nuclear EoS (see [98] and references therein).
Last but most importantly, we have just seen the advent of an exciting time related to multi-messenger astrophysics. An event associated with the merging of neutron stars, GW170817, has led for the first time to the detection of gravitational waves, of the whole electromagnetic signal from gamma-rays to radio waves, and of neutrinos [99][100][101][102]. This has triggered a number of analysis, part of which still ongoing, related to various crucial physics questions like, e.g. the nucleosynthesis process [103]. As far as the EoS is concerned, numerical General Relativity (GR) simulations performed by the LIGO and Virgo (LV) Collaboration have been able to set constraints on the so-called 'chirp' mass associated with the event, that is Þ 1=5 (M i are the masses of the two neutron stars), as well as on the tidal polarizability. This latter quantity reads where G is the gravitational force constant and k 2 the Love number [104,105], and it is associated with the ratio of the induced quadrupole mass moment and the inducing gravitational field. It is very sensitive to the mass and radius of the compact object, as they enter with a large power, and these, in turn, are strongly related to the EoS. In future, new events of the same or other types (like the merging of a neutron star and a black hole) are expected to give a strong boost to our understanding of neutron matter at high density, perhaps changing the current paradigm.

Perspectives
It is hard to survey all possible perspectives in the exciting, and rapidly developing, domain of nuclear DFT. In particular, these perspectives are expected to be shaped by (i) discoveries of new nuclei and new nuclear properties, close to the drip lines and in the superheavy region; (ii) new observations of neutron star properties and/or merging events with the associated multi-messenger signal; (iii) progress in ab initio nuclear theory. Not so much has been done, so far, to ground nuclear DFT on underlying theories, in contrast with the case of electronic DFT. The systematic construction of exchange-correlation functionals based on the state-of-the-art calculations of the electron gas does not have a real counterpart in nuclear physics, despite some attempts to start from Brueckner-Hartree-Fock calculations and build the volume part of an EDF [106], or a Skyrme interaction [107,108]. The strategy proposed in Ref [59] has still to be developed, and also in the case of EDFs built around the idea of symmetry breaking and restoration [109] there are still steps to be undertaken in order to include all correlations from many-body perturbation theory, overcoming the technical difficulties. The most general program for an ab initio-based DFT, that has been outlined in Ref [110], based on the idea of the Legendre transformation of the effective action, is even more in its infancy. Yet, it might be the most consistent path to cast DFT into the EDF scheme [111].
Despite these difficulties one can expect, generally speaking, important steps forward in the next future as far as merging ab initio and DFT are concerned. DFT is a mature field and ab initio techniques are, on the one hand, progressing in a significant fashion and need, on the other hand, to be recast into something less computationally demanding if heavy systems and/or highly excited states are targeted. A different avenue, that may lead to EDFs rooted either in experimental densities or in densities calculated by means of an underlying ab initio theory, consists in exploiting reverse engineering. Basically, one formulates a sort of Kohn-Sham inverse problem, so to infer the Kohn-Sham potential from the densities and, in a further step, the energy functional from the Kohn-Sham potential. Very preliminary steps have been taken along this line recently [112]. Machine learning techniques have been introduced as a tool to attack the manybody problem [113], so one can expect progress in this respect.
It is likely that nuclear structure physics, and DFT-related works in particular, will proceed in several directions rather than in a single one. I believe that hybrid approaches, that combine robust statistical and fitting techniques starting from the ab initio findings with appropriate physics insight, have the best chance to give a new boost to the field. One can optimistically conclude that we expect to be able to write a more conclusive review paper soon, with less open questions, along the lines of the predictions made in Ref [114].