Enhanced dilepton emission from a phase transition in dense matter

Baryon interaction dynamics in our approach consist of two contributions: (i) consecutive two-body scatterings and (ii) density-dependent potential. The latter does include finite-range many-body interactions. While we do not incorporate explicit three-body scatterings, to our knowledge, it has never been shown unambiguously that their introduction leads to significant changes in the dynamics [38, 39].

The treatment of the density-dependent potential in UrQMD as well as the construction and motivation for the phase transitions used in the following has been described in more detail in [40]. In the following, we will simply give a short summary and refer the interested reader to the above paper. As a baseline for the density-dependent EoS we use the CMF model in its most recent version [41]. The CMF model incorporates a realistic description of nuclear matter with a nuclear incompressibility of K₀ = 267 MeV, chiral symmetry breaking in the hadronic and quark sectors, as well as an effective deconfinement transition. Recently it was shown how the CMF EoS can be effectively implemented into the QMD equations of motion of the UrQMD model [34]:

$$\begin{eqnarray}\begin{array}{rcl}{\dot{{\bf{p}}}}_{i} & = & -\displaystyle \frac{\partial {\rm{H}}}{\partial {{\bf{r}}}_{i}}=-\displaystyle \frac{\partial {\rm{V}}}{\partial {{\bf{r}}}_{i}}\\ & = & -\left(\displaystyle \frac{\partial {V}_{i}}{\partial {n}_{i}}\cdot \displaystyle \frac{\partial {n}_{i}}{\partial {{\bf{r}}}_{i}}\right)-\left(\displaystyle \sum _{j\ne i}\displaystyle \frac{\partial {V}_{j}}{\partial {n}_{j}}\cdot \displaystyle \frac{\partial {n}_{j}}{\partial {{\bf{r}}}_{i}}\right)\,,\end{array}\end{eqnarray} \tag{ 1 }$$

where H is the total Hamiltonian and V = ∑_j V[n_B (r_j )] is the total potential energy of the system ⁹ . The potential derivative in equation (1) is calculated by summing over all possible baryon pairs j ≠ i. V_i corresponds to the potential energy of a baryon at position r_i , and the local interaction density n_k at position r_k is calculated assuming that each particle can be treated as a Gaussian wave packet [30, 42].

In addition to the derivative of the mean-field potential energy per baryon ∂V(n_B )/∂n_B that is used as input to UrQMD, only the local densities and gradients need to be calculated to solve the equations of motion for each baryon. The baseline CMF mean-field energy per baryon V(n_B ) as a function of the baryon density is shown as an orange line in the upper panel of figure 1. The minimum around saturation density corresponds to the nuclear liquid–gas transition; at higher density the potential monotonously increases due to the presence of repulsive interactions and the absence of any other phase transitions.

Zoom In Zoom Out Reset image size

A simple augmentation is used in order to implement a phase transition in the CMF model [40]. To provide for another metastable state in the mean-field energy per baryon at large densities (in addition to the bound state from the nuclear liquid–gas transition), the original potential of CMF is cut at density ${n}_{B}^{{cut}}=2.1\,{n}_{0}$ for PT1 and 2.6 n₀ for PT2. V(n_B ) for ${n}_{B}\gt {n}_{B}^{{cut}}$ is then shifted by Δn_B = 1.7 n₀ and 2.6 n₀, respectively.

The specific choice of the value for n_cut is motivated by two factors:

• 1.  
  The transition should be reachable with HIC experiments (which limits the density range to ${n}_{B}^{{cut}}\lt 4-5{n}_{0}$) [34].
• 2.  
  The density is higher than 2n₀ to avoid discrepancies with available constraints from HICs and astrophysical observations [43].

Since the pressure consists not only of the potential but also the kinetic part, which increases with the density, an equal potential depth essentially means a higher pressure. Therefore, the high-density state of matter will coexist with a lower-density system but become unstable as soon as the density keeps decreasing due to the expansion of the low-density matter surrounding it.

The mean-field energy between ${n}_{B}^{{cut}}\lt {n}_{B}\lt {n}_{B}^{{cut}}+{\rm{\Delta }}{n}_{B}$ is then interpolated by a third-order polynomial in order to create a second minimum in energy per particle V(n_B ) and ensure that its derivative is also a continuous function. There are infinite ways to create such a construction, each leading to different properties of the so-constructed transition. In [44], it was shown that the pressure that is obtained with such a density-dependent potential indicates instabilities beyond the spinodal line of a phase transition for which the baryon density is the order parameter. In addition, the derivatives of that pressure with respect to the chemical potential will diverge at the spinodal lines as expected for a first-order phase transition. For infinite systems we therefore expect the system to behave as it would undergo a first-order phase transition exactly. Since we are dealing with systems that are finite in space and time, one does not strictly obtain an exact phase transition, as there is always a finite probability that the system will be in either phase or a metastable state. Thus, we have checked in box simulations that indeed we obtain phase separation in clusters in the unstable region of our phase diagram as indicated in figure 1. It is therefore valid to claim that we are investigating the dynamics of interacting matter which, in the thermodynamic limit, would undergo a phase transition at the indicated densities. Other methods can be employed to construct high-density phase transitions; for example, using a vector density functional [45]. The goal of the present work is to study general signatures of a significant phase transition in the SIS100 energy regime, so we chose to add a significant dip in the potential for densities that are low enough to be reached by HICs at the SIS18/SIS100 accelerators. Since this procedure modifies the CMF EoS only at high densities, it leaves the low-density description consistent with nuclear matter properties and lattice QCD constraints. The resulting average field energy per baryon used for the CMF and PT1 as well as PT2 EoS are shown as orange, green and purple lines in panel (a) of figure 1. Here, the explicitly visible additional minima of the field energy leads to the presence of two coexisting different-density states.

Let us stress that the current investigation does not employ an explicit change of the dynamical degrees of freedom from hadrons to quarks, but incorporates the effects of the phase transition in the EoS on the dynamics of the system. Therefore, any effects that would specifically depend on the deconfinement aspects of the transition will not be included in the transport simulation; however, they can be included in the coarse-graining procedure. In case of a crossover transition in a mean-field approach, a possible solution for changing the degrees of freedom has been presented in [46].

To have a better understanding of the expected influence of the modified potential on the phase structure, phase diagrams of the two CMF-PT models in the thermodynamic limit are constructed and shown in panel (b) of figure 1. The spinodal lines [33], boundaries of regions where matter becomes mechanically unstable, are defined by the densities at which the isothermal speed of sound becomes imaginary. The pressure as well as the chemical potential at given (n_B , T) can be calculated using the single-particle energy:

$$\begin{eqnarray}&&U({n}_{B})=\displaystyle \frac{\partial ({n}_{B}\cdot V({n}_{B}))}{\partial {n}_{B}}\,,\end{eqnarray} \tag{ 2 }$$

$$\begin{eqnarray}&&\mu ({n}_{B})={\mu }_{\mathrm{id}}({n}_{B})-U({n}_{B}),\end{eqnarray} \tag{ 3 }$$

$$\begin{eqnarray}&&P({n}_{B},T)={P}_{\mathrm{id}}({n}_{B},T)+{\int }_{0}^{{n}_{B}}n^{\prime} \displaystyle \frac{\partial U(n^{\prime} )}{\partial n^{\prime} }{dn}^{\prime} .\end{eqnarray} \tag{ 4 }$$

where P_id(n_B , T) and μ_id(n_B ) are the pressure and chemical potential of an ideal hadron resonance gas. Having calculated the pressure one can define the spinodal lines, which are shown as dashed lines in figure 1. The shaded areas are the spinodal regions of the corresponding EoS. Note that the phase transitions implemented in this work are rather strong and since we did not consider any explicit temperature dependence of the potential, there is no critical endpoint present at any finite temperature in both cases.

We would like to clarify that the QMD approach solves an interacting many-body system and is not a mean-field solution. It is true that the particular density-dependent potential used in the solution of the molecular dynamics equations of motion can be related to thermodynamics quantities in the mean-field approximation, which is what we do if we calculate the pressure and/or speed of sound of the system, as shown in [34, 44], for example. The CMF model that was used to construct the density-dependent potential is also based on a mean-field approximation. However, as we have discussed, in principle any density EoS is possible.

Panel (b) of figure 1 also includes the averaged trajectories of the central cell for central collisions at different beam energies, obtained by the coarse-graining method described below. In all three cases, the nature of the nuclear liquid–gas transition is not changed and is shown as the gray shaded area. The trajectories are shown only after the point of highest compression is reached, during the compression of the colliding nuclei, which is a direct result of the stiffness of the EoS. At the lowest beam energies the maximally reached compression in the CMF and PT2 cases is similar, whereas in the PT1 case one reaches much larger densities due to the phase separation. At the highest energy, the effect of the transition for PT1 is already much less significant, since the maximum compression is beyond the spinodal region. Therefore, we would expect to see effects of the phase transition for PT1 already at the lowest beam energy, while effects for PT2 should show up for E_lab = 2 AGeV. In addition, the softening will also lead to a prolonged lifetime of the system and thus larger dilepton emission, as has been shown recently in a fluid dynamic simulation [28]. In the following, we study whether such a signal survives the bulk evolution in a fully microscopic non-equilibrium model. The effects of the bulk viscosity [47] do not change the lifetime of the system significantly, as can be seen from the extended lifetime, which only appears in the simulation with a phase transition. A recent study also shows the formation of density clumps only in the presence of a phase transition [48]. We observe an almost EoS-independent pion multiplicity, which indicates no significant additional entropy production [49].

To calculate the dilepton emission during the collision, a coarse-graining approach will be employed. This approach was used in [50–60] to extract the properties of the medium from microscopic transport simulations. For this procedure, 50 000 events of the most central (b < 2 fm) Au–Au collisions are simulated for several beam energies, using either of the three EoS introduced above. To calculate the medium properties, an ensemble average of the local energy and baryon densities is calculated on a space-time grid with Δx, Δy, Δz, Δct = 0.5 fm. The grid constitutes ${\left(100\right)}^{4}$ cells spanning −25 < x [fm] <25 in space and 0 < t [fm/c] <50 in time. This allows one to incorporate most (if not all) of the system created in Au–Au collisions at energies E_lab = 1 − 10 AGeV.

The medium properties, i.e. the baryon density and temperature, can be extracted from the CMF EoS only in the local rest frame of each gridpoint. Therefore, the local flow velocity $u=\gamma (1,\vec{v})$ has to be calculated using the ideal fluid equilibrium stress-energy tensor:

$$\begin{eqnarray}&&{T}^{\mu \nu }=(\varepsilon +P){u}^{\mu }{u}^{\nu }-{{Pg}}^{\mu \nu },\end{eqnarray} \tag{ 5 }$$

and the 0th component of baryon current:

$$\begin{eqnarray}&&{N}_{B}^{0}=\gamma {n}_{B},\end{eqnarray} \tag{ 6 }$$

where ε, n_B , P are the energy and baryon densities and pressure in the co-moving frame. Using the kinetic energy, momentum and baryon densities extracted from the coarse-grained UrQMD simulations and the pressure from the EoS, the local rest frame properties n_B and T can be found using a fixed-point iterative method used also in ideal fluid dynamical simulations [61]. ¹⁰ There are several options on how to define the local co-moving frame either as the Landau frame or the Eckart frame, which differ in second-order viscous fluid dynamics [62]. Both definitions have been used in previous implementations of the coarse-graining procedure and no significant differences have been observed [50, 51, 53–60]. This is likely due to the fact that small differences in the flow velocity do not translate into significant changes in the overall dilepton emission. In order to self-consistently extract the temperature using the CMF model in the non-ideal hydrodynamics scenario, one would not only need to include the full CMF fields (vector and scalar) in UrQMD simulation but also transport coefficients in order to properly deal with the off-diagonal elements [63]. This is something that requires a proper covariant description of the QMD simulation, and will be discussed in future works. We have estimated that using the energy-momentum tensor extracted from UrQMD, as described in the text, and matching this to the CMF EoS will give temperatures that differ only by less than 2% from the 'correct' equilibrium temperature if the full fields are taken into account.

Note that we do not employ the coarse-graining method for the early stage of the collision when the system is far away from thermal equilibrium. Though the relaxation time is predicted to be smaller than the interpenetration time [28], the emission of dileptons in the very early stages of collisions is often not considered because of this uncertainty. On the other hand, a recent study has shown that shortly after interpenetration, the systems described do reach a state of local isotropy and equilibrium, and the assumption of a near-hydro T^{μ ν} is justified [64]. In addition to extracting the local density and temperature, the coarse-graining procedure allows the local densities of π and ρ mesons and Δ baryons to be determined, which can be used to check whether a system is in local chemical equilibrium.

It was found that the pion multiplicity n_π quickly exceeds the statistical model prediction of an ideal pion gas at a given temperature ${n}_{\pi }^{\mathrm{id}}(T)$ by a factor

$$\begin{eqnarray}&&{\lambda }_{\pi }=\displaystyle \frac{{n}_{\pi }}{{n}_{\pi }^{\mathrm{id}}(T)}\,,\end{eqnarray} \tag{ 7 }$$

which can be referred to as pion fugacity. To take into account the pion number excess, the values of λ_π are used in calculations of the dilepton emissivity.

Once the system's space-time evolution is known from the coarse-graining, the emitted dileptons and their spectra as a function of the pair invariant mass M and pair momentum k can be inferred from the emission rate:

$$\begin{eqnarray}\begin{array}{rcll}\displaystyle \frac{{d}^{8}N}{{d}^{4}{{xd}}^{4}k} & = & -{\lambda }_{\pi }^{1.3}\displaystyle \frac{{\alpha }^{2}}{{\pi }^{3}{M}^{2}}{f}^{{BE}}({k}^{0},T)\displaystyle \frac{1}{3}{g}^{\mu \nu } & \mathrm{Im}\left[{{\rm{\Pi }}}_{{EM}}^{\mu \nu }(M,k,{n}_{B},T)\right]\,,\end{array}\end{eqnarray} \tag{ 8 }$$

where α is the fine structure constant, and the power of 1.3 [55, 57] for the pion fugacity comes from the fact that not only pions and the ρ meson, but also direct decays from other resonances that contain one or more pions will contribute to the dilepton emissivity. The fugacity effects on the ρ and ω mesons are assumed to be similar. Additional investigations are required to determine the precise values of the exponents. However, our main conclusions will involve parts of the spectra outside of the two peaks. The in-medium properties enter through the spectral function ${{\rm{\Pi }}}_{{EM}}^{\mu \nu }(M,k,{n}_{B},T)$, where n_B and T are the co-moving frame baryon density and temperature.

In the present work, we use the most recent version of the Rapp–Wambach–van Hees spectral function as introduced in [65–67]. This spectral function describes the EM properties of hot and dense hadronic matter calculated from the in-medium ρ, ω and ϕ meson spectral functions from hadronic many-body theory [22], supplemented by a four-pion continuum with chiral mixing at masses above 1 GeV. We assume that most other sources of dilepton emission, e.g. late-stage hadronic contributions, like the final decays of the π⁰, η, ω and ϕ mesons (often referred to as 'cocktail' contributions) as well as very early time emission, can be subtracted from the experimental measurements as was done in the NA60 and HADES experiments, for example. Therefore, we only simulate the dilepton yield, which carries information from the hot and dense phase, i.e. the in-medium ρ, ω and ϕ channels. Figure 2 shows the spectral functions of all three included vector mesons in the vacuum (solid line) and in the medium (dashed line). The in-medium spectral function is shown at a temperature of T = 120 MeV and three times saturation density, which is representative of the systems created in the energy range we study. The ρ meson shows the strongest broadening, while ω and ϕ are also broadened in the medium. Interestingly, the in-medium ω appears similar to a vacuum ρ, which means one can easily misinterpret the small peak in the measured dilepton spectrum for a ρ, while it in fact originates from the ω.

Zoom In Zoom Out Reset image size

It should also be noted that this in-medium rate becomes more similar to a pure quark rate as the temperature or density is increased, which is a consequence of the quark–hadron duality near the deconfinement transition. This also means that the spectral function is in principle also valid for densities above the deconfinement transition, and thus allows us to make realistic predictions on the effect of the phase transition. Further explicit effects on the spectral function, e.g. from a chiral critical point [68] or diquark fluctuations [69], are not part of this work.

To understand and interpret the resulting dilepton spectra and yields, it is instructive to look first at the regions in the phase diagram where the dileptons are emitted in the three different scenarios. Figure 3 shows the normalized emission rate as a function of the baryon density and temperature for the CMF, as well as PT1 and PT2 EoS (from left to right). For each EoS, three separate beam energies are presented: E_lab = 1.23, 4 and 10 AGeV, from bottom to top. One can clearly observe that once the system undergoes the phase transition and/or softening of the EoS, the emission at high baryon densities increases. This occurs for PT1 already at 1.23 AGeV and for PT2 at 4 AGeV. Therefore, for the case of a phase transition, we would expect an increased dilepton emission and possible increase of the extracted dilepton temperature at those beam energies.

Zoom In Zoom Out Reset image size