Density fluctuations in baryon-rich quark matter

Density fluctuations in baryon-rich quark matter 1 corresp first-author Corresponding author, ko@comp.tamu.edu Corresponding author, ko@comp.tamu.edu ko@comp.tamu.edu 1 Cyclotron Institute and Department of Physics and Astronomy, Texas A&M University, College Station, TX 77843-3366, USA Cyclotron Institute and Department of Physics and Astronomy, Texas A&M University, College Station, TX 77843-3366, USA At finite baryon chemical potential, the density of a quark matter develops large fluctuations when it undergoes a first-order phase transition. Based on the transport equation derived from the Nambu-Jona-Lasinio (NJL) model, we have studied the density fluctuations in a baryon-rich quark matter that is confined in a finite volume. Allowing the expansion of the quark matter using initial conditions from either a blast wave model or a multiphase transport (AMPT) model, we have further studied the survivability of the density fluctuations as the density and temperature of the quark matter decrease. Possible experimental signatures of the density fluctuations are suggested. NJL model Baryon-rich quark matter Density fluctuations Heavy Ion collisions 1 Introduction 1 Studying the properties of baryon-rich quark-gluon plasma (QGP) is the main focus of the beam energy scan (BES) experiments [ 1 1 - 3 3 1-3 ] at the Relativistic Heavy Ion Collider (RHIC) and the future Facility for Antiproton and Ion Research (FAIR). These experiments are expected to shed light on whether the phase transition from the baryon-rich QGP to the hadronic matter is a first-order one and the location of the critical end point in the QCD phase diagram if the phase transition is first-order. To help understand what could happen in a baryon-rich QGP, we have recently used the Polyakov-Nambu-Jona-Lasinia (PNJL) model [ 4 4 ] to study its spinodal instability [ 5 5 ]. We have found via the linear response theory that the spinodal boundary in the temperature and density plane of the QCD phase diagram shrinks with increasing wave number of the unstable mode. In the small wave number or long wavelength limit, the spinodal boundary coincides with that determined from the isothermal spinodal instability in the thermodynamic approach. We have further found that the quark vector interaction suppresses unstable modes of all wave numbers. For the wave number dependence of the growth rate of unstable modes, it initially increases with the wave number but decreases when the wave number is large. Moreover, we have investigated how unstable modes would grow if one goes beyond the linear response or small amplitude limit by using the transport equation derived from the NJL model to study the time evolution of density fluctuations in a confined as well as in an expanding quark matter [ 6 6 ]. We have found that the unstable modes have significant effects on the time evolution of higher-order density moments in the quark matter, the distribution of quark number in a sub-volume of the quark matter, the quark momentum anisotropy, and dilepton production rate from quark-antiquark annihilation. We have thus suggested to use these observables as signatures for a first-order phase transition of the baryon-rich quark matter produced in heavy ion collisions. In the present paper, we review the results from the transport model study of baryon-rich quark matter. The paper is organized as follows: We first briefly review in Sect. 2 the NJL model and in Sect. 3 the transport equation that is derived from the NJL model. The transport equation is then used in Sect. 4 to study the density fluctuations in a baryon-rich quark matter that is confined in a periodic box and in Sect. 5 to study how density fluctuations are affected by the expansion of the system as in heavy ion collisions. Finally, a summary is given in Sect. 6. 2 The NJL model 1 The NJL Lagrangian for three quark flavors has the form [ 7 7 ]: http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=10344269&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=10344266&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=10344269&type= where ψ =( ψ u , ψ d , ψ s ) T , M =diag( m u , m d , m s ), and λ a is the Gell-Mann matrices with λ 0 being the identity matrix multiplied by <math id="M1"><mrow><msqrt><mrow><mn>2</mn><mo>/</mo><mn>3</mn></mrow></msqrt></mrow></math> . The last term is the Kobayashi-Maskawa-t’Hooft (KMT) interaction that breaks U (1) A symmetry [ 8 8 ] with <math id="M2"><mrow><mtext>d</mtext><mi>e</mi><msub><mi>t</mi><mi>f</mi></msub><mo stretchy="false">(</mo><mover accent="true"><mi>ψ</mi><mo>¯</mo></mover><mi>Γ</mi><mi>ψ</mi><mo stretchy="false">)</mo><mo>=</mo><mstyle displaystyle="true"><msub><mo>∑</mo><mrow><mi>i</mi><mn>,</mn><mi>j</mi><mn>,</mn><mi>k</mi></mrow></msub></mstyle><msub><mi>ε</mi><mrow><mi>i</mi><mi>j</mi><mi>k</mi></mrow></msub><mo stretchy="false">(</mo><mover accent="true"><mi>u</mi><mo>¯</mo></mover><mi>Γ</mi><msub><mi>q</mi><mi>i</mi></msub><mo stretchy="false">)</mo><mo stretchy="false">(</mo><mover accent="true"><mi>d</mi><mo>¯</mo></mover><mi>Γ</mi><msub><mi>q</mi><mi>j</mi></msub><mo stretchy="false">)</mo><mo stretchy="false">(</mo><mover accent="true"><mi>s</mi><mo>¯</mo></mover><mi>Γ</mi><msub><mi>q</mi><mi>k</mi></msub><mo stretchy="false">)</mo></mrow></math> denoting the determinant in flavor space [ 9 9 ]. It gives rise to four-point interactions in two flavors and six-point interactions in three flavors. In the two flavor case, the sum of scalar and pseudo-scalar interactions and the KMT interaction with K =- G S reduces to the original NJL model [ 10 10 , 11 11 ]. We note that for sufficient large G V , such as G V = G S , the first-order phase transition induced by the attractive scalar interaction would disappear [ 5 5 , 12 12 ]. 3 Transport model based on the NJL model 1 For the long-time behavior of a quark matter or an expanding quark matter produced in a heavy ion collision, it can be studied using the Boltzmann or transport equation derived from the NJL Lagrangian for the phase-space distribution functions f a ( X , p) of quarks and antiquarks of flavor a [ 13 13 ], that is [ 14 14 ], <math display="block" id="M3"><mtable columnalign="left"><mtr><mtd><msub><mo>∂</mo><mrow><msup><mi>X</mi><mn>0</mn></msup></mrow></msub><msub><mi>f</mi><mi>a</mi></msub><mo stretchy="false">(</mo><mi>X</mi><mn>,</mn><mstyle mathvariant="bold" mathsize="normal"><mi>p</mi></mstyle><mo stretchy="false">)</mo><mo>+</mo><mfrac><mrow><msup><mi>p</mi><mrow><mi>i</mi><mo>±</mo></mrow></msup></mrow><mrow><msub><mi>E</mi><mrow><msup><mstyle mathvariant="bold" mathsize="normal"><mi>p</mi></mstyle><mo>±</mo></msup></mrow></msub></mrow></mfrac><msub><mo>∂</mo><mrow><msup><mi>X</mi><mi>i</mi></msup></mrow></msub><msub><mi>f</mi><mi>a</mi></msub><mo stretchy="false">(</mo><mi>X</mi><mn>,</mn><mstyle mathvariant="bold" mathsize="normal"><mi>p</mi></mstyle><mo stretchy="false">)</mo></mtd></mtr><mtr><mtd><mo>−</mo><msub><mo>∂</mo><mrow><msup><mi>X</mi><mi>i</mi></msup></mrow></msub><msubsup><mi>V</mi><mi>a</mi><mi>S</mi></msubsup><mo stretchy="false">(</mo><mi>X</mi><mo stretchy="false">)</mo><mfrac><mrow><msub><mi>M</mi><mi>a</mi></msub></mrow><mrow><msub><mi>E</mi><mrow><msup><mstyle mathvariant="bold" mathsize="normal"><mi>p</mi></mstyle><mo>±</mo></msup></mrow></msub></mrow></mfrac><msub><mo>∂</mo><mrow><msub><mi>p</mi><mi>i</mi></msub></mrow></msub><msub><mi>f</mi><mi>a</mi></msub><mo stretchy="false">(</mo><mi>X</mi><mn>,</mn><mstyle mathvariant="bold" mathsize="normal"><mi>p</mi></mstyle><mo stretchy="false">)</mo><mo>∓</mo><msub><mo>∂</mo><mrow><msup><mi>X</mi><mi>i</mi></msup></mrow></msub><msubsup><mi>V</mi><mn>0</mn><mi>V</mi></msubsup><mo stretchy="false">(</mo><mi>X</mi><mo stretchy="false">)</mo><msub><mo>∂</mo><mrow><msub><mi>p</mi><mi>i</mi></msub></mrow></msub><msub><mi>f</mi><mi>a</mi></msub><mo stretchy="false">(</mo><mi>X</mi><mn>,</mn><mstyle mathvariant="bold" mathsize="normal"><mi>p</mi></mstyle><mo stretchy="false">)</mo></mtd></mtr><mtr><mtd><mo>∓</mo><msub><mo>∂</mo><mrow><msup><mi>X</mi><mi>i</mi></msup></mrow></msub><msubsup><mi>V</mi><mi>j</mi><mi>V</mi></msubsup><mo stretchy="false">(</mo><mi>X</mi><mo stretchy="false">)</mo><mfrac><mrow><msup><mi>p</mi><mrow><mi>j</mi><mo>±</mo></mrow></msup></mrow><mrow><msub><mi>E</mi><mrow><msup><mstyle mathvariant="bold" mathsize="normal"><mi>p</mi></mstyle><mo>±</mo></msup></mrow></msub></mrow></mfrac><msub><mo>∂</mo><mrow><msub><mi>p</mi><mi>i</mi></msub></mrow></msub><msub><mi>f</mi><mi>a</mi></msub><mo stretchy="false">(</mo><mi>X</mi><mn>,</mn><mstyle mathvariant="bold" mathsize="normal"><mi>p</mi></mstyle><mo stretchy="false">)</mo><mo>=</mo><mi mathvariant="script">C</mi><mo stretchy="false">[</mo><msub><mi>f</mi><mi>a</mi></msub><mo stretchy="false">]</mo><mn>.</mn></mtd></mtr></mtable></math> In the above, p ± p±V V is the kinetic momentum with the subscript + referring to quarks and - referring to antiquarks, <math id="M4"><mrow><msubsup><mi>V</mi><mi>μ</mi><mi>V</mi></msubsup></mrow></math> is the vector potential, <math id="M5"><mrow><msub><mi>M</mi><mi>a</mi></msub><mo>=</mo><msub><mi>m</mi><mrow><mn>0</mn><mi>a</mi></mrow></msub><mo>−</mo><msubsup><mi>V</mi><mi>a</mi><mi>S</mi></msubsup></mrow></math> is the effective quark mass with <math id="M6"><mrow><msubsup><mi>V</mi><mi>a</mi><mi>S</mi></msubsup></mrow></math> being the scalar potential, and the collision term <math display="block" id="M7"><mtable columnalign="left"><mtr><mtd><mi>C</mi><mo stretchy="false">[</mo><msub><mi>f</mi><mi>a</mi></msub><mo stretchy="false">]</mo><mo>≡</mo><mstyle displaystyle="true"><munder><mo>∑</mo><mrow><mi>b</mi><mi>c</mi><mi>d</mi></mrow></munder></mstyle><mfrac><mn>1</mn><mrow><mn>1</mn><mo>+</mo><msub><mi>δ</mi><mrow><mi>a</mi><mi>b</mi></mrow></msub></mrow></mfrac><mstyle displaystyle="true"><mrow><mo>∫</mo></mrow></mstyle><mfrac><mrow><msup><mi>d</mi><mn>3</mn></msup><msub><mstyle mathvariant="bold" mathsize="normal"><mi>p</mi></mstyle><mi>b</mi></msub></mrow><mrow><msup><mrow><mo stretchy="false">(</mo><mn>2</mn><mi>π</mi><mo stretchy="false">)</mo></mrow><mn>3</mn></msup><mn>2</mn><msub><mi>E</mi><mi>b</mi></msub></mrow></mfrac><mfrac><mrow><msup><mi>d</mi><mn>3</mn></msup><msub><mstyle mathvariant="bold" mathsize="normal"><mi>p</mi></mstyle><mi>c</mi></msub></mrow><mrow><msup><mrow><mo stretchy="false">(</mo><mn>2</mn><mi>π</mi><mo stretchy="false">)</mo></mrow><mn>3</mn></msup><mn>2</mn><msub><mi>E</mi><mi>c</mi></msub></mrow></mfrac><mfrac><mrow><msup><mi>d</mi><mn>3</mn></msup><msub><mstyle mathvariant="bold" mathsize="normal"><mi>p</mi></mstyle><mi>d</mi></msub></mrow><mrow><msup><mrow><mo stretchy="false">(</mo><mn>2</mn><mi>π</mi><mo stretchy="false">)</mo></mrow><mn>3</mn></msup><mn>2</mn><msub><mi>E</mi><mi>d</mi></msub></mrow></mfrac></mtd></mtr><mtr><mtd><mo>×</mo><mfrac><mrow><msup><mrow><mo stretchy="false">(</mo><mn>2</mn><mi>π</mi><mo stretchy="false">)</mo></mrow><mn>4</mn></msup></mrow><mrow><mn>2</mn><msub><mi>E</mi><mi>a</mi></msub></mrow></mfrac><msup><mi>δ</mi><mn>4</mn></msup><mo stretchy="false">(</mo><msub><mi>p</mi><mi>a</mi></msub><mo>+</mo><msub><mi>p</mi><mi>b</mi></msub><mo>−</mo><msub><mi>p</mi><mi>c</mi></msub><mo>−</mo><msub><mi>p</mi><mi>d</mi></msub><mo stretchy="false">)</mo><mo>|</mo><msub><mi mathvariant="script">M</mi><mrow><mi>a</mi><mi>b</mi></mrow></msub><msup><mo>|</mo><mn>2</mn></msup></mtd></mtr><mtr><mtd><mo>×</mo><mrow><mo>[</mo><mrow><msub><mi>f</mi><mi>c</mi></msub><msub><mi>f</mi><mi>d</mi></msub><mo stretchy="false">(</mo><mn>1</mn><mo>−</mo><msub><mi>f</mi><mi>a</mi></msub><mo stretchy="false">)</mo><mo stretchy="false">(</mo><mn>1</mn><mo>−</mo><msub><mi>f</mi><mi>b</mi></msub><mo stretchy="false">)</mo><mo>−</mo><msub><mi>f</mi><mi>a</mi></msub><msub><mi>f</mi><mi>b</mi></msub><mo stretchy="false">(</mo><mn>1</mn><mo>−</mo><msub><mi>f</mi><mi>c</mi></msub><mo stretchy="false">)</mo><mo stretchy="false">(</mo><mn>1</mn><mo>−</mo><msub><mi>f</mi><mi>d</mi></msub><mo stretchy="false">)</mo></mrow><mo>]</mo></mrow></mtd></mtr></mtable></math> that describe scatterings among quarks and antiquarks. The transport equation can be solved using the test particle method [ 15 15 ] by expressing the distribution function in terms of the density of test particles that follow the Newton’s equations of motions determined by the left hand side of Eq.(2) and undergo scatterings according to Eq.(3) [ 14 14 , 16 16 ]. 4 Quark matter in a box 1 We first consider a quark matter that is confined in a cubic box with periodic boundary conditions. The system is prepared by distributing many test particles inside the box according to the density of the system with their momenta given by the Fermi-Dirac distribution at certain temperature. To study how density fluctuations emerge and grow, we compare results from two calculations based on the same initial conditions but with and without the spinodal instability in the equation of state, which is achieved by introducing a vector interaction in the NJL model. In both calculations, 300 test particles are used in solving the transport equation. Figure 1 Figure 1 shows the time evolution of the density distribution of a quark matter with temperature T =20 MeV and net quark density n q =0.5 fm -3 in a box of size 202020 fm 3 for the two cases of G V =0 and G V = G S , with the darker color denoting the high density regions and the lighter color denoting the low density regions. Although the system is initially uniform in space, some dense spots are present due to statistical fluctuations as a result of finite number of test particles used in the calculation. In the case of G V = G S without a first-order phase transition or spinodal instability, the density distribution in the box remains unchanged with time. This changes dramatically, however, for the case of G V =0. Due to the spinodal instability, the initial dense spots act like "seeds", which create several small low pressure regions that attract nearby partons and lead to the formation of many clusters at t =20 fm/ c . These clusters further grow in size by connecting with each other and form stable large structures at t =40 fm/ c , when the system clearly separates into two phases of matter with one of high density and the other of low density as shown in Fig. 2 Fig. 2 . The density fluctuations can be quantified by the scaled density moments ρ N / ρ N > [ 17 17 ], where ρ N d 3 r ρ (r) N+1 / d 3 r ρ (r). The scaled density moments are all equal to one for a uniform density distribution but become greater than one as the density fluctuations grow. In Fig. 3 Fig. 3 , we show by dotted, dashed, and solid lines the scaled density moments for N =2, 4 and 6, respectively, which are obtained by averaging over 1,000 events with 300 test particles in an event. Our results show that the scaled moments increase during the phase separation and reach their saturated values at about t =40 fm/ c , when the phase separation almost ends. Also, moments with larger N increase faster and saturate at larger values. Other quantities of interest are the skewness and kurtosis of the particle multiplicity distribution, which have been proposed as possible signals for the critical phenomena [ 18 18 ] and have been studied in the beam energy scan experiments at RHIC[ 1 1 , 2 2 ]. They are defined as follows: <math display="block" id="M8"><mrow><mtext>skewness</mtext><mo>≡</mo><mfrac><mrow><mo>〈</mo><mi>δ</mi><msubsup><mi>N</mi><mi>q</mi><mn>3</mn></msubsup><mo>〉</mo></mrow><mrow><msup><mrow><mo>〈</mo><mi>δ</mi><msubsup><mi>N</mi><mi>q</mi><mn>2</mn></msubsup><mo>〉</mo></mrow><mrow><mn>3</mn><mo>/</mo><mn>2</mn></mrow></msup></mrow></mfrac><mn>,</mn><mtext> </mtext><mtext>kurtosis</mtext><mo>≡</mo><mfrac><mrow><mo>〈</mo><mi>δ</mi><msubsup><mi>N</mi><mi>q</mi><mn>4</mn></msubsup><mo>〉</mo></mrow><mrow><msup><mrow><mo>〈</mo><mi>δ</mi><msubsup><mi>N</mi><mi>q</mi><mn>2</mn></msubsup><mo>〉</mo></mrow><mn>2</mn></msup></mrow></mfrac><mo>−</mo><mn>3.</mn></mrow></math> Both quantities characterize how far an event-by-event multiplicity distribution deviates from a normal distribution. Theoretical calculations based on the grand canonical picture predict that both quantities diverge with the correlation length when a system approaches its critical point [ 18 18 ], with the kurtosis diverging faster than the skewness. To be consistent with the grand canonical picture, we consider quarks in a sub-volume of the box in our study and treat the remaining part as the reservoir. When the system is initially inside the spinodal instability region, quarks in the reservoir can sometimes move into the sub-volume, but in most of the times quarks would leave from the sub-volume to the reservoir. The number of quarks inside this sub-volume thus varies drastically from event to event, leading to large values for the skewness in its event-by-event distribution. In Fig. 4 Fig. 4 , we show the event-by-event distribution of the number of quarks in a sub-volume from the 1,000 events at t =0, 20, and 40 fm/ c by the solid, dashed and dotted lines, respectively, for the two cases of sub-volume size of 0.6 fm 3 (left window) and 30 fm 3 (right window), respectively. The distribution in the case of small sub-volume clearly becomes asymmetric as time increases, starting with an initial skewness of 0.11 and increasing to 0.60 at 20 fm/ c and 0.75 at 40 fm/ c . This feature is absent in the case of large sub-volume, where the distribution remains essentially symmetric with increasing time, with the skewness changing slowly from -0.001 (t=0) to 0.086 (t=20 fm/ c ) and 0.132 (t=40 fm/ c ), and there is no apparent increase or decrease in the kurtosis. 5 Expanding quark matter 1 To study how large density fluctuations due to the spinodal instability as a result of a first-order phase transition obtained from the box calculation in the previous section are affected by the expansion of the system as in a heavy ion collision, we have carried out a dynamical calculation using the transport equation introduced in Section 3 with initial conditions taking from either a blast wave model or the AMPT model [ 19 19 ]. 5.1 Blast wave initial conditions 2 For the blast wave initial conditions, the positions of quarks and antiquarks are assumed to follow a spherical Wood-Saxon form: <math display="block" id="M9"><mrow><mi>ρ</mi><mo stretchy="false">(</mo><mi>r</mi><mo stretchy="false">)</mo><mo>=</mo><mfrac><mrow><msub><mi>ρ</mi><mn>0</mn></msub></mrow><mrow><mn>1</mn><mo>+</mo><mi>exp</mi><mo stretchy="false">(</mo><mo stretchy="false">(</mo><mi>r</mi><mo>−</mo><mi>R</mi><mo stretchy="false">)</mo><mo>/</mo><mi>a</mi><mo stretchy="false">)</mo></mrow></mfrac></mrow></math> with a radius R =5 fm and a surface thickness parameter a =0.5 fm, similar to that expected from a central Au+Au collisions. The momenta of these partons are taken to be that of a Fermi-Dirac distribution at certain temperature. The density fluctuations obtained from calculations based on 1,000 test particles in solving the transport equation can be seen from the density distribution on a plane such as the one at z =0 shown in Fig. 5 Fig. 5 . The left window shows the density distribution at t =20 fm/ c for the case with a first-order phase transition, while the right window shows that at t =10 fm/ c for the case without a first-order phase transition, when the density of the central cell is about 0.2 fm -3 in both cases. Although density clumps appear in both cases, those in the one with a first-order phase transition are significantly larger. As in the case of quark matter in a box, we can quantify the density fluctuations by the scaled density moments [ 17 17 ]. In both cases with and without a first-order phase transition, the scaled density moments first increase and then decrease with time. In the case without a first-order phase transition, this is caused by the fast increase of the surface of the quark matter and the quick deviation from its initial smooth Wood-Saxon density distribution. To the contrary, the scaled density moments in the case with a first-order phase transition becomes much larger with time and only decreases slightly afterwards, reflecting the effect due to density clumps that distribute randomly inside the expanding quark matter. Therefore, the saturated scaled density moments, which are larger for larger N , can still be regarded as signals for a first-order phase transition in an expanding baryon-rich quark matter. Since density fluctuations can lead to spatial anisotropy even in central heavy ion collisions, it has been suggested that they may affect the anisotropic flows in the transverse plane [ 20 20 , 21 21 ]. The latter are defined by the coefficients v n in the expansion of the transverse momentum distribution f ( p T , ϕ ) as a Fourier series in the azimuthal angle ϕ , <math display="block" id="M10"><mrow><mi>f</mi><mo stretchy="false">(</mo><msub><mi>p</mi><mi>T</mi></msub><mn>,</mn><mi>ϕ</mi><mo stretchy="false">)</mo><mo>=</mo><mfrac><mrow><mi>N</mi><mo stretchy="false">(</mo><msub><mi>p</mi><mi>T</mi></msub><mo stretchy="false">)</mo></mrow><mrow><mn>2</mn><mi>π</mi></mrow></mfrac><mo>{</mo><mn>1</mn><mo>+</mo><mn>2</mn><mstyle displaystyle="true"><munderover><mo>∑</mo><mrow><mi>n</mi><mo>=</mo><mn>1</mn></mrow><mi>∞</mi></munderover></mstyle><msub><mi>v</mi><mi>n</mi></msub><mo stretchy="false">(</mo><msub><mi>p</mi><mi>T</mi></msub><mo stretchy="false">)</mo><mi>cos</mi><mo stretchy="false">[</mo><mi>n</mi><mo stretchy="false">(</mo><mi>ϕ</mi><mo>−</mo><msub><mi>ψ</mi><mi>n</mi></msub><mo stretchy="false">)</mo><mo stretchy="false">]</mo><mo>}</mo><mn>,</mn></mrow></math> where ψ n is the event plane angle [ 22 22 ]. To calculate the anisotropic flow coefficients, we use the two particle cumulant method [ 23 23 , 24 24 ], namely, <math id="M11"><mrow><msub><mi>v</mi><mi>n</mi></msub><mo>{</mo><mn>2</mn><mo>}</mo><mo>=</mo><msqrt><mrow><mo>〈</mo><mi>cos</mi><mo stretchy="false">(</mo><mi>n</mi><mi>Δ</mi><mi>ϕ</mi><mo stretchy="false">)</mo><mo>〉</mo></mrow></msqrt></mrow></math> by averaging over all particle pairs in an event. We have calculated the elliptic flow ( v 2 {2}) and quadrangular flow ( v 4 {2}) of an expanding quark matter with the same blast wave initial conditions for 100 events, and their final event distributions are shown in the left and right windows of Fig. 6 Fig. 6 , respectively, with the solid and dashed lines for the cases with and without first order phase transition, respectively. Both distributions peak at a larger value for the case with a first-order phase transition, particularly for v 4 , thus providing a plausible signal for the first-oder phase transition. However, the values of the fluctuation induced v 2 and v 4 are much smaller than those in non-central heavy ion collisions. We have also studied the effect of density fluctuations on dilepton production from a quark matter. Since the dilepton production rate is proportional to the square of parton density, more dileptons are produced when the density fluctuation is large. Also, a longer partonic phase as a result of a first-order phase transition would increase the depletion yield as well. As usually used in studying dilepton production in heavy ion collisions [ 25 25 ], we adopt the perturbative approach to calculate the dilepton yield from the quark-antiquark scattering by neglecting its effect on the dynamics of the expanding quark matter, and the results from the average over 100 events are shown in Fig. 7 Fig. 7 by the solid and dashed lines for the cases with and without first-order phase transition, respectively. As expected, more dileptions are produced from the quark matter with a first-order phase transition. We note the dilepton invariant mass spectrum peaks at <math id="M12"><mrow><msqrt><mi>s</mi></msqrt><mo>≈</mo><mn>0.5</mn></mrow></math> GeV with the peak value being about 3.5× 10 -4 GeV -1 , which is comparable with the result obtained from a hadronic transport model [ 26 26 ]. This enhancement in dilepton production may thus be detectable in experiments. 5.2 AMPT initial conditions 2 For a more realistic initial conditions, we use the quark and antiquark distributions from the AMPT model with string melting [ 19 19 ] that uses the heavy ion jet interaction generator (HIJING) [ 27 27 ] as the input. This model includes not only the mini-jet partons from initial hard collisions but also hadrons produced from excited strings, which are projectile and target nucleons that have suffered interactions, by converting them to partons according to the flavor and spin structures of their valence quarks. In particular, a meson is converted to a quark and an anti-quark, while a baryon is first converted to a quark and a diquark, and the diquark is then decomposed into two quarks. The quark masses are taken to be m u =5.6, m d =9.9,ss and m s =199 MeV/ c 2 as in the PYTHIA program [ 28 28 ]. The above two-body decomposition is isotropic in the rest frame of the parent hadron or diquark. These partons are produced after a formation time of <math id="M14"><mrow><msub><mi>t</mi><mi>f</mi></msub><mo>=</mo><msub><mi>E</mi><mi>H</mi></msub><mo>/</mo><msubsup><mi>m</mi><mrow><mi>T</mi><mn>,</mn><mi>H</mi></mrow><mn>2</mn></msubsup></mrow></math> , with E H and m T , H denoting, respectively, the energy and transverse mass of the parent hadron. We obtain these partons as the initial conditions for our study of an expanding quark matter by running the AMPT program with vanishing parton scattering cross sections in Zhang’s parton cascade (ZPC) [ 29 29 ] and with the hadronic afterburner based on a relativistic transport (ART) [ 30 30 , 31 31 ] turned off. Using the partons from Au+Au collisions at zero impact parameter and a center-of-mass energy <math id="M15"><mrow><msqrt><mrow><msub><mi>s</mi><mrow><mi>N</mi><mi>N</mi></mrow></msub></mrow></msqrt><mo>=</mo><mn>2.5</mn></mrow></math> GeV as the initial distribution, we have found that some parts of the system go through the spinodal region when the SU(3) NJL model with G V =0 is used in the Boltzmann equation and in constructing the phase diagram. The results obtained with 1,000 test particles are shown in Fig. 8 Fig. 8 . It is seen that the quark matter is initially largely confined in a thin disk of thickness less than 0.5 fm. When it is allowed to free streaming without any interactions (upper row), there appear two high density clumps that fly apart in the opposite directions. This feature becomes less prominent after the inclusion of quark scattering and mean-field potentials but without a phase transition in the quark matter, i.e., taking G V = G S , as shown in the middle row of Fig. 8 Fig. 8 . With a first-order phase transition in the quark matter by setting G V =0, the lower row of Fig. 8 Fig. 8 shows that the initial central disk evolves into three disks of dense matter with one in the middle due to the strong attractions that keep some partons from moving away, besides the two forward and backward moving disks. As the quark matter expands, these disks transform into rings and finally turn into disjointed clumps. Also, the quark matter with a first-order phase transition expands twice as slow as that without a first-order phase transition. Because of the non-trivial spatial distribution even in the case of free-streaming quark matter, the scaled density moments are no longer useful quantities to characterize the density fluctuations of an expanding quark matter. On the other hand, the different density variations along the beam axis affect the parton rapidity distribution. This is because partons in the middle disc, which is present only in the case with a first-order phase transition, have a small rapidity and due to the attractive quark interactions, they attract partons from the other two discs and slow down their expansion in the longitudinal direction, thus restricting their rapidities to a narrow region around the midrapidity as shown by the solid line in Fig. 9 Fig. 9 , which is indeed much narrower than that in the case without a first-order phase transition, shown by the dashed line. This effect can thus be regarded as a possible signal of a first-order phase transition and is worth studying in experiments. We have also studied the dilepton invariant mass spectrum from an expanding quark matter with initial conditions from the AMPT model. Results obtained from averaging over 10 events with 1,000 test particles in each event is shown in Fig. 10 Fig. 10 by the solid and dashed lines for the cases with and without a first-order phase transition, respectively. As in the previous subsection using the blast-wave initial conditions, the presence of a first-order phase transition enhances the dilepton yield as a result of density fluctuations and a longer partonic phase. However, the dilepton yield is lower than that obtained from the calculation with the blast wave initial conditions by two orders of magnitude because there are very few antiquarks in the partonic matter produced in heavy ion collisions at such a low energy and also because we have not included the bremsstrahlung contribution to dilepton production from the quark-quark scattering. 6 Conclusion 1 The spinodal instability is a thermodynamic feature of a first-order phase transition in a many-body system. It occurs when its pressure in some parts decreases with increasing density. This can amplify the density fluctuations and lead to a phase separation in the system. We have studied this phenomenon by solving the Boltzmann equations using the test particle method. The calculations are based on the NJL model, which has been shown to give a good description of the vacuum properties of the hadrons and also predicts the existence of a first-order phase transition in baryon-rich quark matter. We have obtained some intuitive pictures on the phase separation in a quark matter that is either in a static box or undergoes expansion. For the case of a static box, we have found that the higher-order density moments of the quark matter increase and saturate at large values after phase separation, making them possible signals for the first-order phase transition. The skewness of the quark number event-by-event distribution in a small sub-volume of the quark matter is also found to increase, but this feature disappears if the sub-volume is large. As for the expanding quark matter, two cases have been studied. One is based on the blast-wave initial conditions, while the other using the AMPT initial conditions, which are disc-like as a result of the strong correlations between the parton rapidity and longitudinal coordinate. In both cases, we have found that the expansion of the quark matter is slowed down by the presence of a first-order phase transition. Density clumps are found to appear and lead to an anisotropy in the momentum space, which can be characterized by the scaled density moments and the anisotropic flows v 2 and v 4 , respectively. An enhancement in the dilepton yield is also observed. These observables can thus be used as signatures for a first-order phase transition of the baryon-rich quark matter produced in heavy ion collisions. References T. K. Nayak (STAR Collaboration) , Study of the Fluctuations of Net-charge and Net-protons Using Higher Order Moments . Nucl. Phys. A 830 , 555c - 558c ( 2009 ). doi: 10.1016/j.nuclphysa.2009.09.046 http://doi.org/10.1016/j.nuclphysa.2009.09.046 Study of the Fluctuations of Net-charge and Net-protons Using Higher Order Moments M. M. Aggarwal et al . (STAR), Higher moments of net proton multiplicity distributions at RHIC . Phys. Rev. Lett . 105 , 022302 ( 2010 ). doi: 10.1103/PhysRevLett.105.022302 http://doi.org/10.1103/PhysRevLett.105.022302 (STAR), Higher moments of net proton multiplicity distributions at RHIC D. McDonald , Overview of results from phase I of the Beam Energy Scan program at RHIC . EPJ Web Conf . 95 , 01009 ( 2015 ). doi: 10.1051/epjconf/20159501009 http://doi.org/10.1051/epjconf/20159501009 Overview of results from phase I of the Beam Energy Scan program at RHIC K. Fukushima , Phase diagrams in the three-flavor Nambu-Jona-Lasinio model with the Polyakov loop . Phys. Rev. D 77 , 114028 ( 2008 ). doi: 10.1103/PhysRevD.77.114028 http://doi.org/10.1103/PhysRevD.77.114028 Phase diagrams in the three-flavor Nambu-Jona-Lasinio model with the Polyakov loop F. Li and C. M. Ko , Spinodal instabilities of baryon-rich quark-gluon plasma in the Polyakov-Nambu-Jona-Lasinio model . Phys. Rev. C 93 , 035205 ( 2016 ). doi: 10.1103/PhysRevC.93.035205 http://doi.org/10.1103/PhysRevC.93.035205 Spinodal instabilities of baryon-rich quark-gluon plasma in the Polyakov-Nambu-Jona-Lasinio model F. Li and C. M. Ko , arXiv:1606.05012 [nucl-th] . arXiv:1606.05012 [nucl-th] N. M. Bratovic , T. Hatsuda , and W. Weise , Role of vector interaction and axial anomaly in the PNJL modeling of the QCD phase diagram . Phys. Lett. B 719 , 131 - 135 . ( 2013 ). doi: 10.1016/j.physletb.2013.01.003 http://doi.org/10.1016/j.physletb.2013.01.003 Role of vector interaction and axial anomaly in the PNJL modeling of the QCD phase diagram G. ’t Hooft , Computation of the quantum effects due to a four-dimensional pseudoparticle . Phys. Rev. D 14 , 3432 ( 1976 ). doi: 10.1103/PhysRevD.14.3432 http://doi.org/10.1103/PhysRevD.14.3432 Computation of the quantum effects due to a four-dimensional pseudoparticle M. Buballa , NJL-model analysis of dense quark matter . Phys. Rept . 407 , 205 - 376 . ( 2005 ). doi: 10.1016/j.physrep.2004.11.004 http://doi.org/10.1016/j.physrep.2004.11.004 NJL-model analysis of dense quark matter A. Masayuki and Y. Koichi , Chiral restoration at finite density and temperature . Nucl. Phys. A 504 , 668 - 684 . ( 1989 ). doi: 10.1016/0375-9474(89)90002-X http://doi.org/10.1016/0375-9474(89)90002-X Chiral restoration at finite density and temperature Y. Nambu and G. Jona-Lasinio , Dynamical model of elementary particles based on an analogy with superconductivity . Phys. Rev . 122 , 345 ( 1961 ). doi: 10.1103/PhysRev.122.345 http://doi.org/10.1103/PhysRev.122.345 Dynamical model of elementary particles based on an analogy with superconductivity J. Xu , T. Song , C. M. Ko , and F. Li , Elliptic flow splitting as a probe of the QCD phase structure at finite baryon chemical potential . Phys. Rev. Lett . 112 , 012301 ( 2014 ). doi: 10.1103/PhysRevLett.112.012301 http://doi.org/10.1103/PhysRevLett.112.012301 Elliptic flow splitting as a probe of the QCD phase structure at finite baryon chemical potential S. P. Klevansky , A. Ogura , and J. Hufner , Ann. Phys . 261 , 37 ( 1997 ). Ann. Phys C. M. Ko , T. Song , F. Li , V. Greco , and S. Plumari , Partonic mean-field effects on matter and antimatter elliptic flows . Nucl. Phys. A 928 , 234 - 246 ( 2014 ). doi: 10.1016/j.nuclphysa.2014.05.016 http://doi.org/10.1016/j.nuclphysa.2014.05.016 Partonic mean-field effects on matter and antimatter elliptic flows C. Y. Wong , Dynamics of nuclear fluid. VIII. Time-dependent Hartree-Fock approximation from a classical point of view . Phys. Rev. C 25 , 1460 ( 1982 ). doi: 10.1103/PhysRevC.25.1460 http://doi.org/10.1103/PhysRevC.25.1460 Dynamics of nuclear fluid. VIII. Time-dependent Hartree-Fock approximation from a classical point of view G. F. Bertsch and S. Das Gupta , A guide to microscopic models for intermediate energy heavy ion collisions . Phys. Rep . 160 , 189 - 233 ( 1988 ). doi: 10.1016/0370-1573(88)90170-6 http://doi.org/10.1016/0370-1573(88)90170-6 A guide to microscopic models for intermediate energy heavy ion collisions J. Steiheimer and J. Randrup , Spinodal density enhancements in simulations of relativistic nuclear collisions . Phys. Rev. C 87 , 054903 ( 2013 ). doi: 10.1103/PhysRevC.87.054903 http://doi.org/10.1103/PhysRevC.87.054903 Spinodal density enhancements in simulations of relativistic nuclear collisions M. A. Stephanov , Non-Gaussian fluctuations near the QCD critical point . Phys. Rev. Lett . 102 , 032301 ( 2009 ). doi: 10.1103/PhysRevLett.102.032301 http://doi.org/10.1103/PhysRevLett.102.032301 Non-Gaussian fluctuations near the QCD critical point Z. W. Lin , C. M. Ko , B. A. Li , et al ., Multiphase transport model for relativistic heavy ion collisions . Phys. Rev. C 72 , 064901 ( 2005 ). doi: 10.1103/PhysRevC.72.064901 http://doi.org/10.1103/PhysRevC.72.064901 Multiphase transport model for relativistic heavy ion collisions C. Herold , M. Nahrgang , and I. Mishustin , and M. Bleicher , J. Phys. Conf. Ser . 509 , 012065 ( 2014 ). J. Phys. Conf. Ser J. Steiheimer , J. Randrup , and V. Koch , Non-equilibrium phase transition in relativistic nuclear collisions: Importance of the equation of state . Phys. Rev. C 89 , 034901 ( 2014 ). doi: 10.1103/PhysRevC.89.034901 http://doi.org/10.1103/PhysRevC.89.034901 Non-equilibrium phase transition in relativistic nuclear collisions: Importance of the equation of state B. H. Alver , C. Gombeaud , M. Luzum , J. Y. Ollitrault , Triangular flow in hydrodynamics and transport theory . Phys. Rev. C 82 , 034913 ( 2010 ). doi: 10.1103/PhysRevC.82.034913 http://doi.org/10.1103/PhysRevC.82.034913 Triangular flow in hydrodynamics and transport theory S. Wang , Y. Z. Jiang , Y. M. Liu , et al ., Measurement of collective flow in heavy-ion collisions using particle-pair correlations . Phys. Rev. C 44 , 1091 ( 1991 ). doi: 10.1103/PhysRevC.44.1091 http://doi.org/10.1103/PhysRevC.44.1091 Measurement of collective flow in heavy-ion collisions using particle-pair correlations N. Borghini , P. M. Dinh , and J. Y. Ollitrault , Flow analysis from multiparticle azimuthal correlations . Phys. Rev. C 64 , 054901 ( 2001 ). doi: 10.1103/PhysRevC.64.054901 http://doi.org/10.1103/PhysRevC.64.054901 Flow analysis from multiparticle azimuthal correlations L. Xiong , Z. G. Wu , C. M. Ko , and J. Q. Wu , Dielectron production from nucleus-nucleus collisions . Nucl. Phys. A 512 , 772 - 786 ( 1990 ). doi: 10.1016/0375-9474(90)90234-D http://doi.org/10.1016/0375-9474(90)90234-D Dielectron production from nucleus-nucleus collisions T. Galatyuk , P. M. Hohler , R. Rapp , F. Seck , and J. Stroth , Thermal dileptons from coarse-grained transport as fireball probes at SIS energies . Eur. Phys. J. A 52 , 131 ( 2016 ). doi: 10.1140/epja/i2016-16131-1 http://doi.org/10.1140/epja/i2016-16131-1 Thermal dileptons from coarse-grained transport as fireball probes at SIS energies M. Gyulassy and X. N. Wang , HIJING 1.0: A Monte Carlo program for parton and particle production in high energy hadronic and nuclear collisions . Comput. Phys. Commun . 83 , 307 - 331 . ( 1994 ). doi: 10.1016/0010-4655(94)90057-4 http://doi.org/10.1016/0010-4655(94)90057-4 HIJING 1.0: A Monte Carlo program for parton and particle production in high energy hadronic and nuclear collisions T. Sjöstrand , High-energy-physics event generation with PYTHIA 5.7 and JETSET 7.4. Comput . Phys. Commun ., 82 , 74 - 89 . ( 1994 ). doi: 10.1016/0010-4655(94)90132-5 http://doi.org/10.1016/0010-4655(94)90132-5 High-energy-physics event generation with PYTHIA 5.7 and JETSET 7.4. Comput B. Zhang , ZPC 1.0.1: a parton cascade for ultrarelativistic heavy ion collisions . Comput. Phys. Commun ., 109 , 193 - 206 ( 1998 ). doi: 10.1016/S0010-4655(98)00010-1 http://doi.org/10.1016/S0010-4655(98)00010-1 ZPC 1.0.1: a parton cascade for ultrarelativistic heavy ion collisions B. A. Li and C. M. Ko , Formation of superdense hadronic matter in high energy heavy-ion collisions . Phys. Rev. C 52 , 2037 ( 1995 ). doi: 10.1103/PhysRevC.52.2037 http://doi.org/10.1103/PhysRevC.52.2037 Formation of superdense hadronic matter in high energy heavy-ion collisions B. A. Li and A. T. Sustich , B. Zhang , and C. M. Ko , Studies of superdense hadronic matter in a relativistic transport model . Int. J. Mod. Phys. E 10 , 267 ( 2001 ). doi: 10.1142/S0218301301000575 http://doi.org/10.1142/S0218301301000575 Studies of superdense hadronic matter in a relativistic transport model 10.1007/s41365-016-0141-3 This work was supported in part by the US Department of Energy under Contract No. DE-SC0015266 and the Welch Foundation under Grant No. A-1358. 2016-08-19 2016-09-23 2016-09-27 2016-10-17 2016-12-20 2021-04-11 © China Science Publishing & Media Ltd. (Science Press), Shanghai Institute of Applied Physics, the Chinese Academy of Sciences, ChineseNuclear Society and Springer Nature Singapore Pte Ltd. 2016 This is a post-peer-review, pre-copyedit version of an article published in Nuclear Science and Techniques. The final authenticated version is available online at: http://dx.doi.org/10.1007/s41365-016-0141-3 http://dx.doi.org/10.1007/s41365-016-0141-3 . 2016 Nuclear Science and Techniques 06 27