Dynamical formation of center domains in quark-gluon plasma

TL;DR

This paper investigates how center domains form dynamically in quark-gluon plasma during the deconfinement transition, using a phenomenological model fitted to lattice QCD data to understand domain formation times and structures.

Contribution

It introduces a phenomenological model for Polyakov loop propagation to study center domain formation in quenched QCD, connecting theory with lattice data.

Findings

Insights into the dynamical formation of center domains.

Estimation of formation timescales above critical temperature.

Observation of energy band formation along domain walls.

Abstract

We study the formation of domain structures due to spontaneous breakdown of center symmetry at high temperatures in quenched QCD. We develop a phenomenological model for the explicit propagation of the Polyakov loop as the relevant order parameter of the deconfinement phase transition. The surface tension in the equation of motion is fit in comparison with lattice QCD data. Results give insight into the dynamical formation of center domains as well as the formation of energy bands along domain walls and let us estimate the required time to form such structures above the critical temperature.

Figures

Equations

$$L(\vec{x})=\frac{1}{3}\mathrm{tr}P\exp \left[ig{\int }_0^{1/T}{A}_4(\tau ,\vec{x})d\tau \right],$$

$$U(L,T)=\left(-\frac{b_2}{2}{\left|L\right|}^2-\frac{b_3}{6}(L^3+{\bar{L}}^3)+\frac{1}{4}({\left|L\right|}^2{)}^2\right)b_4{T}^4,$$

$$b_2(T)=((1-1.11/x)(1+0.265/x{)}^2(1+0.3/x{)}^3-0.487)/r^2,$$

$$\mathcal{L}(L,T)=\frac{\sigma T_{\mathrm{c}}^2}{2}{\left|{\partial }_{\mu }L\right|}^2-U(L,T),$$

$$\sigma T_{\mathrm{c}}^2{\partial }_{\mu }{\partial }^{\mu }L+\frac{\partial U(L,T)}{\partial L}=0.$$

$$\sigma T_{\mathrm{c}}^2{\partial }_{\mu }{\partial }^{\mu }L+\eta \frac{\partial L}{\partial t}+\frac{\partial U(L,T)}{\partial L}=\xi .$$

$$⟨\xi (x,t)\xi (x^{\prime },t^{\prime })⟩=2\eta T\delta (x-x^{\prime })\delta (t-t^{\prime }).$$

$$\sigma =\{\begin{array}{ll}{\sigma }_0& \text{mbox}ifT\le T_{\mathrm{c}}\\ {\sigma }_0+a{\left(\frac{T-T_{\mathrm{c}}}{T_{\mathrm{c}}}\right)}^2& \text{mbox}ifT>T_{\mathrm{c}}\end{array}.$$

$$C(\left|x-y\right|)=⟨L_2(x)L_2(y){⟩}_T\propto {\mathrm{e}}^{-|x-y|/\xi },$$

$$\sigma (L_2)=\sqrt{⟨(\delta L_2{)}^2{⟩}_V},$$

Full text

I Introduction

The term quark-gluon plasma (QGP) denotes a state of matter, supposedly present in the early universe until about ${10}^{-5}\text{\,}\mathrm{s}$ after the Big Bang. This primordial phase is characterized by the deconfinement of color charges and the restoration of chiral symmetry. Nowadays, the necessary energy densities to recreate a QGP are reached in heavy-ion collision experiments such as the Relativistic Heavy-Ion Collider (RHIC) at BNL or the Large Hadron Collider (LHC) at CERN.

The transition from confinement to deconfinement is fairly well understood for pure gauge quantum chromodynamics (QCD), where it is possible to define an order parameter, the so-called Polyakov loop, related to the spontaneous breakdown of Z(3) center symmetry at high enough temperatures McLerran and Svetitsky (1981). For this case, lattice QCD data finds a first-order phase transition at $T_{\rm c}=$270\text{,}\mathrm{MeV}. With the inclusion of dynamical quarks, the transition temperature is significantly lowered to $T_{\rm c}=$150\text{\,}\mathrm{MeV}160\text{\,}\mathrm{MeV} and the phase transition is smeared out to an analytic crossover Aoki et al. (2006); Borsanyi et al. (2010). At high temperatures, the Z(3) symmetric Polyakov loop potential has three degenerate minima, leading to a domain structure in the deconfined phase, where center symmetry spontaneously breaks into different gauge configurations in different spatial regions. Lattice QCD studies have confirmed the existence of these center domains for the SU(2) gauge group Fortunato and Satz (2000, 2001); Fortunato (2003) and later for SU(3) Gattringer (2010); Danzer et al. (2010); Borsanyi et al. (2011); Endrodi et al. (2014). The formation of such structures necessarily comes with domain walls, interpolating between the different values of Z(3) in neighboring domains Bhattacharya et al. (1991, 1992). Center domains have been claimed to provide a simultaneous explanation for two distinct properties of QGP, namely, the low ratio of shear viscosity over entropy density $\eta/s$ and jet quenching Asakawa et al. (2013), properties that have been experimentally confirmed at the LHC Muller et al. (2012).

Domain formation via bubble nucleation at the confinement-deconfinement transition has been studied within an effective model in Gupta et al. (2010, 2012) both for quenched QCD and including dynamical quarks. In this work, the authors studied the evolution of QGP phase bubbles via the Kibble mechanism Kibble (1976); Zurek (1996) in a confining background in (2+1) dimensions, leading to Z(3) domains separated by domain walls and strings. It was shown that these structures cause inhomogeneities in the energy density even after cooling below $T_{\rm c}$, possibly influencing dilepton or direct photon distributions. It was argued furthermore that the expansion of these energetic fronts in the medium might leave imprints in experimentally detectable flow coefficients.

Our research here focuses on developing a dynamical model for an effective Polyakov loop field to describe the dynamical breaking of center symmetry and the formation of domains in the deconfined phase. For this purpose, we propose an effective Lagrangian consisting of a Polyakov loop potential from fits to lattice QCD and a phenomenological kinetic term, similar to what has been done in Dumitru and Pisarski (2001); Fraga et al. (2007); Herold et al. (2013, 2014). We study the evolution of the Polyakov loop field in an isothermal heat bath in (3+1) dimensions, considering the relaxational dynamics after temperature quenches to the plasma phase. Effects of the heat bath are included by using a Langevin equation of motion. This enables us to include thermal fluctuations and thus compare results to recent lattice QCD studies at finite temperature Borsanyi et al. (2011). Comparison of the correlation length to the one obtained in Borsanyi et al. (2011) allows us to fix the coefficient in front of the kinetic term in the region around $T_{\rm c}$. We are then able to study percolation of domains above the transition temperature and give estimates for formation times to better understand the possible role and relevance of these domain structure for heavy-ion collision experiments.

This paper is organized as follows: We begin with a description of the model in Sec. II, followed by the numerical procedure in Sec. III. In Sec. IV, we fit the surface tension from lattice QCD correlation lengths and discuss formation procedures and estimate formation times in Sec. V. We conclude with a summary and outlook in Sec. VI.

II Effective model

Confinement can be mathematically described by the Polyakov loop potential. The fundamental Polyakov loop is defined as

[TABLE]

where $P$ denotes the path-ordering operator, $g$ is the strong-coupling constant, $T$ is the temperature, and $A_{4}$ is the temporal component of a static gluon background field in Euclidean space-time.

From fits of lattice QCD data in the pure gluon sector, it is possible to obtain a polynomial potential for the Polyakov loop Wozar et al. (2006); Heinzl et al. (2005); Ratti et al. (2006); de Forcrand and von Smekal (2002); Pisarski (2000); Boyd et al. (1996); Roessner et al. (2007). In our work, we use the version from Roessner et al. (2007),

[TABLE]

with the temperature-dependent coefficient

[TABLE]

and the parameters $b_{3}=2/r$, $b_{4}=0.61r^{4}$, where $x=T/T_{c}$ and $r=2.23$. This potential leads to an expectation value of $\langle L\rangle=0$ at temperatures $T<T_{\rm c}$ with $T_{\rm c}=$270\text{,}\mathrm{MeV}$$. Above $T_{\rm c}$, spontaneous symmetry breaking leads to three degenerate states $\langle L\rangle=\mathrm{e}^{i2\nu\pi/3}$ with $\nu=0,1,2$, the three elements of the center subgroup Z(3). The transition at $T_{\rm c}$ is of first-order type. In the presence of dynamical quarks, it becomes an analytic crossover and the Z(3) symmetry is also explicitly broken, preferring the state $\langle L\rangle=1$ at high temperatures. One can account for this by adding a term proportional to $L$ to the potential in Eq. (2) as shown in Dumitru et al. (2004); Meisinger and Ogilvie (1995); Banks and Ukawa (1983); Green and Karsch (1984). Note that as this potential is obtained from lattice QCD fits in equilibrium, it is well determined around the respective minima, but afflicted with uncertainties away from them. As discussed in Asakawa et al. (2013), the logarithmic form would lead to values of $L\approx 0$ along the domain walls representing gauge configurations similar to those in the confined phase. Within our dynamical model, we will demonstrate the emergence of these walls in Sec. V.

Effective potentials for the Polyakov loop are often used in low-energy models such as the Polyakov loop Nambu-Jona-Lasinio (PNJL) model Fukushima (2004, 2008) or the Polyakov-Quark-Meson (PQM) model Schaefer et al. (2007); Herbst et al. (2011).

As the Polyakov loop carries no explicit time dependence, we apply an effective theory based on an effective Lagrangian of the form

[TABLE]

with the parameter $\sigma$ playing the role of a surface tension; cf. Fraga et al. (2007). It is clear that $\sigma$ will influence the domain size and therefore we shall assume it to be temperature-dependent, $\sigma\equiv\sigma(T)$. The factor of $T_{\rm c}^{2}$ has been added to account for the right dimensions in the fluctuation term Dumitru and Pisarski (2001). From the Lagrangian, we can obtain the equation of motion

[TABLE]

To capture the effect of thermal fluctuations, we go beyond this classical equation by using a Langevin equation of motion including dissipation and noise in a thermalized heat bath:

[TABLE]

The dissipation coefficient $\eta$ has been estimated in Fraga et al. (2007); Ananos et al. (2007) by studying the exponential growth of the correlation function of the Polyakov loop within Glauber dynamics of pure SU(3) lattice gauge theory. Relating Monte Carlo time and real time as in Tomboulis and Velytsky (2005) provides us with a value of $\eta=5/\textrm{fm}^{3}$. The stochastic noise field $\xi$ is Gaussian and white, and follows

[TABLE]

from the dissipation-fluctuation theorem. From Eqs. (4) and (5), it is clear that the correlation length $\xi$ depends on the coefficient $\sigma$ as $\xi\sim\sqrt{\sigma}$. It has been shown in Gattringer (2010); Borsanyi et al. (2011) that in quenched QCD, the cluster size $d=2\xi$ is constant below $T_{\rm c}$ at about $0.5\text{\,}\mathrm{fm}$ and then linearly rising above $T_{\rm c}$, thus signaling the phase transition. The authors argued that the value of $0.5\text{\,}\mathrm{fm}$ proves reasonable as it resembles the size of a heavy quark meson. We therefore assume a quadratic increase of $\sigma$ with temperature above the phase transition and use the following ansatz for the kinetic coefficient:

[TABLE]

In a next step, we fix the constants $\sigma_{0}$ and $a$ from the simulation by extracting the correlation length $\xi$ and comparing it with the values obtained in Borsanyi et al. (2011). The ansatz in Eq. (8) is of course an oversimplification as we neglect the curvature of the Polyakov loop potential which also influences the correlation length. It is also not clear if the linear rise in the cluster size is at some point halted as the data provided in Borsanyi et al. (2011) only reaches up to $1.2\,T_{\rm c}$.

III Numerical setup

To study the evolution of the center domains, we numerically solve the equation of motion (6) in (3+1) dimensions. We apply periodic boundary conditions for the spatial coordinates on a cubic lattice with $200^{3}$ sites, large enough to contain several larger domains at higher temperatures. We correlate the noise field over a volume of $(1/T)^{3}$ for each box with given temperature $T$ as we expect the correlation length to be of the order of $1/T$. We chose the lattice spacing to be $\Delta xN_{\rm corr}=1/T$ with $N_{\rm corr}=4$, following Bazavov et al. (2008) such that for each box, the noise field is correlated over a volume of $N_{\rm corr}^{3}$ cells. Finally, we set the time step to $\Delta t=0.01\cdot\Delta x$ to ensure numerical stability. At time $t=$0\text{,}\mathrm{fm}$$, the Polyakov loop $L=|L|\mathrm{e}^{i\phi}$ is initialized with a flat distribution in its argument $\phi$ and a Gaussian with width $0.1$ in the modulus $|L|$, corresponding to small fluctuations around the low-temperature expectation value of $L=0$. Note that Eq. (5) does not include the Hubble term; thus we neglect the effects of an expanding medium and consider a static box only.

IV Domain size

Our first goal is to determine the behavior of the kinetic coefficient $\sigma$ as a function of $T$. We use the temperature values studied in Borsanyi et al. (2011) as quench temperatures above $T_{\rm c}$ and then determine the correlation length $\xi$ from the two-point correlation function of the imaginary part $L_{2}$ of the Polyakov loop $L=L_{1}+iL_{2}$,

[TABLE]

taking the average $\langle\dots\rangle_{T}$ over an ensemble with the same temperature $T$. The value of $L_{2}$ serves to distinguish different types of domains which are characterized by $L_{2}=0$, $L_{2}=\sin(2\pi/3)$, and $L_{2}=-\sin(2\pi/3)$. For the real part $L_{1}$, only two different values occur. Alternatively, one could use the argument $\phi$ of $L$.

After setting the temperature in our box to the desired value above $T_{\rm c}$, we follow the evolution of the system to the deconfined phase, where the small initial fluctuations amplify and form center domains. As soon as the evolution has come to a halt and no more domains form or merge, we determine $\xi$ from the correlation function (9). We tune our input $\sigma$ such that the cluster diameter $d=2\xi$ equals the equally obtained $d$ from the lattice calculation Borsanyi et al. (2011). From a fit of $\sigma$ as function of $T$ according to Eq. (8), we finally obtain the values $\sigma_{0}=0.004$ and $a=4.3$. Next, we use $\sigma(T)$ to calculate the domain size for several values of $T$, ranging from $0.5\,T_{\rm c}$ to $1.3\,T_{\rm c}$, and show the result in Fig. 1. We see a clear resemblance between the data from our dynamical model and the lattice QCD results in the limited range investigated in Borsanyi et al. (2011). In both cases, the phase transition is clearly indicated by the sudden increase in $d$. For temperatures above the provided lattice QCD data, the cluster diameter from our model no longer increases linearly with temperature. This is due to the obvious oversimplification in our ansatz (8), where we neglected the curvature of the potential as an influence on the correlation length. Calculations of the volume of the largest percolating cluster in Endrodi et al. (2014) show that above $T_{\rm c}$, this volume increases less strongly with growing temperature. We may therefore expect a similar behavior for the average diameter, thus at least qualitatively justifying our result.

V Dynamical domain formation

To give an estimate for the formation time of domains, we use the standard deviation of $L_{2}$ over all cells in our box,

[TABLE]

with $\delta L_{2}=L_{2}-\langle L_{2}\rangle_{V}$ and the average $\langle\dots\rangle_{V}$ taken over all cells in our volume $V$. We expect $\sigma(L_{2})$ to saturate at some value as soon as the dynamical formation of domains in the isothermal heat bath has settled. We show this quantity as a function of time in Fig. 2 for several quench temperatures. Note here that the temperature is increased directly at $t=$0\text{,}\mathrm{fm}. We obtain values in the range between $17.0$ and $50.0\text{\,}\mathrm{fm}$ for the formation. This formation time increases with temperature as with increasing surface tension $\sigma$, the velocity $\partial L/\partial t$ decreases. We point out that the obtained curves lie closer and closer with increasing quench temperature, indicating a slowing down in the increase of the formation time. For center domains to play a relevant role in heavy-ion experiments, one needs to consider two things: First, the average domain volume should be smaller than the volume of the medium created after the collision, and second, the formation time has to be smaller than the lifetime of the QGP. Both quantities have been determined by Hanbury-Brown-Twiss (HBT) measurements Aamodt *et al.* ([2011](#bib.bib41)) who found a homogeneity volume at LHC of around 300\text{,}\mathrm{fm}$^{3}$, twice the size as at RHIC. This volume would be enough to contain several smaller domains, given the average domain size does not exceed about $7\text{\,}\mathrm{fm}$ for LHC or $5.5\text{\,}\mathrm{fm}$ for RHIC. The decoupling time has been determined to reach up to $10\text{\,}\mathrm{fm}11\text{\,}\mathrm{fm}$ at LHC energies in comparison to $7\text{\,}\mathrm{fm}8\text{\,}\mathrm{fm}$ at RHIC, well in agreement with previous transport model simulations Bass and Dumitru (2000). Therefore, taking into account our estimates for the required formation time, the possibility of domain formation at the LHC seems doubtful. However, it would be favorable to extend results to higher temperatures and also consider an expanding medium to simulate the situation in a heavy-ion collision. Work in this direction is currently in progress.

It is instructive to follow a visualization of the evolution of center domains. We show plots of the imaginary part of $L$ in Fig. 3 for constant value of $z$, corresponding to the curve of $T=1.1\,T_{\rm c}$ in Fig. 2. Center domains then appear as regions of the same color. We choose four distinct times, namely, at the beginning of the steep increase of $\sigma(L_{2})$ ($t=$20\text{,}\mathrm{fm}), during the increase ($t=$30\text{\,}\mathrm{fm}), at its end ($t=$40\text{,}\mathrm{fm}), and at the end of the simulation ($t=$50\text{\,}\mathrm{fm}). In the first two plots, the shape of the domains becomes more distinct while the field values move towards their respective equilibrium values. In the next two plots, we see the domains further sharpen in the edges and observe some domains merging with neighboring domains of the same type. It is important to note that these domains are all continuously connected and there is no abrupt change in the values, but a smooth transition between volumes with different Z(3) configuration. This border is characterized by a gauge configuration corresponding to the confined low-temperature phase, as already discussed in Sec. II. We demonstrate that this different gauge configuration effectively leads to a difference in the energy (see Fig. 4), which shows the energy density $\Delta U$ above the ground state at time $t=$50\text{,}\mathrm{fm}$$, corresponding to the plot of $L_{2}$ in Fig. 3d. Here, we clearly see bands of high energy at the edges of the center domains.

This is an important observation as these borders might act as potential barriers in two ways: Soft partons with thermal momenta will reflect on the domain walls which effectively limits their free wave length, resulting in a small value of $\eta/s$. Hard partons, on the other hand, may cross these walls under the emission of soft gluon radiation. Reflection of these gluons on the walls then makes the jet energy rapidly isotropic Asakawa et al. (2013). It would be interesting in the future to extend our model such that it can describe these two effects and estimate their impact on shear viscosity and jet quenching.

VI Summary and Outlook

In this article, we have introduced an effective model to describe the formation of Z(3) center domains via dynamical symmetry breaking in pure gluon QCD. We have fixed the surface tension in our equation of motion to obtain correlation lengths corresponding to the latest lattice QCD results. Our model signals the phase transition via a sudden increase in the size of center domains. We observed the mechanism of domain formation after a temperature quench, through the formation of small bubbles and the subsequent merging of bubbles with the same gauge configuration into domain structures. We were able to give estimates for the formation time, ranging from $17\text{\,}\mathrm{fm}50\text{\,}\mathrm{fm}$, larger than the estimated lifetime of the QGP at LHC energies. Although the presented results cover only a limited range up to $1.3\,T_{\rm c}$, we see that the formation time increases less and less with increasing temperature.

Finally, we have demonstrated the occurrence of energy bands along the walls of domains, which have previously been claimed to be a possible explanation for the simultaneous emergence of low $\eta/s$ and jet quenching.

We are going to continue our work by considering the effects of dynamical quarks though a linear term in the effective potential for the Polyakov loop, indicating the explicit symmetry breaking present in full QCD. We expect that in this case, domain formation will also occur, as lattice QCD calculations with dynamical quarks also lead to the formation of similar structures above the transition temperature Deka et al. (2012). Furthermore, dynamical models have also found domain formation if explicit symmetry breaking is taken into account Mohapatra and Srivastava (2013). Without a first-order phase transition, it can be expected that near $T_{c}$, the transition between the two phases proceeds more rapidly in a dynamical setup as there is no barrier to overcome, which significantly decreases the formation time for these scenarios. We are going to continue our work by considering the effects of dynamical quarks though a linear term in the effective potential for the Polyakov loop. To study the situation after the collision of two nuclei, it would furthermore be instructive to couple this model to a full (3+1)-dimensional hydrodynamical expansion to better understand what happens to center domains during expansion and cooling.

Acknowledgements

This work is funded by Suranaree University of Technology (SUT) and by the Office of the Higher Education Commission under NRU project of Thailand. C.H. and C.K. acknowledge support from SUT-CHE-NRU (FtR.15/2559) project.

The computing resources have been provided by the National e-Science Infrastructure Consortium of Thailand, the Center for Computer Services at SUT, and the Frankfurt Center for Scientific Computing.