Microscopic-Macroscopic Approach for Binding Energies with the Wigner-Kirkwood Method - II

The binding energies of deformed even-even nuclei have been analysed within the framework of a recently proposed microscopic-macroscopic model. We have used the semiclassical Wigner - Kirkwood $\hbar$ expansion up to fourth - order, instead of the usual Strutinsky averaging scheme, to compute the shells corrections in a deformed Woods - Saxon potential including the spin-orbit contribution. For a large set of 561 even-even nuclei with $Z\ge 8$ and $N\ge 8$, we find an {\it rms} deviation from the experiment of 610 keV in binding energies, comparable to the one found for the same set of nuclei using the FRDM of M\"oller and Nix (656 keV). As applications of our model, we explore its predictive power near the proton and neutron drip lines as well as in the superheavy mass region. Next, we systematically explore the fourth - order Wigner - Kirkwood corrections to the smooth part of the energy. It is found that the ratio of the fourth - order to the second - order corrections behaves in a very regular manner as a function of the asymmetry parameter $I=(N-Z)/A$. This allows to absorb the fourth - order corrections into the second - order contributions to the binding energy, which enables to simplify and speed up the calculation of deformed nuclei.


I. INTRODUCTION
The models of nuclear masses are continuously challenged by the advances in experimental techniques which nowadays are extending the nuclear chart to previously unexplored regions of exotic isotopes and superheavy elements. The theoretical description of nuclear masses takes place primarily along two main approaches. On the one hand, in the microscopic nuclear models, the nuclear binding energy is obtained from calculations with energy density functionals based on effective nuclear interactions [1][2][3]. In the microscopic-macroscopic (mic-mac) models [2, 4,5], the nuclear binding energy is obtained as the sum of a part that varies smoothly with the number of nucleons plus an oscillatory correction originated by the quantum effects. The smooth part of the mic-mac models is obtained from a liquiddrop model approach, whereas the shell correction is usually evaluated by the Strutinsky averaging method in an external potential well.
In our previous works [6,7], we have demonstrated that the Strutinsky average can be replaced by the semiclassical energy computed by means of the Wigner -Kirkwood (WK) expansion of the one -body partition function [8][9][10][11][12][13][14][15], in order to evaluate the shell corrections of a system of N neutrons and Z protons at zero temperature in an external potential. There are some reasons supporting this choice as we have discussed in Ref. [6].
On the one hand, it has been shown that the Strutinsky level density is an approximation to the WK level density [16]. On the other hand, the WK -expansion of the density matrix has a variational content and it is possible to establish a variational theory based on a strict -expansion [15,17]. We shall point out that the WK expansion is also well suited to deal with nuclei close the drip lines. Although the WK level density exhibits a well known ε −1/2 divergence as ε → 0 for a potential that vanishes at large distances, integrated moments of the level density, such as the energy and the accumulated level density, are well behaved in the ε → 0 limit as it has been demonstrated in Ref. [15]. It has been shown that these shell corrections, along with a simple six parameter liquid drop formula, yield a good description of ground -state masses of spherical nuclei spanning the entire periodic table [6]. The model has also been applied to calculate the binding energies of few deformed nuclei, with a good degree of success [7]. In the present work, we extend the work reported earlier [6] to the deformed nuclei and explore the predictions of the model in exotic scenarios such as drip line nuclei and the superheavy region. In this work, we mainly restrict our attention to the 3 even -even nuclei.
One of the important conclusions of Ref. [6] is that in this model it is necessary to carry out the WK expansion up to the fourth -order in to obtain accurate shell corrections, which implies that in this case one needs to work out derivatives of the single particle potentials (nuclear potential, Coulomb potential as well as the spin -orbit potential) up to the fourth -order, which is a rather cumbersome task. Therefore, this gives rise to an interesting and important question: can the effects of the fourth -order corrections to the binding energy be absorbed into the second -order ones? This question is important from theoretical as well as practical point of view. Theoretically, this would imply that the WK series has been partially re-summed, whereas from a practical point of view, it implies that it is sufficient to expand the one -body partition function up to second -order in to obtain shell corrections with comparable accuracy.
The absorption, if possible, would imply that there is a factor (we denote the factor by α), which may be a function of mass number, charge number, neutron number or combinations thereof, defined as such that where, E ( 2 ) and E ( 4 ), respectively, are second and fourth -order WK corrections to energy. This is an important issue discussed in the present article.
We summarise the essential details of the semiclassical Wigner -Kirkwood expansion of the one -body partition function in the second section. The detailed results and their analysis forms the subject matter of the third section. The parameters of the macroscopic part of our mic-mac model, which also includes curvature correction [5] and the Wigner term [5], have been obtained by minimizing the χ 2 value of the energies using a selected set of 561 even-even deformed and spherical nuclei. The ability of this mic-mac model to describe nuclei in the exotic scenarios is explored in section 4. On the one hand, masses of very proton rich nuclei, measured recently [18], are compared with the predictions of our model. On the other hand, the upper limit of the outer crust in neutron stars is studied, which involves nuclei near the neutron drip line. Finally, we explore the superheavy region, 4 and compare the theoretical alpha decay Q values and the corresponding half lives with the experimental values [19]. The systematic investigation of the absorption factor α as defined above is contained in the fifth section. The summary and conclusions are given in the last section.

II. FORMULATION
For a system of N non interacting Fermions at zero temperature in a given external potential, the quantal one -body partition function is given by: The Hamiltonian of the system (Ĥ) is expressed as: with V ( r) being the one-body central potential andV LS ( r) the spin-orbit interaction. The replacement of the Hamiltonian in the above equations by the corresponding classical Hamiltonian leads to the well -known Thomas -Fermi equations for particle number and total energy. The Wigner -Kirkwood semiclassical expansion amounts to expansion of the quantal one -body partition function in the powers of Planck's constant, , yielding systematic corrections to the Thomas -Fermi energy and particle number [8][9][10][11][12][13].
As stated before, in this work, we carry out the WK expansion up to the fourth -order in . With the spin -orbit interaction, the WK expansion of the partition function can be written schematically as: where, Z (4) (β) (Z SO (β)) is the WK partition function for the central potential (spin -orbit part). The explicit expressions for these partition functions can be found in [6,10].
The level density g W K , the particle number N and the energy E W K are obtained by appropriate Laplace inversion of the WK partition function, as follows: In these expressions, V is the mean field, f is the spin -orbit form factor, κ is the strength of spin -orbit interaction, and λ is the chemical potential.
The shell corrections, which are the difference between the quantum mechanical and the corresponding averaged energies, can now be obtained by subtracting E W K from the quantum mechanical energy. For our calculations we choose a Woods-Saxon potential as 6 mean field and a suitable Woods -Saxon form factor in the spin -orbit sector. These potentials are generalised for taking into account deformation effects and their corresponding parameters are given in Ref. [6]. The Coulomb potential has been obtained by folding the proton density distribution with the Coulomb interaction [6]. In the microscopic part we have also included pairing correlations using the Lipkin -Nogami scheme [20][21][22], as described in details in Ref. [6].

III. CALCULATION OF BINDING ENERGIES
In the present work, we generalise the liquid drop formula employed in [6] by adding a deformation dependent curvature energy term and the Wigner term. The curvature energy term is found to be important in improving the agreement achieved between calculations and the corresponding experimental binding energies [5]. The Wigner term is expected to be important for light nuclei as well as to describe nuclei close to the proton drip line.
Therefore, the modified liquid drop formula used in this work reads: where the terms respectively represent: volume energy, surface energy, curvature energy, Coulomb energy, correction to Coulomb energy due to surface diffuseness of charge distribution and the Wigner energy. The coefficients a v , a s , a cur , k v , k s , k cur , r 0 and C 4 are free parameters; T z is the third component of isospin, and e is the electronic charge.
Several parametrisations of the Wigner term are available in the literature (see, for example, [2, 5,23]). Here, we adopt the following ansatz for the Wigner term with a cut off on charge and mass numbers: where, w 1 and w 2 are free parameters. The cut offs on charge and mass numbers have been introduced since it is expected that the Wigner term will make significant contributions for nuclei with low masses.
The Coulomb, surface and curvature terms appearing in the liquid drop formula, as defined above in Eq.(15), need to be modified for the deformed shapes. In particular, the where, the symbols have their usual meanings. Notice that the integrals have been carried out over nuclear volume, and the lengths have been measured in units of the radius parameter R o of the nucleus with zero deformation. The transformation from six dimensional to four dimensional integrals has been accomplished by following the technique developed by Kurmanov et al. [24]. The surface term, on the other hand, is simply modified by the ratio of deformed to the corresponding spherical surface areas. The curvature energy term, too, needs to be modified to take the deformation effects into account. The modified curvature energy (E cur ) reads: where, E 0 cur is curvature energy at zero deformation; R 1 and R 2 are the principal radii of curvature of the nuclear surface (in the units of R o ), defined by r = r s ; and dS refers to the area element of the nuclear surface. The surface parametrisation assumed in the present work is given by: Here, the Y λ,µ functions are the usual spherical harmonics and the constant C is the volume conservation factor (the volume enclosed by the deformed surface should be equal to the volume enclosed by an equivalent spherical surface of radius R 0 ): The term Z 2 /A, which is the correction to Coulomb energy due to surface diffuseness of the charge distribution, does not have any explicit deformation dependence. This is because the distance function chosen here is such that the surface thickness is the same in all the directions (see discussion about this in Ref. [6]).
The total binding energy of a nucleus with N neutrons, Z protons and deformation parameters β 2 , β 4 and γ is given by: 8 where, δE represents the microscopic part of the binding energy (shell correction plus pairing energy). The microscopic part has been multiplied by a factor η, which is chosen to be 0.85.
One of the reasons for introducing such a factor is that the Coulomb potential used in the present work is less repulsive near r = 0 than the corresponding value obtained by using the hard sphere approximation, used in the fit of proton mean field (see discussion on this point in Ref. [6]).
The free parameters of the liquid drop formula are determined by minimising the χ 2 value in comparison with the experimental binding energies [25]: where E(N j , Z j ) is the calculated total binding energy for the given nucleus, E (j) expt is the corresponding experimental value [25], and ∆E as open shell nuclei, many of which are expected to be deformed. The main task now is to determine the liquid drop parameters as well as the optimal deformation parameters. The calculation proceeds in the following steps: 1. Assuming the previously reported [6] values of the liquid drop parameters, the binding energies of these nuclei are obtained by minimising on a range of β 2 values (β 4 is set to zero in this step). This gives a preliminary estimation of β 2 . Next, keeping this β 2 fixed, β 4 is varied to obtain minimum energy. Thus, we now have preliminary values of both the deformation parameters.
2. In the next step, keeping the deformation parameters fixed as obtained in the earlier step, the liquid drop parameters are fitted by minimising χ 2 .
3. With the new values of liquid drop parameters, the deformation parameters are obtained once again as described in step 1, followed by a final re-fit to the liquid drop parameters.

9
The numerical values of the new constants of the liquid drop formula obtained through this minimisation procedure are: a v = -15.435 MeV, a s = 16.673 MeV, a cur = 3.161 MeV, k v = -1.874, k S = -2.430, k cur = 0 (see discussion below), r 0 = 1.219 fm, C 4 = 0.963 MeV, w 1 = -2.762 MeV and w 2 = 3.725. The values of volume, surface and Coulomb coefficients differ from those reported earlier [6], primarily due to the inclusion of curvature and Wigner terms and the deformation effects. The curvature term, as described earlier, depends on the mean curvature of the nucleus, which is a function of the geometry of the nuclear surface.
Therefore, the curvature energy, a priory, is expected to modify the surface energy term as well as the Z 2 /A term, which is the correction due to the surface diffuseness of the charge density term. The somewhat smaller value of the volume coefficient reported here, is not surprising. The reduction is due to the influence of the curvature term, as it has also been found by Pomorski and Dudek (see Table 1 of Ref. [5]).
It is to be noted that the coefficient of the isospin dependent term in the curvature energy is very difficult to determine with experimental masses. In our case the resulting statistical error in the corresponding parameter turns out to be more than 50% of the numerical value of the coefficient. Further, this term is found to weaken the strength of the isospin dependent term in the surface energy by a factor of 5. The isospin dependence in the curvature term, therefore, has been dropped from the present investigation.
The rms deviation of the calculated binding energies with respect to the experiment obtained is 610 keV. The Möller -Nix calculations [27], for the same set of nuclei, yield a deviation of 656 keV. The explicit values of binding energies of our selected set of 561 even-even nuclei used in the minimisation procedure can be found at [26]. The present calculation establishes that our model is indeed capable of reproducing binding energies of deformed nuclei as well, with excellent accuracy. The difference between the calculated and the corresponding evaluated [25] binding energies is presented in Fig. 1 and that the shell -gaps are also reproduced nicely.
In addition, the systematics of deformation parameters obtained in these calculations turns out to be reasonable. As an illustrative example, we focus on the Sr -Zr region. It is well known from the systematics of experimentally measured charge radii [28] that the charge radii increase dramatically by 2% for 97 Rb, 98 Sr and 100 Zr, in comparison to their respective lighter isotopes. This jump may be attributed to the possibility of onset of highly deformed shapes in the ground -state, around this neutron number (see, for example, [29]).
Our calculations, too, reveal existence of highly deformed ground -states (with β 2 ∼ 0. Further, it is also well known that the ground -states of 72 Kr, 76 Sr and 80 Zr have very large (∼ 0.4) deformation. This is known to be due to population in the intruder 1g 9/2 state. Thus, the ground -state of 80 Zr is a 12 particle -12 hole state, which is manifested again by an extremely large stable deformation in the ground -state of 80 Zr. This has been verified independently, for example, by the relativistic mean field calculation [30], density dependent Hartree Fock calculation with Skyrme interaction [31], as well as by the Hartree Fock band mixing calculation [32]. The deformation parameters reported in the Möller -Nix table [23], too, are consistent with the discussion above. It is gratifying to note that the The masses of 63 Ge, 65 As, 67 Se and 71 Kr have recently been measured [18]. These nuclei are very proton rich, and are expected to be close to the drip -line. Notice that these nuclei are odd -even and even -odd. In this preliminary test of our model near proton drip line, we use the simple uniform filling approach for the calculation of the pairing energy.
The calculated binding energies and one proton separation energies (S p ) for these nuclei, along with the corresponding experimental values [18] and those reported by Möller and Nix [27] are presented in Table 1. The binding energies as well as S p values obtained in the present work are found to be quite close to the experiment. This indicates that the present model extrapolates reliably up to the proton drip lines. The nucleus 65 As is reported to be slightly unbound against proton emission with S p = −90 ± 85 keV [18]. Our calculation, on the other hand, yields a positive value of S p for 65 As, indicating a proton bound nucleus.
However, it should be noted that the separation energies are obtained by taking differences of the relevant binding energies, and hence are very sensitive to the precise details of the same. The fact that the theoretical separation energies obtained in this work differ from the corresponding experimental values only by a few hundred keV's is quite remarkable.

B. Composition of the Outer Crust of Neutron Stars
The masses of very neutron-rich nuclei are particularly interesting for some astrophysical calculations. We next compute the composition of the outer crust of a neutron star as a further application of our present mass model. As one moves from the surface of a neutron star to its interior, the outer crust is the region comprising matter at densities between ∼ 10 4 g/cm 3 and ∼ 10 11 g/cm 3 . Matter at those densities consists of fully-ionised, neutron-rich atomic nuclei that arrange themselves in the lattice sites of a Coulomb crystal embedded in a degenerate electron gas [33]. inner crust of the neutron star, where the atomic nuclei are immersed in an electron gas and a neutron gas.
In order to compute the composition of the outer crust we follow the usual formalism as described in Refs. [34][35][36] and references quoted therein. That is, we consider cold and electrically neutral matter which is assumed to be in thermodynamic equilibrium and in its absolute ground -state. We calculate the Gibbs free energy of this system by adding the contributions of the nuclear, electronic, and lattice terms [34][35][36] and, finally, we evaluate the equilibrium composition (Z,N) at a certain pressure by minimising the obtained Gibbs free energy per nucleon.
We display our predictions for the equilibrium nuclear species present in the outer crust in Fig. 7. We perform the calculations within the range ρ = 10 7 g/cm 3 to ρ = 3 × 10 11 g/cm 3 . The variation of the neutron and proton numbers with increasing crustal density shows a structure of plateaus that are interrupted by abrupt jumps in the composition. As exemplified by the N = 50 plateau, the prevalence of a given nucleon number over a large 16 range of densities is related with the shell effect due to the filling of a nuclear shell. The N = 50 neutron plateau also is very illustrative of the fact that, with increasing density, it is energetically favorable for the nuclei of the crust to capture electrons from the degenerate electron gas. This results in increasingly neutron-rich nuclides along the neutron plateau.
Eventually, the mismatch between the neutron and proton numbers is too large and the jump to the next neutron plateau takes place in an effort to reduce the penalty imposed on the system by the nuclear symmetry energy [35,36].  Fig. 7 for comparison. Though the overall pattern is quite similar to the results obtained with our calculated masses, the Möller-Nix mass table predicts more structure in the variation of the neutron and proton numbers with the crustal density, and the jump to the N = 50 plateau is delayed to a little higher density. This fact suggests that in the present mass region the shell effects due to the filling of nuclear shells and sub-shells are somewhat weaker in the Möller-Nix mass formula than in our model.

C. Superheavy Nuclei
Production and study of superheavy nuclei is of current interest from both theoretical [37][38][39][40] and experimental [19,41] aspects. With the advent of increasingly sensitive detection methods, it is possible to identify the superheavy elements, and measure α decay Q values precisely. The elements with Z = 118 have been produced so far [19]. Here, we apply our mic-mac model to a few recently reported superheavy nuclei [19]. In particular, we focus on the α decay Q values (Q α ). The binding energies of the parent as well as the daughter nuclei, necessary to obtain the Q α values, are obtained within our mic-mac model by minimising over the deformation (β 2 , β 4 ) mesh. The binding energy of the α particle is adopted from the Audi -Wapstra compilation [25]. The calculated (Calc.) as well as the experimental Q values [19] are presented in Table 2. We find that the calculated Q α values are very close to the experiment. This is quite encouraging, since as in the case of the separation energies, the Q values as well are obtained by taking differences between two large quantities. The α decay Q values can be related to the half lives through the Viola -Seaborg relation [42]. In particular, following Oganessian [19], we adopt: where, Z is the charge number of the parent nucleus; Q α is the α decay Q value, and a, b, c and d are parameters, taken to be [19]: Large scale calculations using the proposed mic-mac model can be cumbersome and highly time consuming. Therefore, it may be very useful to look for simplifications that allow to speed up the calculations without loss of accuracy. To this end, we explore the possibility of absorbing the fourth -order correction into the net second -order contribution: Here, (n) and (p) stand for neutronic and protonic contributions. See Eqs. (11)-(14) for the definitions of the different terms appearing in these two equations. This absorption is expected to have two major effects. Clearly, if such an absorption is possible, the factor α (see Eqs. (1) and (2) for definition), should be expressible as a function of neutron number, proton number, or some combinations thereof. Before discussing the possibility of absorbing fourth -order terms into second -order terms for a Woods -Saxon potential, we demonstrate the existence of such a functional form for the simple Harmonic Oscillator potential.

A. The Harmonic Oscillator Potential
The harmonic oscillator (HO) potential provides a unique opportunity to investigate the details of the WK expansions analytically. Therefore, first we consider the simplest form of the HO potential, without spin -orbit interaction. It can be shown that for the HO potential, the different WK corrections are given by [10]: where, λ is the chemical potential, determined as described earlier, and ω is the oscillator frequency. For the HO potential, assuming degeneracy of 2, the particle number (see Eq. (7)) is given by: we take an alternative and physically more transparent approach, wherein, we express λ as [15,43] where, λ j is correct up to order j . Starting from the Thomas Fermi expression for the chemical potential, and noticing that the normalisation is true order by order, we get the following expression for chemical potential, correct up to 4 : This, along with the second and fourth -order WK corrections to energy (see Eqs. (11), (12)), where, λ p and λ n are chemical potentials for Z protons and N neutrons respectively. Further, notice that the neutron and proton numbers can be written as: where, the terms up to the order A −2/3 are retained, and the expansion in I has been carried out only up to second -order in I. It can be therefore seen that the factor α can indeed be written as a function of mass number and I, implying that it is in principle possible, at least in the case of HO potential, to absorb the fourth -order WK corrections to the energy into the second -order WK corrections.
To understand the behaviour of α with respect to I, we plot the factor α as a function of I in Fig. 8. It is seen that the factor α has a very regular behaviour with respect to

B. Woods -Saxon potential
Next, we investigate the factor α for the Woods -Saxon potential. In order to achieve this, we choose a set of 2171 known nuclei [25] with Z > 5. Spherical symmetry is assumed.
The nuclear, spin -orbit and Coulomb potentials have been taken as defined in Ref. [6].
The full Wigner -Kirkwood calculations up to the fourth -order in are carried out for these nuclei, and the exact values of the factor α are obtained. These are then plotted as a function of the asymmetry parameter I in Fig. 9. The figure exhibits that the factor α has a very regular behaviour as a function of asymmetry. In order to understand the detailed structure of the factor α, we plot the same results with a greater resolution in Fig. 10.
A remarkable and regular pattern emerges from the plots. In comparison with the case of the HO potential, the pattern is inverted. The pattern consists of 'fan like' structures. There are groups of points stacked exactly along vertical lines, as indicated in Fig. 10 accompanied by symmetrically placed, slanting groups of points. All these groups of points constitute nearly perfect straight lines. This is in contrast with the case of HO potential, where the lines were curved.
A closer examination of the behaviour of the factor α reveals several interesting features. To understand them better, we shall first enlist the nuclei appearing in a particular 'fan' structure. We shall designate the slanting lines appearing in the 'fan' structure as 'rays'. to the mass number. Thus, it is expected that in the limit of A → ∞, the α values will approach some constant value, say, α 0 , which is approximately 1.125, according to the figure above.
Considering these observations, we propose the following parametrisation for the factor α: where, α j 's are adjustable parameters. Considering all the 2171 nuclei (see above), we carry out a least squares fit to determine these parameters. The fit turns out to be exceptionally good, with rms deviation 1.09 × 10 α 1 =2.26744; α 2 =-0.02659 and α 3 =0.29987. The difference between the exact and the corresponding fitted α values is plotted in Fig. 11, indicating that the agreement is almost perfect, and that the phenomenological formula that has been proposed here is indeed robust, for all the mass regions.
We shall now investigate the deformation effects particularly with reference to the factor α. In order to achieve that, we once again consider the set of 561 even -even nuclei (see The shell corrections require averaged energies, which are calculated here using the WK expansion. Here, we consider the WK expansion only up to second -order, and simulate the effects of fourth -order through the factor α (Eq. (34)). This defines the averaged energies and hence the shell corrections completely. The difference between the shell corrections thus obtained and the 'exact' shell corrections is found to be indeed small, the maximum deviation being of the order 150 keV, implying that the factor α obtained merely by using the spherical nuclei works very well for deformed systems as well (with both deformation parameters β 2 and β 4 ). This observation is indeed of great practical importance.
With these approximate shell corrections, we make a re-fit to the liquid drop parameters.
Comparison between the liquid drop parameters as reported in Section 3 and the ones obtained with the approximate shell corrections is presented in Table 3. It is indeed gratifying to note that the liquid drop parameters obtained in the two cases are almost identical, and so is the rms deviation of the calculated binding energies with respect to experiment [25].
This substantiates the validity of the parametrisation of α.
To test the robustness of the parametrisation of α further, we calculate the constants α j 's in Eq. (34) using just four nuclei ( 40 Ca, 100 Sn, 146 Gd and 208 Pb) instead of 2171 nuclei as described above. It is found that the numerical values of the constants practically remain the same. To test the validity of these parameters, the liquid drop parameters are re-worked employing the new values of α j 's. It is found that the liquid drop parameters thus obtained are practically equal to the ones reported in the right most column of Table 3. We close this section, by concluding that the absorption of fourth -order Wigner -Kirkwood corrections into the second contributions is reliable, and can be used in large scale mic-mac calculations. The absorption also has the advantage of reducing the numerical noise that might arise in the higher order derivatives of the potentials.

VI. SUMMARY AND CONCLUSIONS
The semiclassical Wigner -Kirkwood expansion of the one -body partition function has been employed instead of the Strutinsky averaging scheme to calculate the shell corrections within the framework of a mic-mac model. The microscopic part of the energy also contains pairing contributions that are obtained using the Lipkin -Nogami scheme. We have improved the macroscopic part of the model as compared with the one used in our previous work [6,7] by including the curvature term as well as the Wigner contribution. With just ten adjustable parameters, our model reproduces the binding energies of 561 even -even spherical and deformed nuclei with rms deviation of 610 keV. We have tested this new mic-mac model near the proton and neutron drip lines as well as in the superheavy region. Our present 25 calculations show that the mic-mac model proposed in this paper reproduces remarkably well the recent experimental results in these exotic scenarios.
Further, a systematic study of the ratio of the fourth -order and second -order Wigner -Kirkwood energies has been carried out. We find that the ratio of these two energies behaves in a very systematic manner. We have shown that this ratio can be parametrised accurately by a simple expression, implying that the fourth -order corrections can be absorbed into the second -order contributions in a very simple way. We have checked that using this simple procedure, we recover practically the same parameters of the macroscopic part, without deteriorating the quality of agreement achieved with the full Wigner Kirkwood calculation including explicitly the fourth -order contributions. Therefore, this simplified calculation of shell corrections can be used confidently in the large scale mic-mac calculations that we plan to carry out as the next step.
Finally, we point out that there is still some room for improving our model particularly in two specific directions. On the one hand, the full blocking procedure in the pairing calculations of odd -odd, odd -even and even -odd nuclei, that may be particularly relevant for spherical nuclei, has to be introduced. On the other hand, refinements in the mean field Woods -Saxon potential and in the distance function are still needed to study with our model not only neutron rich nuclei, but also fission barriers. This would require large scale calculations with the model, for which, the simplification proposed above may be very useful.