Empirical information on nuclear matter fourth-order symmetry energy from an extended nuclear mass formula

Abstract We establish a relation between the equation of state of nuclear matter and the fourth-order symmetry energy a s y m , 4 ( A ) of finite nuclei in a semi-empirical nuclear mass formula by self-consistently considering the bulk, surface and Coulomb contributions to the nuclear mass. Such a relation allows us to extract information on nuclear matter fourth-order symmetry energy E s y m , 4 ( ρ 0 ) at normal nuclear density ρ 0 from analyzing nuclear mass data. Based on the recent precise extraction of a s y m , 4 ( A ) via the double difference of the “experimental” symmetry energy extracted from nuclear masses, for the first time, we estimate a value of E s y m , 4 ( ρ 0 ) = 20.0 ± 4.6 MeV . Such a value of E s y m , 4 ( ρ 0 ) is significantly larger than the predictions from mean-field models and thus suggests the importance of considering the effects of beyond the mean-field approximation in nuclear matter calculations.

To the best of our knowledge, unfortunately, there is so far essentially no experimental information on the magnitude of E sym,4 (ρ), even at normal nuclear density ρ 0 . Theoretically, the mean-field models generally predict the magnitude of E sym,4 (ρ 0 ) is less than 2 MeV [16,[18][19][20]. A value of E sym,4 (ρ 0 ) = 1.5 MeV is obtained from chiral pion-nucleon dynamics [21]. The recent study [22] within the quantum molecular dynamics (QMD) model indicates that the E sym,4 (ρ 0 ) could be * Corresponding author (email: lwchen@sjtu.edu.cn) as large as 3.27 ∼ 12.7 MeV depending on the interactions used. Based on an interacting Fermi gas model, a significant value of 7.18 ± 2.52 MeV [23] is predicted for the kinetic part of E sym,4 (ρ 0 ) by considering the highmomentum tail [24] in the single-nucleon momentum distributions that could be due to short-range correlations of nucleon-nucleon interactions. In addition, the divergence of the isospin-asymmetry expansion of nuclear matter EOS in many-body perturbation theory is discussed in Refs. [21,25]. Therefore, the magnitude of E sym,4 (ρ 0 ) is currently largely uncertain and it is of critical importance to obtain some experimental or empirical information on E sym,4 (ρ 0 ).
Conventionally nuclear matter EOS is quantitatively characterized in terms of a few characteristic coefficients through Taylor expansion in density at ρ 0 , e.g., The density in the interior of heavy nuclei is believed to nicely approximate to saturation density of symmetric nuclear matter (nuclear normal density) ρ 0 and the empirical value of ρ 0 ≈ 0.16 fm −3 has been obtained from measurements on electron or nucleon scattering off heavy nuclei [26]. Our knowledge on nuclear matter EOS largely stems from nuclear masses based on various nuclear mass formulae. By analyzing the data on nuclear masses with various nuclear mass formulae (see, e.g., Ref. [27]), consensus has been reached on E 0 (ρ 0 ) and E sym (ρ 0 ) with E 0 (ρ 0 ) ≈ −16.0 MeV and E sym (ρ 0 ) ≈ 32.0 MeV. These empirical values on E 0 (ρ 0 ) and E sym (ρ 0 ) are of critical importance for our understanding on nuclear matter EOS.
Generally speaking, it is very hard to determine the higher-order parameter E sym,4 (ρ 0 ) and the fourth-order symmetry energy a sym,4 (A) of finite nuclei from simply fitting nuclear masses within nuclear mass formulae since the term a sym,4 (A)I 4 (I = N −Z A with N and Z being the neutron and proton number, respectively, and A = N + Z is mass number) is considerably small compared to other lower-order terms in the mass formula for known nuclei, even for the predicted dripline nuclei [28]. Recently, however, by approximating a sym,4 (A) to a constant c sym,4 in the mass formula, several studies [29][30][31][32] have been performed to extract c sym,4 from analyzing the double difference of the "experimental" symmetry energy extracted from nuclear mass data, and robust results with high precision have been obtained, i.e., a sizable positive value of c sym,4 = 3.28 ± 0.50 MeV or 8.47 ± 0.49 MeV is obtained in Ref. [29], depending on the Wigner term form in the mass formula. More recently, a value of c sym,4 = 8.33 ± 1.21 MeV is extracted in Ref. [32] using similar analysis on nuclear masses. These results provide the possibility to extract information on E sym,4 (ρ 0 ). In this work, by self-consistently considering the bulk, surface and Coulomb contributions to the nuclear mass, we extend the mass formula of Ref. [33] to additionally include the corrections due to central density variation of finite nuclei and the higher-order fourth-order symmetry energy term a sym,4 (A)I 4 . In this extended mass formula, a explicit relation between a sym,4 (A) and E sym,4 (ρ 0 ) is obtained. We demonstrate for the first time that the precise value of c sym,4 obtained recently from nuclear mass analysis allows us to estimate a value of E sym,4 (ρ 0 ) = 20.0 ± 4.6 MeV.

II. NUCLEAR MASS FORMULA
There have been a number of nuclear mass models which aim to describe the experimental nuclear mass database and predict the mass of unknown nuclei. Nowadays, some sophisticated mass formulae [27,[35][36][37] (with shell and pairing corrections) can reproduce the measured masses of more than 2000 nuclei with a root-meansquare deviation of merely several hundred keVs. These mass formulae provide us empirical information about the EOS of nuclear matter, especially its lower-order characteristic parameters E 0 (ρ 0 ), E sym (ρ 0 ) and so forth.
To relate the coefficients in the mass formula to the EOS of nuclear matter, one can express the binding energy B(N, Z) of a nucleus with N neutrons and Z protons in terms of the bulk energy of nuclear matter in the interior of the nucleus plus surface corrections and Coulomb energy. Based on such an argument, Danielewicz [33] developed a mass formula with a selfconsistent A-dependent symmetry energy a sym (A) of finite nuclei. Considering that the central density ρ cen in nuclei generally depends on N and Z and deviates from ρ 0 , we here extend the mass formula of Ref. [33] by considering the deviation of ρ cen from ρ 0 , and additionally including the higher-order I 4 terms. In such a framework, a nucleus with N neutrons and Z protons is assumed to localize inside an effective sharp radius R, i.e., where r 0 is a constant satisfying 4 3 πρ 0 r 3 0 = 1 and χ cen = (ρ cen − ρ 0 )/3ρ 0 is a dimensionless variable characterizing the deviation of ρ cen from ρ 0 . Furthermore, we denote the volume (surface) neutron excess as and N s (Z s ) represent the neutron (proton) number in the volume and surface regions of the nucleus, respectively, with N v + N s = N and Z v + Z s = Z. Generally, χ cen and ∆ v (∆ s ) depend on N and Z of the nucleus and can be determined from equilibrium conditions, and this is consistent with the argument of the droplet model (see, e.g., Ref. [34]).
In the present work, the nuclear binding energy consists of volume term B v , surface term B s and Coulomb term B c . The volume part of the binding energy can be treated in nuclear matter approximation, i.e., The surface term comes from surface tension and symmetry potential (detailed argument can be found in Ref. [33]), and it can be expressed as where σ 0 (σ I ) represents the isospin independent (dependent) surface tension, S = 4πR 2 is the surface area of the nucleus, and we define E s0 = 4πr 2 0 σ 0 and β = σI (1) has been used to obtain the second line in Eq. (3). For Coulomb energy, for simplicity we adopt the following simple form without exchange term, i.e., with a c = 3 5 e 2 4πǫ0r0 . The equilibrium condition of nuclei can be obtained from variations of the binding energy B(N, Z) of the nucleus with respect to χ cen and ∆ v , i.e., from which we can obtain χ cen and ∆ v (∆ s ) for different A and Z. The first equation means the mechanical equilibrium and tells us how the surface energy, Coulomb energy and the isospin dependent part of volume energy affect the value of ρ cen , while the second equation represents the balance of the isospin asymmetry chemical potential between the volume and surface regions. To solve Eq. (5), we expand χ cen in terms of ∆v A , and then expand ( ∆v A ) 2 in terms of I, i.e., where the expansion coefficients χ 0 , χ 2 and D 2 might depend on A or Z, consistent with calculations from the droplet model [34] and the Thomas-Fermi approxima-tion [38]. Using Eqs. (2), (3) and (4) and substituting B(N, Z) = B v + B s + B c into Eq. (5) leads to the following two equations By eliminating ∆s A in Eq. (8) and Eq. (9), we obtain From Eqs. (10) and (6), one can obtain χ 0 and χ 2 , i.e., which represents the modification of the central density of finite nuclei due to the surface and Coulomb energy, and which determines the modification of nuclear central density due to the isospin dependent part of nuclear matter EOS. On the other hand, one can figure out the relation between ( ∆v A ) 2 and I from Eq. (9), i.e., from which, together with Eq. (7), one then obtains As can be seen from Eq. (7), the D 2 (A) reflects the fraction of volume neutron excess in the total neutron excess of a nucleus with mass number A. The ratio E sym (ρ 0 )/β is the so-called symmetry volume-surface ratio [33], which can be considered as an independent parameter. From Eqs. (11), (12) and (14), χ cen and ∆ v in Eqs.
where c 00 and c 01 characterize the isospin independent parts of B(A, Z) with while the symmetry energy a sym (A, Z) and the fourthorder symmetry energy a sym,4 (A) of finite nuclei can be expressed, respectively, as Noting that a sym (A, Z) in Eq. (18) includes a small Zdependent term a c A −4/3 Z 2 χ 2 D 2 , which comes from the modification of ρ cen due to the Coulomb energy. The mass formula Eq. (15) is an extended form of the mass formula of Ref. [33] and the latter can be obtained from the former by omitting the corrections due to the central density variation of finite nuclei (i.e., setting χ 0 = χ 2 = 0) and the higher-order I 4 term. In particular, by setting χ 2 = 0, a sym (A, Z) is then reduced to a simpler form of a sym (A) [33], i.e., In addition, in the limit of A → ∞, both the χ 0 (see Eq. (11)) and the D 2 (see Eq. (14)) become to zero, Eq. (15) is then reduced to the binding energy per nucleon of asymmetric nuclear matter at saturation point where the B(A, Z)/A reaches its minimum value, i.e., In the above, the Coulomb interaction is neglected and the I is replaced by δ. It is seen that Eq. (21) is exactly the same as the expression obtained in Ref. [18] for the binding energy per nucleon of asymmetric nuclear matter at saturation point. It should be pointed out that the A-dependence of a sym (A, Z) and a sym,4 (A) mainly comes from D 2 (A) in Eq. (14), and the same form of D 2 (A) has been obtained in Ref. [33]. In nuclear mass formula, the nuclear binding energy per nucleon is usually expanded in two small quantities A −1/3 and I 2 [39,40] with the coefficient of each term being determined by fitting nuclear masses with optimization. However, since the value of Esym(ρ0) β is around 2.4 [33,41], (14) is thus not small enough to obtain a rapid converging A −1/3 expansion. Several studies indicate the convergence of A −1/3 expansion is not as good as that of I 2 (see, e.g., Ref. [34]). For example, the value of a sym (A) (neglecting the small Z-dependent term) in the range of known nuclei is found to be very different from its asymptotic value at infinite A (or a sym (∞) in some literatures) [42][43][44]. Considering the slow convergence of the A −1/3 expansion, in our mass formula we do not expand a sym (A, Z) and a sym,4 (A) in terms of A −1/3 and exactly retain their A dependence.
We would like to point out when the nuclear binding energy B(N, Z) is expanded in terms of A −1/3 and I 2 , the obtained expressions for the coefficients of the five leading-order terms (i.e., A, A A × I 4 ) in the mass formula Eq. (15) are complete and self-contained, which means all other characteristic parameters of nuclear matter EOS and surface energy that are not show up in the mass formula (e.g., the higherorder J 0 and K sym ) will not contribute to the coefficients of the five leading-order terms.

III. SYMMETRY ENERGY OF FINITE NUCLEI AND Esym(ρ)
When the semi-empirical mass formula was first introduced, the symmetry energy term has the simple form of c sym (N − Z) 2 /A, with a constant symmetry energy coefficient c sym . By fitting the newly released nuclear mass table AME2012 [45] (all nuclei with A > 20 are considered) with the following simple Bethe-Weizsäcker mass formula, i.e., Here the constant c sym is just a parameter in mass formula and cannot be simply considered as E sym (ρ 0 ) in the EOS of nuclear matter since the symmetry energy coefficient a sym (A, Z) in the mass formula is sensitively dependent of the mass number A in the mass region of known nuclei as shown in Eq. (18).
It is constructive to figure out the relation between c sym and a sym (A, Z). Since each nucleus was considered equally in our simple fitting, the c sym can be treated as arithmetic average of a sym (A, Z), i.e., where N MN = 2348 is the number of measured nuclei we used in our simple fitting (i.e., the measured nuclei in AME2012 with 20 < A < 270) and the sum runs over all these nuclei. Substituting Eq. (18) into Eq. (23), one can then obtain where the summations are the same as that in Eq. (23) and the average X(A, Z) is defined as X(A, Z)/N MN . Similarly, one can obtain the following relations and Using E sym (ρ 0 )/β = 2.4 [32,33,41], we find , we finally obtain c sym ≈ 21.8 MeV, which is in good agreement with the value 22.2 MeV obtained by simple fitting. Our results also indicate that the E s0 and a c can be nicely approximated, respectively, by c sur and c cou in the simple Bethe-Weizsäcker mass formula. The above demonstration suggests that information on E sym (ρ 0 ) can be extracted inversely from the constant c sym . This feature is rather valuable in the case of a sym,4 and it provides an approach to extract E sym,4 (ρ 0 ) through the obtained constant fourth-order symmetry energy coefficient c sym,4 in the mass formula.

IV. FOURTH-ORDER SYMMETRY ENERGY OF FINITE NUCLEI AND Esym,4(ρ)
Using the empirical values of K 0 = 240 MeV, E sym (ρ 0 )/β = 2.4 and L = 45.2 MeV, we show in Fig. 1 the fourth-order symmetry energy a sym,4 (A) of finite nuclei as a function of mass number A for E sym,4 (ρ 0 ) = 0.0 MeV, 10.0 MeV, 20.0 MeV and 30.0 MeV. Also included in Fig. 1 is the constant fourth-order symmetry energy coefficient c sym,4 = 3.28±0.50 MeV obtained in Ref. [29]. The inset of Fig. 1 displays the asymptotic behavior of the obtained a sym,4 (A) as a function of the mass number A for E sym,4 (ρ 0 ) = 20.0 MeV. As shown in the inset, the value of a sym,4 (A) in the region of A < 270 is only about one fourth of a sym,4 (∞), indicating the very slow convergence of the A −1/3 expansion for a sym,4 (A). Similarly as in our simple fitting of c sym in section III, when extracting the constant fourth-order symmetry energy coefficient c sym,4 from analyzing nuclear masses in Ref. [29][30][31][32], each nucleus in AME2012 with 20 < A < 270 is also considered equally. Therefore, c sym,4 can be also treated as arithmetic average value of the Adependent a sym,4 (A) (i.e., Eq. (19)). Following the discussion in Section III, one can then extract information on E sym,4 (ρ 0 ) from c sym,4 . It is interesting to see from Fig. 1 that although the a sym,4 (A) exhibits Adependence, it changes not so rapidly in the mass range of A < 270. Such a feature makes it more reliable for us to estimate the value of E sym,4 (ρ 0 ) through the obtained constrains on c sym,4 .
We would like to point out that the detailed value of E sym,4 (ρ 0 ) = 20.0 ± 4.6 MeV estimated above relies on the empirical values of the lower-order parameters L and K 0 of nuclear matter EOS as well as the extracted values of E sym (ρ 0 )/β and a sym,4 from analyzing nuclear mass data. For example, if we choose L = 58.7 ± 28.1 MeV [11] or L = 58.9 ± 16.5 MeV [48], the obtained E sym,4 (ρ 0 ) changes to 23.0 ± 8.1 MeV or 23.0 ± 5.9 MeV, respectively. Nevertheless, the choice of L does not change the constrain on E sym,4 (ρ 0 ) much, and our results clearly indicate that a sizable positive E sym,4 (ρ 0 ) is necessary to describe the value of c sym,4 obtained from the double difference of the "experimental" symmetry energy extracted from the nuclear mass data. Considering the fact that the majority of nuclear energy density functionals based on mean-field models give a fairly small magnitude of E sym,4 (ρ 0 ) with its value less than 2 MeV [16,[18][19][20] (Note: E sym,4 (ρ 0 ) = 2 MeV leads to c sym,4 ≈ a sym,4 = −0.45 MeV), the effects beyond the mean-field approximation, such as the shortrange correlation effects, might be needed to explain such a sizable E sym,4 (ρ 0 ). On the other hand, if E sym,4 (ρ 0 ) is indeed very small, then a novel mechanism is called for to explain the large value of c sym,4 from analyzing the data on nuclear masses.

V. CONCLUSION AND OUTLOOK
By self-consistently considering the bulk, surface and Coulomb contributions to the nuclear mass, we have obtained an extended nuclear mass formula. In this mass formula, the symmetry energy a sym (A, Z) and the fourthorder symmetry energy a sym,4 (A) of finite nuclei are related explicitly to the characteristic parameters of nuclear matter EOS. In particular, using the recently extracted constant fourth-order symmetry energy coefficient c sym,4 from analyzing the double difference of the "experimental" symmetry energy extracted from nuclear masses, we have estimated for the first time a value of E sym,4 (ρ 0 ) = 20.0 ± 4.6 MeV for nuclear matter fourthorder symmetry energy at nuclear normal density ρ 0 .
The significant value of E sym,4 (ρ 0 ) = 20.0 ± 4.6 MeV challenges the mean-field models which generally predict E sym,4 (ρ 0 ) 2 MeV. Therefore, it will be interesting to explore E sym,4 (ρ 0 ) within the framework of beyond the mean-field approximation (e.g., by considering the short range correlation effects). This would substantially improve our understanding on the properties of nuclear matter systems at extreme isospin, such as neutron stars.