{'fraction_non_alphanumeric': 0.026806300175004863, 'fraction_numerical': 0.011333648156893247, 'mean_word_length': 4.547850208044383, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Many-body correlations in Semiclassical Molecular Dynamics and Skyrme forces for symmetric Nuclear Matter

Constraint Molecular dynamics CoMD calculations have been performed for symmetric nuclear matter (NM) by using a simple effective interactions of the Skyrme type. The set of parameter values reproducing common accepted saturation properties of nuclear matter have been obtained for different degree of stiffness characterizing the iso-vectorial potential density dependence. A comparison with results obtained in the limit of the Semi-Classical Mean Field approximation performed using the same kind of interaction put in evidence the role played by the many-body correlations present in the model explaining also the noticeable differences obtained in the parameter values in the two cases


Introduction
The description of many-body systems is one of the most difficult problems in nuclear physics due to the complexity of this kind of systems which are quantum objects described by a large number of degrees of freedom. A large variety of theoretical models have been developed using mean-field based and beyond-mean field approaches like Density Functional Theory [1] and Energy Density Function Theory [2,3]. In these approaches which use the independent particle approximation as a starting point, phenomenological effective interaction like Skyrme and Gogny forces are widely used taking advantage of their simple form [4]. In nuclear structure modeling we can quote the Skyrme-Hartree-Fock (SHF) methods, the relativistic mean-field (RMF) approach (as well as their Bogoliubov extensions), and the Hartree-Fock-Bogoliubov method with the finite-range Gogny force [5,6].
Complexity still becomes higher when nuclear dynamics, triggered by nuclear collisions, is studied to understand the property of nuclear forces far from the stability. At the Fermi energy and beyond semiclassical methods become necessary to describe the produced processes in which practically all the degrees of freedom are involved. Many-body correlations are responsible for the reorganization of the hot Nuclear Matter(NM) in to clusters. Two main classes of approaches have been developed to handle this complex scenario. The first one is based on the Boltzman transport equation including, at the higher energy, also three body collisions term [7,8,9]. In this case the theoretical approach is able to describe the time evolution of the one-body distribution function in phase-space. Apart from the collision term, the Boltzman equation corresponds to the semi-classical limit of time-dependent Hartree-Fock equations in phase-space. Different models have been implemented on this ground. They essentially differ from each other for the strategy adopted to simulate the collision terms and including the method used to represent the phase-space distribution starting from the nucleonic degree of freedom (test particle methods). As in Harthrre-Fock mean field methods the main ingredient concerning the interaction is the energy density which for a Skyrme type phenomenological interaction has a rather manageable form. In particular the local part of the interaction produces a simple functional for the potential energy expressed through an algebraic function of the density (see also Sec.II). In the early 90 a further development included stochastic forces (typical of Langevin processes) [10] in the so called Stochastic Mean Field [11] model to produce fluctuations capable of describing associated phenomena like mean-field instabilities leading to the cluster production. The second line of development of theoretical semiclassical approaches to deal the nuclear many-body problem is represented by the so called Molecular Dynamics models. The point of view of these methods is some sense antithetic compared to the ones based on the mean-field concept. However also in these cases, due to their simplicity the same kind of phenomenological effective interactions Skyrme and Gogny are widely used. The starting point of all these models is a basic assumption on the wave function describing the single nucleon which is represented through a wave packet well localized in phase-space with uncertainty satisfying the uncertainty principle. The manybody wave fucntion can be expressed through a direct product, leading to the so called Quantum Molecular Dynamics (QMD)-like models [12,13,14,15] or an anti-symmetrized product like in FMD and AMD [16,17]. Variational principles have been used to obtain the equation of motion for the the many-body functions. Within their own approximation schemes these models treat the many-body problem without the usage of the mean-field approximation. Many-body correlations are spontaneously produced and in general are able to describe the main features concerning the multi-breakup processes [18] observed in heavy ion collisions at intermediate energy and beyond. Due to these deep conceptual differences between the mean-field based models and the molecular dynamics ones it becomes interesting, and necessary in our opinion, to investigate the differences between the energy-density functionals which are produced with these two classes of approaches when one uses, in both the cases, the same kind of effective elementary interaction between nucleons. The study presented in this work is performed by using the Constraint Molecular Dynamical Model CoMD [13,14]. It is a part of a forthcoming more extended study involving also asymmetric NM. As it will be shown the presence of isovectorial forces strongly affects the results of this comparison (see also [19]). The quantitative results obviously will be related to some specificity of the model and to the used effective interaction, but as due to the common features which links most of the molecular dynamics models (semi classical wave packet dynamics), the obtained results could assume a more wide meaning at a qualitative level. The work is organized as follows: in Sec. II we illustrate the choice of the effective interaction and the related total energy functional are introduced. In Sec.III we describe the NM simulations. The results are presented in Sec. IV. Sec.V contains the summary and the concluding remarks.

The effective interaction and the total energy
Before illustrating with some detail the model and the NM simulations, in this section we briefly introduce the effective interaction which we use in our calculations. We also present the way in which the related total energy is obtained in the case of the semiclassical mean-field approximation (Se-MFA) and in the case of the semiclassical-wave packet dynamics.
The elementary interaction between two nucleons with spacial coordinates r,r ′ and third isospin component τ ,τ ′ is of the Skyrme type and has the following form: The first term of eq.(1) represents the iso-scalar contribution to the two-body interaction, the second one is the usual 3-body effective interaction. For this term the ρ density dependence is modeled through the parameter σ. This term appears to be a generalization of the Brink-Voutherin [3] three-body term corresponding to σ = 2. This generalization is sustained both on the basis of phenomenological parametrization used in mean-field approaches [8] and on more complex microscopic calculations based on G matrix calculations with realistic interactions. The third term represents instead the Iso-Vectorial contribution. The form factors F k have the following expression: ρ 0 is the saturation density value. Also in this case the F k form factors are introduced to reproduce the results of more complex microscopic calculations based on realistic interaction [20,21,22,23,24] concerning the symmetry energy and reflecting effects beyond the twobody interaction . These form factors have been widely used in Hartree-Fock calculations, in semiclassical BUU calculations and in others molecular dynamcs models. In particular, F 4 will be used, in the limit ρ/ρ 0 ∼ 1, to perform calculation with stiffness parameter γ different from the ones related to the others form factors (in the same limit F 1 , F 2 , F 3 correspond to γ values 1.5, 1 and 0.5 respectively). Finally we observe that for simplicity reason we do not add to this interaction non local effects. Regarding this assumption we note that some of the suggested Skyrme parametrizations have a vanishing or very small terms related to the non locality at the saturation density in the limit of Nuclear Matter (see for example Ref. [25]). Moreover as it will be shown in the following, in CoMD calculations the effective potential which determines the effective force able to govern the motion of the wave packet centroid is already non-local (sum of Gaussian contributions depending on the intra nucleons distances) being the results of the convolution of the well localized wave packets with δ function in space.
Even if well-known, we shortly recall in the next section the main ingredients of the Se-MFA. We think this section useful to present in the following the way in which this assumption is instead broken in the semiclassical wave-packets dynamics.

The Semiclassical Mean Field Approximation
At a fixed time starting from eq.(1) (we are considering here a stationary problem) , in a rather general way, we can obtain the expression for the total energy W related to the potential interaction by folding the elementary interaction with the two-body density distribution in phase-space D(r, r ′ , p, p ′ ). By taking into account the usual truncation condition on D typical of the mean-field approximation, we can write: where D 1 is the one-body distribution. Therefore by using eq.(1) for W we get: Taking in to account that by definition D 1 (r, p)dp = ρ(r) we obtain: 1 2 V (2) ρ 2 appear to be the energy density associated to the potential energy from which we can obtain the related binding energy per nucleon due to the potential interaction E pot represents the charge/mass asymmetry parameter evaluated from the neutron and proton density ρ n , ρ p respectively. The total binding energy E is obtained by adding the kinetic contribution coming from the Fermi motion In E F kin β terms of order greater than two are neglected. In particular we see how the isovectorial vectorial forces with strength proportional to T 4 contribute to the symmetry energy E sym depending on the asymmetry parameter β. For quadratic form in β we get: The common accepted bulk value of e sym at the saturation density is about 30 MeV even if relativistic Hartree-Fock models can predict higher values up to about 40 MeV [26,27]. However, still under study is the density dependence of this quantity which is able to affect neutron-skin thickness in nuclei and the Giant Monopole Resonances [26,27,28]. Finally we note that the structure of the obtained energy-density functional in Se-MFA makes independent the choice of the parameter describing the main properties of the symmetry energy from the other ones which instead are fixed from the saturation properties of the symmetric nuclear matter. We will show in the next section that this is no longer true in the CoMD approach.
In molecular dynamics approaches a strong assumption is done on the wave functions describing the nucleonic degree of freedom. It is commonly assumed that the wave function is a gaussian wave packet with fixed width σ r in coordinate space. The centroid in phase-space The Wigner transform of φ i is the widths in momentum and space satisfy the minimum uncertainty principle condition σ r σ p = 1 2h . Another assumption concerns the N-body Wigner distribution that in CoMD model [13](like in QMD approach [12]) is a direct product of the single particle distributions. In this case therefore for a system formed by A nucleons the one-body and 2-body distributions above introduced have a multi-component structure that is: From the above relations we can see that in general D = D 1 (r, p)D 1 (r ′ , p ′ ). In particular, for the case in which these distributions are expressed as a sum of different localized components inside a volume V g , we can indicate with a g the number of components which give non negligible contributions in V g . The relative difference associated to eq.(18 can be estimated to be of the order of 1/a g which corresponds to the ratio between diagonal (i = j) and off-diagonal elements (i = j) within the ensemble of a g (a g −1) terms. In semiclassical meanfield models a g can be enough high to make the difference negligible, in fact the single particle distribution usually spreads over the whole system (test particles methods) and therefore the truncation condition eq.(18) can be retained valid [7,8,9]. On the contrary this is not surely the case for the molecular dynamics approaches for which the typical spreading volume is of the order of 2-10 f m 3 . Localization and therefore coherence of the wave-packed used to describe the single-particle wave-functions allows to keep memory of the two-body nature of the interparticles interaction, and at same time, allows for the spontaneous appearance of the clustering phenomena in simulations concerning low-density and excited portion of nuclear matter. With these assumptions on the 2-body phase-space distribution, taking in to account the properties of the δ function, we can obtain the explicity expression for the different terms concerning the total energy; for the two-body isoscalar contribution W twb we get: In the above expression S i v is the normalized sum of the Gaussian terms and it represents just a measure of the overlap between the nucleonic wave-packets. Its two body character is quite explicit. In the calculations concerning the NM simulation that we will illustrate in the next sections , the large number of particles A involved in the system allows us to write the above quantity in a simpler way by introdcing the average overlap per nucleon S i v = S v whose dependence on the particle index in the ideal case can be omitted: By comparing the expression obtained for E twb in the two different appraches (eq.(10) and eq.(23)) we note a formal analogy where the variable ρ is substituted by the overlapp integral S v per nucleon. We however observe that this analogy is only formal, this aspect will be discussed in some detail in the next subsection. Concerning the three-body term, according to the evaluations reported in reference [12] and taking in to account the previous obervations we get: For the term related to the iso-vectorial interaction and for the most simple case F ′ = 1 in the limit A, N, Z >> 1 (N and Z represent the number of neutron and protons) we obtain: whereρ cc ′ with cc ′ equal to nn, pp and np represents the overlap integral per couples of neutrons, protons and neutron-proton. A more convenient form for the above expression is obtained by introducing the two followng quantities: for symmetric NM N = Z and we get: whereρ A ≡ Aρ. The expression in eq.(32) also contains a generalization to the cases in which we use the generic form factors F ′ k . Here F ′ k keep the same functional form in S v using the formal analogy ρ ρ 0 → Sv S v,0 above discussed. From eq.(32) we obtain that for α = 0 (as CoMD calculations predicts, see next sections) the iso-vectorial force produces a term which we name E bias . The related effect is not-negligible if compared to the balance of the different terms appearing in the expression of the total energy. The correlation coefficient α by definition (see eq.(30)) represents the difference in percentage of the overlap between the neutron-proton couples from the one related to the couples formed by homonym nucleons. It also depends on the strength of the iso-vectorial forces (T 4 parameter). Finally, the kinetic contribution E kin is obtained in a self-consistent way for ground state configurations by applying the cooling-warming procedure coupled with the constraint on the occupation numbers (see Ref. [14]) given by the Pauli principle. The total energy per nucleon is therefore obtained by adding all the contributions discussed in this section and it will be indicated as E C = E pot (S v , α,ρ A ) + E kin (ρ). In particular we note that, at this level, it depends on the new defined primary quantity S v , α,ρ A . In the following section we will try to relate these quantities to more fundamental ones that are the density ρ and the spatial correlation function ν between nucleon pairs.

Overlap integrals and spacial correlations
In a rather general way and as suggested from CoMD calculations (see Fig.1), for an uniform many-body system at a fixed density we can introduce the probability p to have a particle in the volume dV 1 localized in r 1 and a second one in the volume dV 2 localized in r 2 . Due to the uniformity condition p depends only on r = |r 1 − r 2 |. p can be expressed in the following form: with ν(0) = 1 and k 0 ≥ 0. Moreover lim r→∞ ν(r) = 0, that is: no spacial correlations can be expected for very distant particles. Non zero values of ν can be instead expected at relatively small distances due both to the interaction (for an attractive interaction the positive sign must be considered in eq.(33)) and to the symmetry of the many-body wave function for quantal systems of identical particles (in the case of identical Fermions we must consider the sign minus and k 0 = 1 ). For a classical system of non interacting particles we have a vanishing correlation effect p = 1 (k 0 = 0). For our aims, we need to evaluate a normalized probability P = cp in such a way: V c represents the volume in which the well localized correlation function ν is different from zero. This volume is always finite and of the order of the σ r 3 in our model calculations. So that in the limit V → ∞ we get c = ρ and : In Fig.1, just as an example, we display results of calculations for a given set of parameters for the Skyrme interaction (see the next section). The calculations show the probability to find two nucleons at a distance r in case of Pauli blocked couples P 1 (red points, neutron and proton couples with same spins), for the case of un-blocked proton and neutron couples P 0 (black point) and finally for neutron-proton couples P np (blue points). The calculations are referred to a symmetric portion of nuclear matter consisting of a sphere containing 2000 nucleons at a density of about 0.17f m −3 . The calculations include also the iso-vectorial potential energy (T 4 = 59 MeV in this example ). The probability distribution has been estimated by taking in to account only couples of nucleons in which at least one of the two nucleons is localized inside a smaller sphere, centered at the origin, with a radius of 4 fm. In this way we can neglect surface effects on this quantity. The error bars indicates the uncertainty related to the statistics of the simulations. We remark that in all the calculations concerning the present work in the range of density explored no cluster formation has been detected by our numerical procedure which was applied at each time step to check for the existence of aggregation processes. Only in the lower limit ρ = 0.7ρ 0 and for only part of the time, the formation of a big cluster with a mass about 98 % the total one and light particles have been observed. From the figure it is clearly seen the depletion of the P 1 distribution at small distances due to Pauli constraining in CoMD. For the other distributions an enhancement is produced due to the attractive effect of the two-body interaction. In particular we note that this enhancement is more pronounced for P np as compared to P 0 ; this is clearly an effect due to the iso-vectorial forces which generate a stronger attraction for neutron-proton couples. Now from eq.(35) we can obtain an expression for S v as: Figure 1. Probability distribution P 1 to find two identical nucleons at a relative distance r. In the same figure P 0 represents the probability evaluated for neutron and proton couples with opposite spin. Finally P np represents the probability obtained for neutron-proton couples. (Color) From the above equations we see that even if the overlap integral per nucleon is closely related to the density ρ, it admits correction terms related to the spacial correlation through the integral I. In general, for function ν(r) localized within a distance of the order of l the value of I decreases with the ratio x = σr l : for example, for an exponential profile e − r l the correction terms go to zero like e − σ 2 r 2l 2 . In other words, according to what observed in the previous section, when the single particles space distribution is enough large an averaging of the spacial correlations effect is obtained and the CoMD functional tends to the one represented in eq.(10) which is typical of the semi-classical mean-field approximation in transport theory. On the contrary for well localized wave-packets with σ r of the order of 1-2 f m, I is different from zero and reflects the behavior of ν(r) at small distances. The above relations can be also generalized for a multi-component system characterized by a charge/mass parameter different from zero.
The determination of the correlation functions ν(r) and of the related integral I is a problem which can not be solved in a general way. Some special cases are well known. They concern the study of the density fluctuations in the hydrodynamical limit valid for large distances compared to the mean free-path in non-equilibrium cases (see also [10,11]).However in this limit simultaneous spacial correlations are supposed to be zero at small distances. At short distances, comparable with the range of the effective interaction, as in our case, approximations schemes can be developed actually corresponding to a power decomposition in ρ as we will perform in the next section.

The nuclear matter simulation
As mentioned in the introduction the main goal of the present work is to understand also at a quantitative level and for a given simple form of the effective interaction, what are the consequences produced by the correlations above discussed by taking as a reference point the results obtainable in the framework of the Se-MFA approach. The study will be performed in a narrow region of density around the saturation density ρ 0 . To this aim in the following we will try to find the set of parameters for the effective interaction which reproduce in both the cases some of the commonly accepted saturation properties of symmetric nuclear matter. In particular we will refer to ρ 0 =0.165 f m −3 , binding energy E(ρ 0 ) = −16 MeV, compressibility modulus for symmetric nuclear mater K N M (ρ 0 ) = 220 MeV. The strength of the iso-vectorial interaction T 4 will be chosen in such way to produce a symmetry energy value e sym (ρ 0 ) ≃ 28.6 MeV in the case of the Se-MFA, i.e. T 4 = 32 MeV. From the functional E reported in eqs. (10)(11)(12) we can find easily the parameter values satisfying the above conditions by solving the system formed by the following set of equations for symmetric NM: In concrete cases due to the finite steps with which we perform the variation on the parameters the values of E(ρ 0 ) and K N M (ρ 0 ) are obtained within ±0.5%. Within the above specified uncertainty the solution gives the following values for the parameters: T 0 = −263 ± 1.3 MeV, T 3 = 208 ± 1 MeV and σ = 1.25 ± 0.02.

NM calculations and CoMD model
The evaluation of the total energy per nucleon E C related to the CoMD calculations requires the solution of the many-body problem using the equation of motion regulating the wave-packet dynamics. Details on the numerical procedure for the CoMD calculations are described in ref. [13,14]. At the different densities changing between 0.7-1.2ρ 0 with steps equal to 0.1ρ 0 , the calculations have been performed by enclosing a relatively large number of particles A 1 = 1600 and A 2 = 3560 in spherical volumes of radii R = r 0 A 1/3 with r 0 = ( 3 4πρ 0 ) 1/3 . Particles trying to escape from the spheres have been re-scattered inside through an elastic reflection at the surface. For symmetrical configuration, starting from the parameter values minimizing the functional E in eq.(10-12) we have searched for the stationary minimum energy conditions by applying the cooling-warming procedure coupled with the constraint related to the Pauli principle as described in Ref. [13]. Calculations have been performed for the two systems having number of particles equal to A 1 and A 2 . The value of T 4 has been fixed to 32 MeV and the calculations have been performed for different values of the stiffness parameter γ. From the minimum energy configurations we have evaluated the related intensive quantities S v (ρ), α(ρ),ρ A (ρ) and E kin . Corrections due to the surface effects, which are necessary to estimate the associate bulk values have been evaluated using the following relation: Q i indicates the quantity valuated for the system with mass A i (i=1,2). Q b and Q s are the bulk and surface coefficients. Effects related to the curvature are represented by the last term of eq.(42). The coefficient 0.45 has been deduced performing a couple of calculations in boxes having the same volume as the considered spheres.
As an example, for β = 0 and γ = 1 (form factor F 2 ) in Fig.2 we show as a function of the reduced density the values of S v (ρ), α(ρ), ρ A (ρ) evaluated for the systems 1 (marked points with star), for the system 2 (marked points with square) and the bulk estimated values (dot points). Corrections of the same order are evaluated also for the kinetic energy contribution. In the following for simplicity we will refer to the bulk quantities without using the subscript b. After this first step, we fit the evaluated bulk quantities with a polynomial function of the reduced density ρ ρ 0 . The above quantities show in fact deviations from a simple linear behavior as a function of ρ. This can be seen by looking at the red line in the figure (ρ as a function of ρ/ρ 0 ). A second order polynomial reproduces very well the behavior in the range of explored densities. The results of the fit are shown in the figures with lines.
The obtained functions are then substituted in E C and therefore the total binding energy can be now evaluated with continuity as a function of ρ and we can search for the parameter values solving the CoMD functional E C obtained in Sec. II.1 and satisfying the conditions expressed in eqs. (38-41). We have to note that the numerical solution of this system of coupled equations in general can not be obtained with the same precision as the one involving the functional E. In most of the cases, depending on the stiffness parameter γ it is possible to obtain solutions reproducing the E(ρ 0 ) and ρ 0 values within 10% while a larger spread is obtained for the K N M (ρ 0 ). The chosen solutions will be the one which minimize the total relative difference from the reference values. Having found the best solution for the functional E C , in the sense above specified, with the new set of parameter values for T 0 , T 3 and σ we perform another series of microscopic NM simulations on the systems of A 1 ,A 2 particles. After having included the correction for surface effects, we do the polynomial fit. Then using the new calculated quantities

Results and discussion
In this section we illustrate the results obtained from the recursive procedure previously described.

Results on the primary quantities
As an example in Fig.3 we show the final values of S v as a function of the reduced density obtained in the case of symmetric NM and for different form factors F k (T 4 = 32 MeV). The lines represent the results of the fit with a second order polynomial. The red line represents instead the linear relation corresponding to the density ρ. As we can see S v shows deviations from ρ depending also on the used form factor.
Under the same conditions in Fig.4 with black line and points, we show the value of α as a function of the reduced density. The red line (and points) represents the values of α obtained in the case of T 4 = 0 MeV using the form factor F 2 . The blue line represents the values obtained for T 4 = 59 MeV which corresponds in the case of the Se-MFA to a value of e sym (ρ 0 ) equal to 42 MeV.
As can be seen the value of α decreases with the density and increases with T 4 . For T 4 = 0 finite values of α are essentially due to correlations imposed by the Pauli principle in the system of interacting particles.

Total binding energy
In Fig.5, upper panel, we show with a solid line the total energy per nucleon E as a function of the reduced density (eqs. (11,12)) satisfying the requested conditions on NM saturation properties. This result represent our reference point. The curves with broken lines represent the results of our NM simulations with CoMD model at the first step of the iteration for the different indicated form factors F k . In this case the parameter used for the effective interaction are the same like the ones obtained from the minimization of the energy functional related to the Se-MFA. As can be seen the curves in this cases are rather different. In particular for F 1 and F 2 the minimum energies are shifted to lower values and the compressibility modulus also shows large deviations compared to the chosen reference value. F 3 shows instead even a maximum at the saturation density. In the lower panel we show the corresponding values of E bias (see eq.(32)). Together with the density dependence of the primary quantities above described, this term is the main co-responsible of the observed deviations. The figure shows that the density behavior of E bias strongly depends at a quantitative and qualitative level on the used form factors F k . In particular we note an increasing slope from positive to negative values with the increasing of the degree of stiffness γ characterizing the behavior of the iso-vectorial forces. We can expect therefore that the necessary correction on the parameters of the Iso-scalar effective interaction to reproduce the reference properties of the symmetric NM will show a dependence on γ.
In Fig.6 we show analogues results after that the self consistent iteration procedure has been completed as described in Sec. III . In the interval γ ≃ 0.85 − 1.5 the saturation density and binding energies are reproduced within some percent. The compressibility is instead obtained within about 20%. In the figure we show some example of these solutions with black lines.  In general we note a more pronounced asymmetry of the curves around the saturation density with a reduced slope of the lower density branch and an increased slope for reduced density larger than 1. This behaviors is a direct consequence of the density dependence of the average overlap integral (see Fig.5). For stiffness parameter values out of the indicated interval we observe a fast increasing of the compressibility (well beyond 400 MeV) and of the corresponding binding energy at the saturation density. In Fig.6 an example showing this trend is represented by colored lines obtained for γ = 0.75 and γ = 1.75. All these circumstances shows that, in the frame work of the present molecular dynamics approach, the F 4 form factor with stiffness parameter values external to the range γ ≈ 1.5 − 0.85, due to the correlations discussed, is not able to reproduce the chosen commonly accepted saturation properties of the symmetric nuclear matter. These results are consistent with recent findings on the study of the 40,48 Ca + 40,48 Ca systems at 25 MeV/nucleon [29,30] concerning the balance between the yields of incompletefusion and multi-breakup processes. Finally in Fig.7 we show as a function of the γ parameter the values of the T 0 , T 3 and σ obtained from the iterative procedure described in Sec. III.1. Apart from the extremal plotted values, corresponding to high values of the compressibility ( see Fig.6), the internal ones correspond to saturation density, binding energies and K N M values well within 20% of the value obtained in the case of the Se-MFA. From the figure we observe a dependence of the parameter values describing the iso-scalar forces on the stiffness parameter γ associated to the iso-vectorial interaction. In the internal region of the explored interval, even if the dependence is moderate, maximum changes of the order of 16% and 20% are obtained for the T 3 and σ values respectively. However it is remarkable that the average values of the isoscalar interaction parameters show large differences compared to the reference values obtained in the case of the Se-MFA (T 0 = −263 MeV,T 3 = 208 MeV and σ = 1.25. see Sec.III) which are instead independent on γ.

Summary and Concluding Remarks
Many-body correlations produced in Molecular Dynamics approach based on the CoMD model have been discussed and their connection with the used effective interaction have been analyzed in the case of symmetric nuclear matter simulations. This study has been performed by comparing the results obtained for the total energy functional and the nuclear matter (NM) main saturation properties with the ones obtainable in the case of a semiclassical mean field approximations (Se-MFA). This comparison has been performed by using the same kind of simple effective Skyrme interaction. While we can expect the two approaches to produce large differences at low density due to the cluster formation process, the following study shows that noticeable differences are obtained also in a narrow range of densities around the saturation one where no cluster production has been observed. The effective Skyrme interaction include two-body and three-body effective iso-scalar interactions plus a two-body iso-vectorial one with form factors commonly used also in Se-MFA. In CoMD model calculations, around the saturation density, the effects related to the spacial correlations generated by the usage of the localized wavepackets and the ones associated to the multi-particles correlations produced through the Pauli principle constraint can be well described through a second order polynomial decomposition of the total energy as a function of the density. The obtained results show that, contrary to the case of the Se-MFA, the discussed correlations produce an interdependence between the parameters describing the iso-scalar forces and the ones related to the iso-vectorial interaction. The usage of an iterative procedure tuned to obtain in both the cases (Se-MFA and in CoMD approaches) very similar saturation density, total energy and compressibility values for symmetric NM allowed to extract the "good" set of parameters used for CoMD calculations. The obtained values differ under many aspects from the ones obtained from the Se-MFA: in particular the density dependence of the used form factors describing the iso-vectorial forces has to change now in a more restricted range of stiffness values. Moreover the values of the coefficients describing the Figure 7. The values of the parameters T 0 T 3 and σ (panel a,b,c respectively) obtained through the iterative procedure applied to CoMD calculations (see the Sec. III.1)are shown as a function of the stiffness parameter γ. The values of γ equal to 1 and 1.5 are associated to the functional form F 2 and F 1 respectively (see Sec.II). The others ones are instead associated to F 4 . The bar errors represent the global uncertainty on the parameter values related to the numerical procedures iso-scalar interactions are rather different in the two cases. Work is in progress to extend these studies to asymmetric NM. Finally we conclude by observing that even if from a numerical point of view the obtained results are strictly valid for the CoMD model, the performed study shows that the observed differences in the parameter values describing the chosen effective interaction can have a more wide meaning. They are indeed strictly linked to some general properties of the semiclassical wave packets dynamics.