I. INTRODUCTION

The symmetry energy term in the nuclear equation of state (EOS) intimately relates to a wealth of the dynamical process of nuclear reactions, the structure character of nuclei and astrophysical phenomena [1]. The investigations of the symmetry energy, especially focusing on its density dependence, have been conducted through many observables such as isotopic ratio [2], isospin diffusion [3], neutron-proton emission ratio [4, 5], giant monopole resonance [6], pygmy dipole resonance [7], giant dipole resonance [8], collective flow [9] and isoscaling [10-12]. In our recent works, we utilized the isotopic yields to extract the symmetry energy coefficient relative to the temperature, a_sym/T, in the framework of the nuclear phase transition theory [13-17]. Along this scenario, a_sym/T was extracted, using m-scaling [18] and isobaric yield ratios [19]. Ono et al. also independently introduced a generalized free energy, K(N,Z), and extracted the a_sym/T values from their quadratic distributions [20]. In this method, for a given Z, all contributions from the volume, surface, Coulomb and pairing terms in the free energy are squeezed into an N proportional term and a constant, ξ(Z)N + η(Z), as a course approximation.

While temperature is one of the key variables in characterizing nuclear reactions, it is very difficult to determine the temperature of hot nuclear matter in a dynamical process. Several nuclear thermometers have been proposed [21]. These include the slope of energy spectra [22, 23], momentum fluctuations [24, 25], double isotope yield ratios [26] and excited state distributions [27]. However they may not be generally applicable in all circumstances and even for a given system the extracted temperature values from these thermometers may be quite different from each other [21]. In our recent work, temperature is evaluated from a (μn-μp)/T analysis, but no density determination was possible [28]. In our another recent work, the isotopic yield ratio method was applied to extract a_sym/T values from the experimentally reconstructed primary fragment yields. These ratios were compared to those calculated from AMD primary generated fragment yields obtained using Gogny interactions with different density dependencies of the symmetry energy [29, 30]. In the analysis, we found that the extracted a_sym/T values change according to the interactions used. From the dependence on interaction, the density, symmetry energy and temperature at the time of fragment formation were determined in a self-consistent manner.

In this article, an improved method is used to extract the a_sym/T values. In the improved method, all available isotope yields are employed. The improvement is made possible by taking into account the mass dependent temperature [31] in the free energy in an iterative manner. Using the self-consistent manner, the density and the temperature at the time of fragment formation are carefully determined from the obtained a_sym/T values from different interactions with different density dependencies of the symmetry energy term, i.e., the standard Gogny interaction which has an asymptotic soft symmetry energy (g0), the Gogny interaction with an asymptotic stiff symmetry energy (g0AS) and the Gogny interaction with an asymptotic super-stiff symmetry energy (g0ASS) [12, 32]. Our present analyses are conducted in the framework of AMD. There are three major reasons to use AMD as the event generator for this work. One is its capability to reproduce the experimental isotope yields. AMD results, such as multiplicity, angular distribution and energy spectra, have often been compared with those from the experimental data for intermediate energy heavy ion collisions and reproduce them reasonably well [32-39]. In one of our recent works in Ref.[30], the yields of the experimentally reconstructed primary hot isotopes are well reproduced by those of the AMD simulations. Second is to eliminate the secondary cooling effect. As shown in Ref.[19, 30], the sequential decay of the primary hot isotopes significantly alters the yield distribution and distorts the information inherent in the primary hot fragment yields. Third is to simplify the initial conditions using zero impact parameter to eliminate the effects of transverse flow and neck emission among others. Thus the AMD events are generated for central collisions (b=0 fm) of ⁴⁰Ca + ⁴⁰Ca at 35 MeV/nucleon.

II. IMPROVED MFM MODEL AND EXTRACTION OF A_SYM/T

In the framework of MFM, the yield of an isotope with mass A and I=N-Z (N neutrons and Z protons) produced in a multifragmentation reaction, can be given as [14, 16-19, 40]

Using the generalized Weiszäcker-Bethe semiclassical mass formula [41, 42], W(I,A) can be approximated as

In Eq.(1), A^-τ and Nln(N/A)+Zln(Z/A) originate from the increases of the entropy and the mixing entropy at the time of the fragment formation, respectively. μn (μp) is the neutron (proton) chemical potential. τ is the critical exponent. In this work, the value of τ=2.3 is adopted from the previous studies [40]. In general coefficients, av, as, a_sym, ap and the chemical potentials are temperature and density dependent. In this formulation a constant volume process at equilibrium is assumed in the free energy and therefore the term "symmetry energy" is used throughout this work along Ref. [43]. If one assumes a constant pressure at the equilibrium process [44], the term "symmetry enthalpy" should be used. Experimentally, whether the equilibrium process takes place at constant pressure or volume can not be determined, and thus we use "symmetry energy" throughout the paper, keeping in mind the ambiguity [43].

In the previous analyses [19, 30, 45, 46], the temperature in Eq.(1) was assumed to be identical to the temperature of the fragmenting source and treated as a constant for all isotopes. However as seen below, this temperature turns out to be (fragment) mass dependent. This mass dependence on the temperature was not recognized in the previous analyses, just because the mass dependence was masked by the larger error bars. In this improved method, the error bars become small and the mass dependence becomes evident. In order to take into account the mass dependence of the temperature in Eq.(1), the temperature T is replaced by an apparent temperature T(A). We attribute this mass dependence to the system size effect as discussed in Section IV. In the improved MFM formulation, therefore, this system size effect is empirically realized by reducing the apparent temperature as A increases as T(A)=T₀(1-kA). T₀ is the temperature of the fragmenting source and k is a constant quantifying the mass dependence.

In order to study the density, temperature and symmetry energy in the fragmenting source, the improved MFM of Eq.(1) is utilized to extract the a_sym/T₀ value from the available isotope yields. Since the a_sym/T₀ value in Eqs.(1) and (2) depends on 5 parameters, av, as, ac, ap and Δμ (defined by Δμ = μn-μp), the optimization process of these parameters is divided into the following three steps to minimize the ambiguity for each parameter. For a given k value

1. Optimize Δμ/T₀ and ac/T₀ values from mirror isobars and fix these parameter values.

2. Optimize av/T₀, as/T₀ and ap/T₀ values from N=Z isotopes.

3. Using extracted parameters in step (1) and step (2), a_sym/T₀ values are extracted from all available isotopes. Comparing the extracted a_sym/T₀ values with three different interactions, the density of the fragmenting source is extracted. Using this density, the value of the symmetry energy coefficient, a_sym, for each interaction is determined. The temperature is then calculated as the ratio of a_sym to a_sym/T₀.

It is expectable that if the k value is properly selected which means the mass dependence is well considered, a constant T₀ is obtained. Since the k value is small as seen below, we perform the optimization of the parameter k in an iterative manner in the following analysis, that is, in the first round k=k₁=0 is set in T(A)=T₀(1-kA) and calculate the temperature as a function of A, using steps (1)–(3). From this plot a new k value, $k=k_1\text{'}$, is extracted from the slope. In the second round, $k=k_2=k_1+\frac{1}{2}{k}_1\text{'}$ is used for the steps (1)–(3) and a new k value, $k=k_2\text{'}$, is extracted. If $k_2\text{'}$ is 0 within a given error range, the iteration stops and the k₂ value is fixed as the mass dependent parameter of the apparent temperature and the T₀ value is determined. Otherwise the iteration continues.

The details of steps (1)–(3) are first described below for a given k value. In the step (1), following Ref. [19], the isotope yield ratio between isobars with I+2 and I, R(I+2,I,A), is utilized, which is

where Y(I,A) is the yield of isotopes with I and A, and Δ(I+2,I,A)=S_mix(I+2,A) - S_mix(I,A). When the above equation is applied for a pair of mirror nuclei of odd mass isotopes with I = -I and I, the symmetry energy term, pairing term and mixing entropy terms drop out and the following equation is obtained.

For all available mirror isobars, Δμ/T₀ and ac/T₀ are optimized in Eq.(4). The ln[R(I, -I, A)]/I values and the fit result for k=0 is shown in Fig.1. In the step (2) we apply Eq.(1) to the isotopes with N=Z with the extracted Δμ/T₀ and ac/T₀ values in the step (1). For the N = Z = A/2 isotopes, the ratio of the free energy relative to the temperature can be calculated from Eq.(1) and Eq.(2) without the symmetry energy term as

Fig. 1.Full-screen View

(Color online) ln[R(I,-I,A)]/I versus A for I=1 and I=3 from the events with the g0 interaction. The curve is the fit result of Eq.(4) for k=0. The extracted values of Δμ/T₀ and ac/T₀ are given in the third and fifth columns of Table 1.

where $\widetilde{a_v}=a_v+\frac{1}{2}({\mu }_n+{\mu }_p)$. The value of $ln[\frac{Y(A/\mathrm{2,}A/2)A^{\tau }}{Y_0}]$ on the right of the second equation can be calculated from the simulated or experimental values when the τ value is fixed. Non-zero values show the deviation of the mass distribution of N=Z isotopes from the power law distribution determined by the critical exponent [40]. When other τ values are used, the parameter values change accordingly. In order to eliminate the constant Y₀, all isotope yields are normalized to the yield of ¹²C [40, 18, 19]. For the first round with k=0, the renormalized values of $-\frac{F(A/\mathrm{2,}A/2)}{T_0}$ from the AMD events with the g0 interaction are plotted as a function of the isotope mass A using solid points in Fig. 2(a). The values of $\widetilde{a_v}/T_0$, as/T₀ and ap/T₀ are used as free parameters to fit the given $-\frac{F(A/\mathrm{2,}A/2)}{T_0}$ values, employing Eq.(5). A typical search result is shown by open circles in Fig. 2(a) for the case of the g0 interaction. Similar quality results are obtained for the events generated using the g0AS and g0ASS interactions. One should note that the value of ap/T₀ makes a small contribution and that the contribution is evident as staggering in the -F(A/2,A/2)/T₀ vs. A plot. Therefore the essential free parameters in this step are $\widetilde{a_v}/T_0$ and as/T₀. The extracted parameter values are summarized in Table 1 for the first round (k=0) and the final round (k=0.007).

Fig. 2.Full-screen View

(Color online) (a) Calculated ratio of free energy to T₀ for N=Z isotopes from the AMD events with g0 (solid points). Circles represent the fit using Eq.(5). (b) Calculated $\frac{\Delta F(N\mathrm{,}Z)}{T_0}$ values for the g0 interaction and quadratic fits using Eq.(6) for Z=2 to 18. The same symbols are used for isotopes with a given Z. (c) Extracted a_sym/T₀ values from (b) for g0 (dots), g0AS (squares) and g0ASS (triangles). Errors are from the quadratic fits. All results are from the first round (k=0).

In step (3), Eq.(1) is applied to the yields of all isotopes with N=Z and N Z. From Eq.(1), a_sym/T₀ and Δμ/T₀ values can be related to the modified free energy, $\frac{\Delta F(N\mathrm{,}Z)}{T_0}$ as

where $\frac{\Delta F(N\mathrm{,}Z)}{T_0}$ is the free energy relative to the temperature, $\frac{F(N\mathrm{,}Z)}{T_0}$, subtracted by the calculated contributions of the volume, surface, Coulomb and pairing terms, using the parameters in Table 1. Resultant $\frac{\Delta F(N\mathrm{,}Z)}{T_0}$ values are shown by symbols in Fig. 2(b). They exhibit quadratic relationships with minimum values close to zero. The minimum values are at or near N=Z isotopes and therefore reflect approximately the difference between the data and fitting points in Fig. 2(a). In this step, the a_sym/T₀ and the Δμ/T₀ values are optimized. Since the Δμ/T₀ values are extracted from the step(1), the optimization is made around the values in the fifth column of Table 1 in a small margin. The a_sym/T₀ values are extracted from the quadratic curvature of the isotope distribution for each given Z and plotted in Fig. 2(c) separately for the g0, g0AS and g0ASS interactions. As one can see for the first round with k=0, the extracted a_sym/T₀ values increase as Z increases in all cases, and they more or less parallel each other.

III. SELF-CONSISTENT DETERMINATION OF DENSITY AND TEMPERATURE

In order to determine the density and temperature at the time of the fragment formation, the parallel behavior of the observed a_sym/T₀ values in Fig. 2(c) is utilized. As suggested in Ref. [12], the observed differences are attributed to the difference of the symmetry energy at the density at the time of the fragment formation. The ratios between g0/g0AS and g0/g0ASS of a_sym/T₀ for the first round are shown in Fig. 3(a). The ratios show flat distributions as a function of Z for both cases. The extracted average ratio values are shown by lines in the figure and the values are given in the first column of Table 2. In Fig. 3(b) the symmetry energy coefficient is plotted as a function of the density for the three interactions used in the calculations and in Fig. 3(c) their ratios, R_sym = a_sym(g0)/a_sym(g0AS) and R_sym = a_sym(g0)/a_sym(g0ASS), are plotted. Using the ratio values determined from Fig. 3(a) and the density dependence of the R_sym values in Fig. 3(c), the implied densities of the fragmenting sources are indicated by the shaded vertical areas shown in Fig. 3(c). The extracted density values for each case are given in the second column of Table 2. Assuming that the nucleon density should be same for the three different interactions used, the nucleon density of the fragmenting source is determined from the overlap of the extracted values. This assumption is reasonable for the central collisions because the nucleon density is mainly determined by the stiffness of the EOS and not by the density dependence of the symmetry energy term. From the overlapped density area in Figs. 3(c), ρ/ρ₀ = 0.67 ± 0.02 is extracted as the density at the time of the fragment formation. Using this density value, the corresponding symmetry energy values at that density are extracted for the three different interactions from Fig. 3(b). They are given in the third column of Table 2.

Fig. 3.Full-screen View

(Color online) (a) The ratios of the a_sym/T₀ values shown in Fig. 2(c), dots for g0/g0AS and circles for g0/g0ASS. (b) Symmetry energy coefficient vs density used in the AMD simulations. In the plot, solid, dashed and dotted lines represent the g0, g0AS and g0ASS interactions, respectively. (c) The ratio of the symmetry energy coefficient in (b). The shaded horizontal lines indicate the ratios extracted in (a) and the vertical shaded area shows the density region corresponding to these ratios. Two different shadings are used for the two ratio values. All results are from the first round.

Once the symmetry energy value is determined for a given interaction, the temperature, T₀, can be calculated as T₀ = a_sym / (a_sym/T₀). The extracted T₀ values are shown as a function of Z by open symbols for the first round in Fig. 4. The larger errors of T₀, comparing to those in Fig. 2(c), originate from the errors of a_sym and a_sym/T₀ extracted for each interaction which are shown in the third column of Table 2 and Fig. 2(c), respectively. The temperature values extracted from the three different interactions agree with each other very well and show a monotonic decrease as Z increases from 5 MeV at Z = 4 to 3 MeV at Z = 18. From this slope the extracted temperature as a function of A, T₀ = 5.5 (1 - 0.012A), is determined for the first round, assuming A 2Z.

Fig. 4.Full-screen View

(Color online) Extracted T₀ values as a function of Z. Open symbols are for the first round (k=0) and closed symbols are for the final round (k=0.007). Different symbols represent results with g0(circles), g0AS(squares) and g0ASS(triangles) interactions.

The iteration is repeated four times in this work. The same plots as Fig.2, but with the k value for the fourth (final) round, k=0.007, are shown in Fig.5 and the extracted parameters are also given in Table 1. A very similar quality of results to those of the first round with k=0 were obtained, even though the optimized parameter values are quite different between those of the first round (k=0) and of the fourth round (k=0.007). The extracted a_sym/T₀ values parallel each other and show a rather flat distribution as a function of Z for Z up to 15 in Fig. 5(c). As seen in Fig. 4, in which the extracted T₀ values are shown by closed symbols as a function of Z for the fourth round with k=0.007, the extracted T₀ values are consistent with 5.5 MeV within the error bars. Since T₀ values show a flat distribution as a function of Z, the iteration is stopped at this round and T₀=(5.5±0.2) MeV is taken as the temperature of the emitting source.

Fig. 5.Full-screen View

(Color online) Same plots as Fig.2, but for the final round (k=0.007). See also the figure caption of Fig.2.

The extracted parameter values of R_sym, density and symmetry energy for the fourth round are very similar to those of the first round as shown in Table 2. The values and errors of these parameters are essentially determined by the ratios and their errors of the a_sym/T₀ values between different interactions as discussed in Fig.3(a). These ratio values are stable between the first and fourth rounds, even though the optimized parameter values in Table 1 are quite different between these rounds. The extracted value of the mass dependent factor k has some errors, but the above fact ensures that the extracted density, temperature and symmetry values and their errors in Table 2 are rather stable independent of the choice of the k values within its error bar, when the parameters in Table 1 are optimized for the given k value.

The decreasing trend of the temperature as A increases is often observed in heavy ion collisions and normally attributed to variations in impact parameter [47]. The heavier fragments tend to be produced in more peripheral collisions and therefore show a lower temperature. However in this study all events analyzed are generated in the same class of events, central collisions with b=0 fm. Therefore we attribute the decreasing trend may originate to different fragment formation processes rather than to the centrality of the events.

IV. SUMMARY

An improved method is proposed for the extraction of the ratio of the symmetry energy coefficient relative to the temperature, a_sym/T(A), taking into account the mass dependence of the apparent temperature, based on the MFM model. This method is applied for the central collisions of the AMD events generated for ⁴⁰Ca + ⁴⁰Ca at 35 MeV/nucleon. The Gogny interactions, g0, g0AS and g0ASS, with three different density dependencies of the symmetry energy are employed. As a function of IMF charge Z, the ratios of the extracted a_sym/T₀ values from different interactions are essentially constant and reflect the differences of the symmetry energy at that density at the time of the fragment formation. Using this correlation, ρ/ρ₀ = 0.67 ± 0.02 is evaluated as the density at the time of fragmenting source and the symmetry energy value at that density are extracted for each interaction. The temperature values are then determined as T₀=5.5 MeV. The apparent temperatures show a monotonic decrease as the fragment Z increases, changing from 5 MeV to 3 MeV when Z increases from 4 to 18 due to the system size effect.