# Influence of relativistic rotation on the confinement/deconfinement transition in gluodynamics

IInﬂuence of relativisic rotation on the conﬁnement/deconﬁnement transition ingluodynamics

V. V. Braguta,

1, 2, 3, ∗ A. Yu. Kotov, † D. D. Kuznedelev, ‡ and A. A. Roenko § Bogoliubov Laboratory of Theoretical Physics, Joint Institute for Nuclear Research, Dubna, 141980 Russia National University of Science and Technology MISIS, Leninsky Prospect 4, Moscow, 119049 Russia Moscow Institute of Physics and Technology, Dolgoprudny, 141700 Russia Jülich Supercomputing Centre, Forschungszentrum Jülich, D-52428 Jülich, Germany (Dated: February 11, 2021)In this paper we consider the inﬂuence of relativistic rotation on the conﬁnement/deconﬁnementtransition in gluodynamics within lattice simulation. We perform the simulation in the referenceframe which rotates with the system under investigation, where rotation is reduced to externalgravitational ﬁeld. To study the conﬁnement/deconﬁnement transition the Polyakov loop and itssusceptibility are calculated for various lattice parameters and the values of angular velocities whichare characteristic for heavy-ion collision experiments. Diﬀerent types of boundary conditions (open,periodic, Dirichlet) are imposed in directions, orthogonal to rotation axis. Our data for the criticaltemperature are well described by a simple quadratic function T c (Ω) /T c (0) = 1 + C Ω with C > for all boundary conditions and all lattice parameters used in the simulations. From this we concludethat the critical temperature of the conﬁnement/deconﬁnement transition in gluodynamics increaseswith increasing angular velocity. This conclusion does not depend on the boundary conditions usedin our study and we believe that this is universal property of gluodynamics. PACS numbers:Keywords: Lattice QCD, conﬁnement/deconﬁnement phase transition, Polyakov loop, rotation, heavy-ioncollisions, open boundary conditions, Dirichlet boundary conditions

I. INTRODUCTION

Recently the study of various physical systems under rotation has become relevant and extremely interesting researcharea. Rotating physical systems frequently appear in astrophysics [1, 2]. Relativistic fermions with angular momentumcan be realized in condensed matter physics [3, 4]. It is believed that rapidly rotating quark-gluon plasma is createdin heavy-ion collision experiments [5–8]. In the last example non-central heavy ion collisions generate nonzero angularmomentum. Partly this angular momentum is taken away by spectator partons, but considerable part is transferredto quark-gluon plasma, created in the collision. The experimental results for Λ , ¯Λ baryons polarization conﬁrm thisexpectation and give the following average value for the angular velocity Ω ∼ MeV [8]. Hydrodynamic simulationsof heavy-ion collisions predict even larger magnitudes of the angular velocity Ω ∼ (20 − MeV [5]. These values ofangular velocity lead to relativistic rotation of quark-gluon plasma.Rotation gives rise to lots of interesting phenomena which can be observed in heavy-ion collision experiments. Forinstance, chiral vortical eﬀect [9–12] and polarization of diﬀerent particles [13, 14] are examples of such phenomena.In addition, relativistic rotation is believed to inﬂuence phase transitions in QCD what also can be observed in theexperiments.There are a lot of theoretical papers dedicated to the phase transitions in rotating QCD matter (see, for instance, [15–22]). Mostly these studies are carried out within Nambu–Jona-Lasinio model (NJL) [23, 24] and they are focusedon inﬂuence of rotation on the chiral symmetry breaking/restoration transition in QCD. Despite variety of detailsand results of there papers, there is one common result: rotation suppresses the chiral condensate, which leads tothe decrease of the critical temperature with rotation. An interesting physical explanation of this suppression wasproposed in paper [17]. The idea is that rotation induces a polarization eﬀect which aligns all microscopic angularmomenta along the total angular momentum. So the states with nonzero angular momentum/spin are more preferredthan the states with zero angular momentum/spin which results into the suppression of the chiral condensate exposedto the rotation. ∗ Electronic address: [email protected] † Electronic address: [email protected] ‡ Electronic address: [email protected] § Electronic address: [email protected] a r X i v : . [ h e p - l a t ] F e b An important disadvantage of all works based on the NJL model is that they take into account only the quarkdegrees of freedom whereas the gluon sector of the theory, which is responsible for the conﬁnement, is integrated out.For this reason it is diﬃcult to study conﬁning properties as well as conﬁnement/deconﬁnement transition withinsuch approaches. For QCD matter under rotation this disadvantage might be crucial since gluons have spin-1 andone can expect that relativistic rotation induces polarization of gluons. The properties of QCD matter with po-larized gluons might be diﬀerent as compared to that without polarization. In particular, rotation can aﬀect theconﬁnement/deconﬁnement transition. So, in order to understand the impact of rotation to the phase transitionsin QCD matter, one needs to apply an approach which takes into account both quark and gluon degrees of free-dom. The papers [20–22] are focused on the conﬁnement/deconﬁnement transition in rotating QCD. The authorsof [21] applied holographic approach. The author of [20] studied rotating compact QED which also possesses con-ﬁnement/deconﬁnement transition. The hadron resonance gas model was adopted in paper [22]. The results of theseworks indicate that the critical temperature decreases with the angular velocity. Although the results of diﬀerentphenomenological studies give interesting insights into QCD properties, we believe that lattice simulation of QCD isthe most appropriate to study how rotation inﬂuences the conﬁnement/deconﬁnement transition.This paper is devoted to the study of SU(3) gluodynamics properties under rotation within lattice simulation. Inparticular, we are going to address the question how rotation aﬀects the conﬁnement/deconﬁnement transition ingluodynamics. It is worth mentioning that the ﬁrst lattice study of rotating QCD matter was carried out in [25], butthe impact of rotation on QCD phase transition was not considered in this paper. In our paper we will use latticeformulation of rotating gluodynamics developed in [25]. We would like also to note that our ﬁrst results devoted tothermodynamic properties of rotation gluodynamics were published in [26]. In this paper we continue our study.This paper is organized as follows. Next section is devoted to the theoretical background of our study. In partic-ular, we discuss Yang-Mills theory in external gravitational ﬁeld, write its discretized action and describe boundaryconditions used in our study. In Section III the results of our study obtained with diﬀerent boundary conditions arepresented. In last section we discuss the results and draw a conclusion. In addition in Appendix A the inﬂuenceof ﬁnite volume eﬀects on thermodynamic properties of gluodynamics is studied. In Appendix B we show latticeparameters used in our simulations.

II. THEORETICAL BACKGROUNDA. Thermodynamical ensemble in presence of rotation

To study the inﬂuence of rotation on gluodynamics properties we are going to use the approach proposed in papers[15–19, 25]. The idea is to carry out the study in the reference frame which rotates with the system. Below it will beassumed that the system rotates around z axis. In this reference frame there appears an external gravitational ﬁeldwith well known metric tensor g µν = − r Ω Ω y − Ω x y − − Ω x − − , (1)where r = (cid:112) x + y is the distance to the axis of rotation.All components of the metric tensor (1) do not depend on time coordinate t . As the result, the Hamiltonian of thesystem which is given by the expression H = (cid:90) dV √ g (cid:15) ( (cid:126)r ) = 12 g (cid:90) d x √ γ √ g (cid:20) F atx F atx + F aty F aty + F atz F atz + F axy F axy + F axz F axz + F ayz F ayz −− Ω (cid:0) xF atx F axy − xF atz F ayz + yF aty F axy + yF atz F axz (cid:1)(cid:21) , (2)is conserved. Here the Greek letters correspond to the Lorentz indices while Latin letters correspond to the colorones. The dV = d x √ γ is the three dimensional volume in our coordinate system, g is strong coupling constant and (cid:15) ( (cid:126)r ) is energy density. Notice that to write conserved quantity it is important to multiply the energy density (cid:15) ( (cid:126)r ) byadditional contribution of the gravitational ﬁeld √ g .Given the Hamiltonian (2), it is straightforward to construct the partition function of the system under study Z = Tr exp (cid:20) − β ˆ H (cid:21) . (3)Here one comment is in order. One can introduce the notation T ( r ) = 1 /β √ g which brings the expression (3) tothe from Z = Tr exp (cid:20) − (cid:90) dV ˆ (cid:15) ( r ) T ( r ) (cid:21) . (4)The eﬀective T ( r ) is the temperature which depends on the position in space and satisﬁes the expression T ( r ) √ g =1 /β = const . Last equation describes Ehrenfest–Tolman eﬀect which states that in gravitational ﬁeld the temperatureis not constant in space in thermal equilibrium. For the rotation one has T ( r ) √ − Ω r = 1 /β = T ( r = 0) . So, onecan conclude that rotation eﬀectively heats the system from the rotation axis to the boundaries T ( r ) > T ( r = 0) . Inwhat follows the temperature at the rotation axis T ( r = 0) = 1 /β will be referred to as T .The calculation of the Tr exp[ ... ] in formula (3) can be carried out through the standard procedure. Applying it,the partition function of the gluodynamics in external gravitational ﬁeld can be written in the following form [25] Z = (cid:90) DA exp ( − S G ) . (5)In last formula the integration over gluon degrees of freedom is carried out. The S G is the Euclidean action of thegluodynamics in external gravitational ﬁeld, which can be written as S G = 14 g (cid:90) d x √ g E g µνE g αβE F aµα F aνβ . (6)The Euclidean metric tensor ( g E ) µν in last formula can be obtained from (1) by Wick rotation t → iτ . As in usualpath integral for the partition function the Euclidean time τ varies in the region τ ∈ (0 , β ) and the gluon degrees offreedom satisfy periodic boundary conditions in temporal direction A µ (0 , x ) = A µ ( β, x ) .Substituting the ( g E ) µν to formula (6) we get the following expression for the Euclidean action S G = 12 g (cid:90) d x (cid:2) (1 − r Ω ) F axy F axy + (1 − y Ω ) F axz F axz + (1 − x Ω ) F ayz F ayz + F axτ F axτ + F ayτ F ayτ ++ F azτ F azτ − iy Ω( F axy F ayτ + F axz F azτ ) + 2 ix Ω( F ayx F axτ + F ayz F azτ ) − xy Ω F axz F azy (cid:3) . (7)It is seen from this formula that the action is a complex function what leads to the sign problem. Unfortunately,direct Monte Carlo simulation of this system is impossible today. To overcome this problem instead of the real angularvelocity Ω , we are going to conduct Monte Carlo simulations with imaginary angular velocity Ω I = − i Ω . The resultsobtained in this way will be expanded in the Ω I and then analytically continued to real angular velocity. B. Lattice formulation for rotating gluodynamics

In order to conduct lattice simulation of rotating gluodynamics, one has to discretize action (7). In this paper weare going use the lattice action proposed in [25], which can be written as S G = 2 N c g (cid:88) x (cid:16) (1 + r Ω I )(1 − N c Re Tr ¯ U xy ) + (1 + y Ω I )(1 − N c Re Tr ¯ U xz ) ++ (1 + x Ω I )(1 − N c Re Tr ¯ U yz ) + 3 − N c Re Tr ( ¯ U xτ + ¯ U yτ + ¯ U zτ ) −− N c Re Tr (cid:0) y Ω I ( ¯ V xyτ + ¯ V xzτ ) − x Ω I ( ¯ V yxτ + ¯ V yzτ ) + xy Ω I ¯ V xzy (cid:1)(cid:17) , (8)where the ¯ U µν denotes clover-type average of four plaquettes (see Fig. 1). In the ﬂat metric case, the use of cloversinstead of plaquettes would lead to the same action after summation. However, in the case of non-uniform metric,the clover-type average allows one to build local expressions for the F µν and reduce discretization errors.In order to implement the terms F µα F αν on the lattice one uses the antisymmetric chair-type average ¯ V µνρ of 8chairs (see Fig. 2). The intuition for such a term in a discretized action may be understood as follows: we have termsof the form F µν F νρ in the action, and loop in the µν ( νρ ) plane gives respectively F µν ( F νρ ) . The simplest gaugeinvariant object has to reside in µν and νρ plane and thus involves at least 6 links.We have performed numerical simulation on lattices N t × N z × N x × N y = N t × N z × N s ( N s = N x = N y ), suchthat sizes in spatial directions, orthogonal to the rotation axis z , are equal and may diﬀer from size along the rotationaxis and temporal direction. The axis z is in the middle of xy -plane. For the all sizes considered, the restriction,following from causality: ( Ω N s a/ √ < ) holds. And in the most cases we have Ω N s a/ √ (cid:28) . µ ν U µν = 14 Figure 1: The clover-type average of four plaquettes. µρ ν µρ ν V µνρ = 18 Figure 2: The antisymmetric chair-type average of eight chairs.

C. Boundary conditions

It is clear that rotating reference frame cannot be extended to arbitrary large distances from the rotation axis, sinceat distances Ω r ≥ the g becomes negative and such rotating system cannot be realized by real bodies. For thisreason in the simulations one has to impose boundary conditions (BC) on our system. Here we would like to stressthat BC are important part of all approaches aimed at studying of rotating quark-gluon plasma rather than a latticeartifact. The results obtained within any approach depends on BC. In order to study this dependence in our paperwe applied a series of diﬀerent BC.For all BC used in our paper we apply periodic boundary conditions for gluon ﬁelds in z - and τ -directions. Whatconcerns the x - and y -directions we used the following BC1. Open boundary conditions (OBC).

For this type of boundary conditions the sum in action (8) is taken overplaquettes and chair-type products of link variables which belong to the lattice volume under investigation. Ifeither plaquette or chair goes beyond the lattice volume it is excluded from the action. No additional restrictionsare put on the link variable located on the boundary. Similar conditions were used in Ref. [27, 28] for temporaldirection and in Ref. [29] for one of spatial directions, where topological susceptibility was investigated. Here we would like to mention that OBC do not violate the Z symmetry and we do not see explicit incompat-ibility of OBC with the ﬁeld of velocities of rotating body. For this reason, we believe that OBC are the mostappropriate for the lattice simulation of rotating gluodynamics.The exclusion of the plaquettes outside the studied lattice volume from the action can be considered as if oneputs these plaquettes to unity what leads to zero lattice action. So, physically, OBC can be interpreted as ifthe studied volume is attached to classical (without quantum ﬂuctuations) zero temperature Yang-Mills theory.Lattice study of non-rotating lattice with OBC conﬁrms this physical picture (see Appendix A).2. Periodic boundary conditions (PBC) for gluon ﬁelds in x - and y -directions. This type of BC is not compatiblewith the velocity distribution in the rotating body but it respects the Z symmetry of the action. Notice that in papers [27–29] the plaquettes on the boundary were accounted in the action with the weight 1/2. We studied this variantof OBC for non-rotating lattice and did not found meaningful diﬀerence between these two approaches. Dirichlet boundary conditions (DBC).

In this case we set all links which lie on the boundary to unit matrix: U µ ( x ) = 1 . The links sticking out from the lattice volume are not included in the lattice action. DBC were usedin Ref. [25] to study rotating QCD. In the continuum limit this corresponds to A µ ( (cid:126)x, τ ) | x = ± R/ = 0 , µ = 0 , , , (9a) A µ ( (cid:126)x, τ ) | y = ± R/ = 0 , µ = 0 , , , (9b)where R is the size of lattice in the directions x and y .We do not see incompatibility of DBC with the ﬁeld of velocities of rotating body, but DBC break the Z symmetry of the action. In principle, without rotation this explicit violation of the Z symmetry is ﬁnite-volume eﬀect which disappears in the inﬁnite volume limit (see Appendix A). However, as explained above, onecannot take inﬁnite volume limit for the x and y directions in the rotating reference frame. So, in the rotatingreference frame the violation of the Z symmetry might play important role.It is interesting to note that physically DBC can be considered as if one ﬁxes high temperature on the boundaryof the volume under investigation. This can be understood as follows. In DBC the Polyakov loop on theboundary equals three, i.e. the Z symmetry of lattice gluodynamics is explicitly broken on the boundary. Ingluodynamics breaking of Z symmetry is the property of high temperature phase. For this reason in DBC theboundary plays a role of a seed of high temperature phase or, in other words, this can be considered as if thereis a high temperature on the boundary. From this perspective DBC is opposite to OBC with low temperatureoutside the lattice volume and it is particularly interesting to compare the results obtained with OBC and DBC.Numerical simulations of non-rotating lattice conﬁrms the physical picture of high eﬀective temperature on theboundary (see Appendix A).One can expect that BC incompatible with the properties of rotating body might lead to unphysical behaviour of thesystem. Consequently, the approach applied in this paper might become inapplicable to study rotating gluodynamics.The results obtained in this paper allow us to state that for suﬃciently large lattice volumes and small angularvelocities BC do not change bulk properties of rotating gluodynamics considerably. We believe this fact can beexplained as follows. For suﬃciently large volume and small angular velocity it becomes energetically favorable forthe system to screen BC incompatible with the bulk properties of rotating body and the bulk properties start todominate over the screened boundary. The screening of the boundary is well seen in Fig. 6 and Fig. 14. D. Measurement of the critical temperature

The main question to be addressed in this paper is how rotation inﬂuences conﬁnement/deconﬁnement phasetransition in gluodynamics. Commonly in non-rotating gluodynamics one exploits the

Polyakov loop as an orderparameter for this transition. Lattice expression for the Polyakov loop can be written as L ( (cid:126)x ) = Tr (cid:34) N t − (cid:89) τ =0 U ( (cid:126)x, τ ) (cid:35) , (10)where U ( (cid:126)x, τ ) is the gauge link variable in the temporal direction.In non-rotating gluodynamics the Polyakov loop can be used as an order parameter since lattice action of gluo-dynamics is invariant under the multiplication of the U link elements in some time slice by center elements of the SU (3) group ( Z symmetry). At the same time the Polyakov loop is not invariant under this transformation. Fromthe symmetry perspective, the conﬁnement is Z symmetric phase with (cid:104) L (cid:105) = 0 , whereas in the deconﬁnement Z symmetry is broken and (cid:104) L (cid:105) (cid:54) = 0 . If we now turn to the rotating gluodynamics, it is clear that the action (8) possesses Z symmetry and above arguments on treating the Polyakov loop as an order parameter persist.In our paper we are going to use the Polyakov loop to label the phases of rotating gluodynamics. In addition tolocal Polyakov loop at point (cid:126)x , we are going use spatially averaged Polyakov loop L = 1 N s N z (cid:88) (cid:126)x L ( (cid:126)x ) , (11)The critical temperature of the conﬁnement/deconﬁnement transition will be determined from the Polyakov loopsusceptibility χ = N s N z (cid:0) (cid:104)| L | (cid:105) − (cid:104)| L |(cid:105) (cid:1) , (12)which has a peak at the critical temperature. We use the Gaussian ﬁt over a set of points near the susceptibility peakto calculate the critical temperature from obtained data χ ( T ) = A + B exp (cid:18) − ( T − T c ) δT (cid:19) . (13)One might assume that due to boundary conditions and inhomogeneity of rotating gluodynamics the Polyakov loop andits susceptibility are not appropriate for ﬁnding the critical temperature. However we believe that these observablescan be used to study conﬁnement/deconﬁnement transition for the following reasons:• In the gluodynamics with OBC (both with and without rotation) observables can depend on the coordinate,because these boundary conditions break translational symmetry. However, since OBC respect Z symmetry,in the conﬁnement phase the Polyakov loop is zero and does not depend on the spatial coordinate (see Fig. 6).In the deconﬁnement phase the Polyakov loop is not zero and the dependence on spatial coordinate appears.From these facts it is clear that the spatially averaged Polyakov loop (11) and its susceptibility (12) can beused to detect the critical temperature. In Appendix A it is shown that in the inﬁnite volume limit the criticaltemperatures in non-rotating gluodynamics with OBC and PBC agree with each other as it should be.• DBC also break translational symmetry. At the same time these boundary conditions, contrary to OBC andPBC, break Z symmetry even at zero temperature: the L ( (cid:126)x ) = 3 on the boundary, and, as the result, it is notzero in the bulk. Under these circumstances the Polyakov loop becomes an approximate order parameter. Stillit is possible to detect the pseudo-critical temperature of this crossover through the peak of the susceptibility(12). For non-rotating gluodynamics the ﬁrst order phase transition and its critical temperature is recovered inthe inﬁnite volume limit N z → ∞ , N s → ∞ (see Appendix A). Unfortunately for rotating gluodynamics thelimit N s → ∞ cannot be taken and the conﬁnement/deconﬁnement transition remains a crossover even in thelimit N z → ∞ .• For all BC in rotating gluodynamics the Polyakov loop acquires additional dependence on spatial coordinatedue to the rotation. Because of Z symmetry, for the PBC and OBC the Polyakov loop is zero and independenton the space coordinate in the conﬁnement phase (see Fig. 10 and Fig. 6). In the deconﬁnement the Polyakovloop is non-zero and depends on the spatial coordinate. So, the spatially averaged Polyakov loop (11) andits susceptibility (12) can be used to detect the critical temperature even in rotating gluodynamics. As wasexplained above, for DBC the Polyakov loop is an approximate order parameter, but its susceptibility can beused to ﬁnd the pseudo-critical temperature for rotating gluodynamics.At the end of this section we would like to note that for all BC used in this paper the value of the critical temperature T c (Ω I ) contains ﬁnite volume eﬀects which depend on the lattice size N s (see the discussion in the Appendix A). Inorder to account for these eﬀects below our results for the critical temperature will be presented in terms of the ratio T (Ω I ) /T c (0) . Finally, we would like to mention that the detailed description of used lattice parameters for all the BCis presented in the Appendix B. III. THE RESULTS OF THE CALCULATIONA. Open boundary conditions

We believe that OBC are the most appropriate for the lattice study of the rotating matter. For this reason we startthe discussion of the lattice results with OBC.The Polyakov loop and the Polyakov loop susceptibility as functions of the ratio

T /T c (0) for various values of (imag-inary) angular velocity Ω I for the lattice size × × are shown in Fig. 3. The conﬁnement/deconﬁnement phasetransition manifests itself as a rapid growth of the Polyakov loop and correspondingly as a peak in the susceptibility.One can easily read from Fig. 3 that the phase transition is shifted to lower temperatures when the (imaginary)angular velocity Ω I grows. To make quantitative predictions, we ﬁt several points in the transition region for thePolyakov loop susceptibility with Gaussian function (13). The χ / ndof of the ﬁt is ∼ . − for all angular velocities Ω I . The ratios T c (Ω I ) /T c (0) as functions of Ω I are presented in the Fig. 4(a). It is also instructive to introduce the(imaginary) linear velocity v I at the points with the coordinates x = ± R/ , y = 0 which are located on the boundary: v I = Ω I ( N s − a/ and to present the ratios T c ( v I ) /T c (0) as functions of v I (see Fig. 4(b)). In order to assess Note that to determine v I we used the lattice spacing a = a ( β c ) at the critical temperature. . . . . . T/T c (0)0 . . . . . . . h | L | i × × , OBCΩ I = 0 MeVΩ I = 10 MeVΩ I = 20 MeVΩ I = 24 MeVΩ I = 30 MeV (a) .

75 0 .

80 0 .

85 0 .

90 0 .

95 1 .

00 1 .

05 1 . T/T c (0)0510152025 χ × × , OBCΩ I = 0 MeVΩ I = 10 MeVΩ I = 20 MeVΩ I = 24 MeVΩ I = 30 MeV (b) Figure 3: The Polyakov loop (a) and the Polyakov loop susceptibility (b) as a function of temperature for diﬀerent values ofimaginary angular velocity Ω I . The results are obtained on the lattice × × with OBC. The lines for the Polyakovloop (a) are drawn to guide the eye. The Polyakov loop susceptibilities (b) are ﬁtted in the vicinity of the phase transition bya Gaussian function (13). I (MeV )0 . . . . . T c ( Ω I ) / T c ( ) N s /N t ’ N s /N t ’ N s /N t ’ OBC8 × × × × × × × × × × × × × × × × × × (a) .

00 0 .

02 0 .

04 0 .

06 0 .

08 0 .

10 0 .

12 0 . v I /c . . . . . T c ( v I ) / T c ( ) N s /N t ’ N s /N t ’ N s /N t ’ OBC8 × × × × × × × × × × × × × × × × × × (b) Figure 4: The ratios T c /T c (0) determined from the peak of the Polyakov loop susceptibility as a function of the imaginaryangular velocity squared Ω I (a) and the linear boundary velocity squared v I (b). Results are presented for several lattice sizes N t × N z × N s with OBC. Lines correspond to simple quadratic ﬁts T c (Ω I ) /T c (0) = 1 − C Ω I and T c ( v I ) /T c (0) = 1 − B v I /c systematic eﬀects, in Fig. 4 we present the results for various lattice sizes.Based on the results, presented in Fig. 3 and Fig. 4, one can draw the following conclusions:• The T c decreases with increasing Ω I . We have found that for the studied parameters the dependence of the T c (Ω I ) on the imaginary angular velocity Ω I can be described by a simple quadratic function ( χ / ndof ∼ . − ) T c (Ω I ) T c (0) = 1 − C Ω I . (14)This conﬁrms that the angular velocities used in the simulation are indeed small and one can expand the criticaltemperature in a series over Ω I . Upon analytical continuation to real angular velocity Ω I → i Ω , which islegitimate for small angular velocities, one gets the following dependence of the critical temperature on the . . . . . . . . . N s /N t . . . . . . . B OBC8 × × N s × × N s × × N s × × N s × × N s × × N s × × N s Figure 5: The coeﬃcient B in Eq. (16) versus the ratio N s /N t for several lattice sizes with OBC. − −

10 0 10 20 x . . . . . . . . | h L ( x , y = ) i | OBC, Ω I = 0 MeV T/T c (0) = 0 . T/T c (0) = 1 . (a) − −

10 0 10 20 x . . . . . . . . | h L ( x , y = ) i | OBC, Ω I = 24 MeV T/T c (0) = 0 . T/T c (0) = 1 . (b) Figure 6: The Polyakov loop |(cid:104) L ( x, y = 0) (cid:105)| as a function of coordinate x for OBC and Ω I = 0 MeV (a), Ω I = 24 MeV (b).The results were obtained on the lattice × × for two temperatures: T /T c (0) = 0 . in the conﬁnement phase and T /T c (0) = 1 . in the deconﬁnement phase. value of the Ω : T c (Ω) T c (0) = 1 + C Ω . (15)Our results indicate that the C > , which leads to the conclusion: With OBC the critical temperature of theconﬁnement/deconﬁnement phase transition grows with increasing angular velocity. • In order to study the dependence of our results on the N z lattice size we calculated the critical temperature onthe lattices × N z × , N z = 20 , , . The results obtained on these lattices agree within the uncertainty (seeFig. 4(a)). In order to study discretization eﬀects, we conducted our study on the lattices × × , × × , × × where the physical sizes are kept ﬁxed. As can be seen from Fig. 4(a), the ratio T c (Ω I ) /T c (0) shows almost no dependence on the lattice spacing a . Next we proceeded to the dependence of the results onsize in the transverse directions N s . To do this we ﬁxed the N t and N z sizes and varied the N s . It is seen fromFig. 4(a) that our data are split into lines with diﬀerent slopes. The dependence of these slopes (diﬀerent C constants) on the lattice sizes N s is quite signiﬁcant.This phenomenon can be understood in the following way. Increasing lattice size in directions, orthogonal to therotation axis, leads to the increase of the linear velocity on the boundary of the rotating lattice. Since this linearvelocity ∼ r Ω enters the metric tensor, the action and the expressions for local temperature T ( r ) = 1 /β √ g ,it is reasonable to assume that the T c is a function of some “velocity-like” parameter ∼ r Ω , not an (imaginary)angular velocity itself. In order to check this assumption, in Fig. 4(b) we present the ratio T ( v I ) /T c (0) as afunction of linear velocity v I on the boundary. Quite remarkably on this ﬁgure all points show a clear tendencyto lie on one line. Data can be well described by a simple quadratic function T c ( v I ) T c (0) = 1 − B v I c , (16)which corresponds to the following relation for real rotation: T c ( v ) T c (0) = 1 + B v c . (17)We present the values of B for several sets of parameters in Fig. 5. The coeﬃcient B has a mild dependenceon the parameters of the system. It does not change upon changing the lattice spacing a and slightly growswith increasing lattice extent N s in x - and y -directions. It is reasonable to assume that for suﬃciently large N s and small angular velocity the bulk dominates over the boundary, i.e. the role of the boundary becomesless important. We believe that this property manifests itself when the B goes to the plateau in Fig. 5 for N s /N t > . So, our second conclusion is: The dependence of the critical temperature on the linear velocity atthe boundary v has the form (17), with weak dependence of the B on the lattice parameters. For lattices withsuﬃciently large N s and OBC the coeﬃcient is B ∼ . . • It is instructive to study how the Polyakov loop depends on the spatial coordinate. Since BC and rotation breaktranslational invariance in x - and y - directions, but preserve the invariance in z -direction, we introduce the localPolyakov loop in x, y -plane L ( x, y ) = 1 N z (cid:88) z L ( x, y, z ) , (18)where L ( x, y, z ) = L ( (cid:126)x ) is deﬁned by Eq. (10) and study its ensemble average. In Fig. 6 we present the Polyakovloop |(cid:104) L ( x, y = 0) (cid:105)| as a function of the coordinate x for the lattice × × . The results are shown for twotemperatures: T /T c (0) = 0 . in the conﬁnement phase and T /T c (0) = 1 . in the deconﬁnement phase. Inaddition we plot data for the lattice without rotation Ω I = 0 MeV (Fig. 6(a)) and with Ω I = 24 MeV (Fig. 6(b)).Fig. 6 illustrates the features of Polyakov loop discussed in Section II D. One sees that Polyakov loop |(cid:104) L ( x, y ) (cid:105)| iszero for all spatial points in the conﬁnement phase, both without rotation and with nonzero angular velocity. Itconﬁrms, that for OBC the average Polyakov loop still acts as the order parameter of conﬁnement–deconﬁnementphase transition. In the deconﬁnement phase one sees nontrivial coordinate dependence of Polyakov loop. Mainlythis dependence can be attributed to the inﬂuence of OBC. When one moves from the boundary to the bulk,one observes that in the deconﬁnement phase the boundary is screened. B. Periodic boundary conditions

In this section we present the results of our study of the conﬁnement/deconﬁnement phase transition for rotatinggluodynamics with PBC. Although PBC are commonly used in lattice simulations, they might not be physical forthe simulation of rotating medium. Nonetheless, we believe that they can be used to check the robustness of ourpredictions against changing the boundary conditions.In Fig. 7 we present the Polyakov loop and the Polyakov loop susceptibility with respect to the temperature forvarious values of (imaginary) angular velocity Ω I for the lattice size × × . It is clearly seen that the behaviourof the critical temperature is the same as for OBC: it decreases with growing imaginary angular velocity.To determine the critical temperature we ﬁtted the susceptibility in the vicinity of the peak by the Gaussianfunction (13). In Fig. 8 we present the obtained dependence of the ratios T c /T c (0) on the imaginary angular velocity Ω I (Fig. 8(a)) and the corresponding linear boundary velocity v I (Fig. 8(b)) for several lattice sizes.In general, the behaviour of the conﬁnement/deconﬁnement phase transition in the rotating gluodynamics withPBC is very similar to the case with OBC. In particular, it is worth to mention the following:0 . . . . . . T/T c (0)0 . . . . . . . . . h | L | i × × , PBCΩ I = 0 MeVΩ I = 10 MeVΩ I = 15 MeVΩ I = 20 MeVΩ I = 24 MeV (a) .

75 0 .

80 0 .

85 0 .

90 0 .

95 1 .

00 1 .

05 1 . T/T c (0)010203040 χ × × , PBCΩ I = 0 MeVΩ I = 10 MeVΩ I = 15 MeVΩ I = 20 MeVΩ I = 24 MeV (b) Figure 7: The Polyakov loop (a) and the Polyakov loop susceptibility (b) as a function of temperature for diﬀerent values ofimaginary angular velocity Ω I . The results are obtained on the lattice × × with PBC. The lines for the Polyakovloop (a) are drawn to guide the eye. The Polyakov loop susceptibilities (b) are ﬁtted in the vicinity of the phase transition bythe Gaussian function (13). I (MeV )0 . . . . . . . . . T c ( Ω I ) / T c ( ) N s /N t ’ N s /N t ’ N s /N t ’ PBC8 × × × × × × × × × × × × × × × × × × × × (a) .

00 0 .

02 0 .

04 0 .

06 0 .

08 0 . v I /c . . . . . . . . . T c ( v I ) / T c ( ) PBC8 × × × × × × × × × × × × × × × × × × × × (b) Figure 8: The ratios T c /T c (0) determined from the peak of the Polyakov loop susceptibility as a function of the imaginaryangular velocity squared Ω I (a) and the linear boundary velocity squared v I (b). Results are presented for several lattice sizes N t × N z × N s with PBC. Lines correspond to simple quadratic ﬁts T c (Ω I ) /T c (0) = 1 − C Ω I and T c ( v I ) /T c (0) = 1 − B v I /c • In Fig. 8(a) we present the dependence of the ratio T c (Ω I ) /T c (0) on the imaginary angular velocity Ω I . One caneasily see that with good accuracy it can be described by the same formula, as for OBC: T c (Ω I ) /T c (0) = 1 − C Ω I , with C > . After analytical continuation we draw a conclusion that with PBC the critical temperature of theconﬁnement/deconﬁnement transition increases with angular velocity .• Analogously to OBC, with PBC we have performed simulations for several lattice sizes, in order to estimatesystematic uncertainties. It is seen that uncertainties of the calculation for PBC are smaller than that for OBC.By looking at the results of the simulations with lattice sizes × N z × , N z = 24 , , we conclude thateﬀects of ﬁnite lattice size in z -direction are small. From the results obtained on the lattice sizes × × , × × and × × with ﬁxed physical volumes one can read oﬀ, that the dependence on thelattice spacing is rather weak. However, it can be easily seen that when one varies lattice size N s , the ratio1 . . . . . . . N s /N t . . . . . . . B PBC8 × × N s × × N s × × N s × × N s × × N s Figure 9: The coeﬃcient B in Eq. (16) versus the ratio N s /N t for several lattice sizes with PBC. − −

10 0 10 20 x . . . . . . . | h L ( x , y = ) i | PBC, Ω I = 0 MeV T/T c (0) = 0 . T/T c (0) = 1 . (a) − −

10 0 10 20 x . . . . . . . . | h L ( x , y = ) i | PBC, Ω I = 24 MeV T/T c (0) = 0 . T/T c (0) = 1 . (b) Figure 10: The Polyakov loop |(cid:104) L ( x, y = 0) (cid:105)| as a function of coordinate x for PBC and Ω I = 0 MeV (a), Ω I = 24 MeV (b).The results were obtained on the lattice × × for two temperatures: T /T c (0) = 0 . in the conﬁnement phase and T /T c (0) = 1 . in the deconﬁnement phase. T c (Ω I ) /T c (0) changes signiﬁcantly. Again, similarly to OBC, this dependence can be absorbed by looking atthe T c ( v I ) /T c (0) versus the linear boundary velocity squared v I (see Fig. 8(b)): T c ( v I ) /T c (0) = 1 − B v I /c .We present the results for the coeﬃcient B in Fig. 9. It is worth noting, that the values of the coeﬃcient B for N t = 8 are slightly larger, then for N t = 10 and , which almost do not diﬀer with each other within errorbars.It may be attributed to ﬁnite lattice spacing eﬀects. For the N t = 8 lattices the dependence of the B on N s is either very slowly rising with N s or constant within the uncertainty. While for the N t = 10 , lattices thereis no such dependence within the uncertainty of the calculation. We thus conclude, that for PBC the relationbetween the critical temperature and the linear boundary velocity also has the form (17) with B ∼ . .• Similarly to Fig. 6, in Fig. 10 we present the Polyakov loop |(cid:104) L ( x, y = 0) (cid:105)| as a function of the coordinate x forthe lattice × × with PBC. The results are shown for two temperatures: T /T c (0) = 0 . in the conﬁnementphase and T /T c (0) = 1 . in the deconﬁnement phase. In addition we plot data for the lattice without rotation Ω I = 0 MeV (Fig. 10(a)) and with Ω I = 24 MeV (Fig. 10(b)). With and without rotation Polyakov loop is zeroin the conﬁnement phase, whereas it develops nonzero values in the deconﬁnement phase. Comparing Fig. 10(a)and Fig. 10(b) it is seen that the Polyakov loop acquires weak dependence on the coordinate due to the rotation.2

C. Dirichlet boundary conditions . . . . . . . T/T c (0)0 . . . . . . . . h | L | i × × , DBCΩ I = 0 MeVΩ I = 10 MeVΩ I = 20 MeVΩ I = 24 MeVΩ I = 30 MeV (a) . . . . . T/T c (0)1 . . . . . . . χ × × , DBCΩ I = 0 MeVΩ I = 10 MeVΩ I = 20 MeVΩ I = 24 MeVΩ I = 30 MeV (b) Figure 11: The Polyakov loop (a) and the Polyakov loop susceptibility (b) as a function of temperature for diﬀerent valuesof imaginary angular velocity Ω I . The results are obtained on the lattice × × with DBC. The lines for the Polyakovloop (a) are drawn to guide the eye. The Polyakov loop susceptibilities (b) are ﬁtted in the vicinity of the phase transition bythe Gaussian function (13). In this section we present the results for the phase diagram of the rotating gluodynamics with DBC. These bound-ary conditions explicitly break Z symmetry what leads to additional lattice artifacts, which disappear only in thethermodynamical limit (see Appendix A). As a consequence of the explicit center symmetry breaking, it is natural toexpect that they make the phase transition smoother. In order observe good peak in the susceptibility of the Polyakovloop one has to conduct lattice simulations on the lattices with much larger spatial volumes as compared to OBC andDBC. For this reason DBC are more expensive from computational point of view and we made an exploratory studyof these BC with few investigated lattice sizes.In Fig. 11 we present the Polyakov loop and its susceptibility as functions of the ratio T /T c (0) for the lattice size × × . One important diﬀerence between DBC and other two boundary conditions is that the Polyakov loopdoes not go to zero in the conﬁnement phase: it is a consequence of the explicit symmetry breaking. Nevertheless,there is a clear inﬂection point for the Polyakov loop, as well as the peak for the susceptibility. Using the standardGaussian ﬁt (13) we determine the critical temperature from the Polyakov loop susceptibility peak, which is shownin Fig. 12(a). In Fig. 12(b) we present, similarly to other BC, the ratio T c ( v I ) /T c (0) as a function of the (imaginary)boundary velocity v I .The behaviour of the critical temperature T c in rotating gluodynamics with DBC is completely analogous to OBCand PBC:• The ratio T c (Ω I ) /T c (0) is described with good accuracy by a function T c (Ω I ) /T c (0) = 1 − C Ω I , with C > .After analytical continuation we draw a conclusion that with DBC the critical temperature of the conﬁne-ment/deconﬁnement transition increases with angular velocity .• The coeﬃcient C has very weak dependence on the lattice spacing and lattice size N z but signiﬁcantly changeswith the size N s . If one takes instead Ω I the linear boundary velocity v I : T c ( v I ) /T c (0) = 1 − B v I , theratio T c ( v I ) /T c (0) exhibits signiﬁcantly smaller dependence on the lattice size N s (see Fig. 12(b)). Looking atFig. 12(b) one can see that due to large uncertainties the divergence of lines with diﬀerent N s is quite signiﬁcant.The coeﬃcient B is presented in Fig. 13. From this ﬁgure it is seen that similarly to OBC, for suﬃciently large N s the contribution of the boundary is suppressed and the B goes to plateau with the value B ∼ . .• In Fig. 14 we show the dependence of the Polyakov loop on the spatial coordinate for the angular velocities Ω I = 0 MeV , MeV. The results are shown for two temperatures:

T /T c (0) = 0 . in the conﬁnement phaseand T /T c (0) = 1 . in the deconﬁnement phase.DBC ﬁx the value of the Polyakov loop L ( x, y ) = 3 on the boundary. Similarly to OBC, the eﬀect of boundaryconditions is screened and at suﬃciently large distances from the boundary the Polyakov loop tends to a constant,3 I (MeV )0 . . . . T c ( Ω I ) / T c ( ) N s /N t ’ N s /N t ’ N s /N t ’ DBC8 × × × × × × × × × × × × × × (a) .

00 0 .

02 0 .

04 0 .

06 0 .

08 0 .

10 0 .

12 0 . v I /c . . . . T c ( v I ) / T c ( ) N s /N t ’ N s /N t ’ N s /N t ’ DBC8 × × × × × × × × × × × × × × (b) Figure 12: The ratios T c /T c (0) determined from the peak of the Polyakov loop susceptibility as a function of the imaginaryangular velocity squared Ω I (a) and the linear boundary velocity squared v I (b). Results are presented for several lattice sizes N t × N z × N s with DBC. Lines correspond to simple quadratic ﬁts T c (Ω I ) /T c (0) = 1 − C Ω I and T c ( v I ) /T c (0) = 1 − B v I /c which is zero in the conﬁnement phase and nonzero for the deconﬁnement. Comparing Fig. 14(a) and Fig. 14(b)one can notice that rotation induces additional inhomogeneity of Polyakov loop but its inﬂuence is much weakeras compared to that due to the BC. . . . . . . . N s /N t . . . . . . B DBC8 × × N s × × N s × × N s × × N s Figure 13: The coeﬃcient B in Eq. (16) versus the ratio N s /N t for several lattice sizes with DBC. IV. DISCUSSION AND CONCLUSION

In this paper we addressed the question how rotation inﬂuences the conﬁnement/deconﬁnement transition in gluo-dynamics within lattice simulation. We preform the simulation in the reference frame which rotates with the systemunder investigation. In this reference frame rotation is reduced to external gravitational ﬁeld. Having constructed theaction of lattice gluodynamics in external gravitational ﬁeld we found that this action is spoiled by sign problem anddirect application of Monte Carlo importance sampling is not possible. To overcome this problem we conducted ourstudy for suﬃciently small imaginary angular velocities Ω I and the results were analytically continued to real valuesof angular velocity Ω . Our results suggest that this approach is applicable for the values of angular velocity which are4 . . − −

10 0 10 20 x . . . . . . . | h L ( x , y = ) i | DBC, Ω I = 0 MeV T/T c (0) = 0 . T/T c (0) = 1 . (a) . . − −

10 0 10 20 x . . . . . . . | h L ( x , y = ) i | DBC, Ω I = 24 MeV T/T c (0) = 0 . T/T c (0) = 1 . (b) Figure 14: The Polyakov loop |(cid:104) L ( x, y = 0) (cid:105)| as a function of coordinate x for DBC and Ω I = 0 MeV (a), Ω I = 24 MeV (b).The results were obtained on the lattice × × for two temperatures: T /T c (0) = 0 . in the conﬁnement phase and T /T c (0) = 1 . in the deconﬁnement phase. characteristic for heavy-ion collision experiments.Because of the causality, the simulation of rotating gluodynamics has to be carried out with boundary conditions(BC). It is important to stress that BC are important part of all approaches aimed at studying of rotating quark-gluonplasma rather than a lattice artifact. The results obtained within any approach depend on BC. In order to study thisdependence in our paper we applied various BC. In particular, our simulations were carried out with open boundaryconditions (OBC), periodic boundary conditions (PBC), Dirichlet boundary conditions (DBC). In our paper we aremainly focused on the inﬂuence of rotation on the critical temperature of the conﬁnement/deconﬁnement transition.The critical temperature was determined from the peak of Polyakov loop susceptibility.The results obtained in our work allow us to state that after the analytical continuation the T c (Ω) can be welldescribed by a simple quadratic function T c (Ω) T c (0) = 1 + C Ω , (19)with C > for all BC and all lattice parameters used in the simulations. From this result we draw the main conclusionof our paper the critical temperature of the conﬁnement/deconﬁnement phase transition grows with increasing angularvelocity. This conclusion does not depend on BC and we believe that this is universal property of gluodynamics.The magnitude of the coeﬃcient C depends on BC. For each boundary condition used in the simulations, the C does not depend on the lattice size along the rotation axis N z , has weak dependence on the lattice spacing, but ithas strong dependence on the lattice size perpendicular to the rotation axis N s = N x = N y . The leading dependenceof formula (19) on the N s can be captured if one rewrites it in terms of the linear velocity v on the boundary v = Ω ( N s − a/ as follows T c ( v ) T c (0) = 1 + B v c . (20)In last formula the coeﬃcient B weakly depends on N s . We believe the possibility to describe our results for all N s by universal formula (20) rather than formula (19) originates from the following fact. For the thermodynamics ofrotating system the most important physical object is the ﬁeld of velocities and how close this ﬁeld approaches thespeed of light. This property is determined by the product Ω a ( N s − / but not the Ω alone. In addition we alsofound that for OBC B ∼ . , for PBC B ∼ . and for DBC B ∼ . . So, although the values of the B are closeto each other, we still see the dependence of our quantitative results on BC.In paper [20] it was proposed that rotation might lead to mixed inhomogeneous phase when the matter is inthe conﬁnement phase close to the center of rotation whereas the deconﬁnement takes place close to the boundary.Unfortunately for lattice parameters used in our study we have not found such state. Probably a more thorough studyis required to answer the question whether such a phase is realized in rotating gluodynamics.5The authors of papers [20–22] studied the conﬁnement/deconﬁnement transition in rotating QCD within phe-nomenological models. The results of these works indicate that the critical temperature decreases with the angularvelocity, which is in disagreement with the results of our work. The origin of this disagreement is not yet clear, in thefuture we plan further investigation of this problem. Acknowledgments

We would like to thank Valentin Zakharov, Oleg Teryev, Maxim Chernodub, Vitaly Bornyakov for useful discussions.This work was supported by RFBR grants 18-02-40126. This work has been carried out using computing resources ofthe Federal collective usage center Complex for Simulation and Data Processing for Mega-science Facilities at NRC“Kurchatov Institute”, http://ckp.nrcki.ru/ ; the cluster of Institute for Theoretical and Experimental Physics andthe Supercomputer “Govorun” of Joint Institute for Nuclear Research.

Appendix A: Finite volume eﬀects on the lattices with Dirichlet and open boundary conditions

In this section we are going to address the question how boundary conditions (BC) considered in this paper inﬂuencethe critical temperature for gluodynamics without rotation. For all BC used in our paper we apply periodic boundaryconditions for gluon ﬁelds in z - and t -directions. What concerns the x - and y -directions we used periodic (PBC),Dirichlet (DBC) and open (OBC) boundary conditions (see Section II B for details). Since PBC are common forlattice simulations, there are a lot of papers where the volume dependence of critical temperature within PBC wasstudied (see, for instance, [30]). For this reason, in Appendix we are going to focus on DBC and OBC only.

1. Open boundary conditions

In this section we consider OBC. The critical temperature is determined from the peak of the Polyakov loopsusceptibility (12). The results of this calculation for diﬀerent lattices can be found in Table I. It is seen that forall lattices presented in this table the critical temperature for OBC is larger than that in gluodynamics with PBC T c / √ σ = 0 . [30]. This feature of OBC can be understood as follows. In OBC the plaquettes outside thelattice volume are excluded. This can be considered as if one puts these plaquettes to unity what leads to zero latticeaction. So, physically this can be interpreted as if the studied volume is attached to classical (without quantumﬂuctuations) zero temperature Yang-Mills theory. From this perspective the regions near the boundary have lowertemperature than the regions remote from the boundary. So, in order to observe conﬁnement/deconﬁnement transitionwith OBC one has to heat the system to larger temperature as compared to homogeneous gluodynamics.To study the inﬁnite volume limit we ﬁt our data for the lattices × N s by the function T c ( N s /N t ) = T + A ( N t /N s ) (A1)where T corresponds to the inﬁnite volume limit. The ﬁt gives T / √ σ = 0 . . This value is in reasonableagreement with that obtained in the inﬁnite volume limit for the PBC: T / √ σ = 0 . [30].From these results one can draw a conclusion that ﬁnite volume eﬀects for OBC enhance the critical temperature.As one increases the volume, the eﬀect of BC becomes screened and the critical temperature for OBC approaches toits acknowledged value [30]. The screening of the boundary is well seen in Fig. 6(a). Notice that in this and next section we compare our results for the critical temperatures with the result obtained in the inﬁnite volumelimit on the lattices with PBC and N t = 8 [30]. This is because in Appendix we do not take continuum limit and conduct the simulationson N t = 8 lattices. Table I: The critical temperature for non-rotating lattices with OBC and DBC. These results to be compared with the criticaltemperature of gluodynamics with PBC: T c / √ σ = 0 . [30].OBC DBCLattice T c / √ σ T c / √ σ × × × ×

2. Dirichlet boundary conditions

In this section we consider DBC. The results for the critical temperature calculation with DBC for diﬀerent latticescan be found in Table I. It is seen that for all lattices presented in this table the critical temperature for DBCis smaller than that in gluodynamics with PBC [30]. In this sense OBC and DBC inﬂuence to the system in theopposite way. The decrease of the critical temperature in DBC as compared to PBC gluodynamics can be understoodas follows. The link variables on the boundary in DBC equal unity, i.e. after taking the trace over colors thePolyakov loop on the boundary equals 3. So, the Z symmetry of lattice gluodynamics is explicitly broken on theboundary. In gluodynamics breaking of Z symmetry is the property of high temperature phase. For this reason inDBC the boundary can be considered as a seed of high temperature phase. On these grounds, one can expect thatthe conﬁnement/deconﬁnement transition takes place at smaller critical temperature.To study the inﬁnite volume limit we ﬁt our data for the lattices × N s by the same function (A1),where T corresponds to the inﬁnite volume limit. The ﬁt gives T / √ σ = 0 . . This value is in reasonableagreement with that obtained in the inﬁnite volume limit for PBC: T / √ σ = 0 . [30].So, ﬁnite volume eﬀects in DBC decrease the critical temperature which approaches to that of PBC gluodynam-ics [30] for suﬃciently large volume. Similarly to OBC this behaviour can be explained by screening of the boundarywhich is well seen on Fig. 14(a). From this ﬁgure one can note that the Polyakov loop increases as one approaches tothe boundary. This can be interpreted as increase of eﬀective temperature of the regions close to the boundary.To summarize of both sections of Appendix, our results allows us to state that for suﬃciently large volume bulkproperties of the system cannot be considerably modiﬁed by either OBC or DBC. We believe that this conclusionremains to be true even for rotating gluodynamics if the angular velocity is not too large. Appendix B: Lattice parameters used in the simulations

In our paper we performed numerical simulations for the lattice sizes and values of imaginary angular velocitylisted in Tab. II. For each lattice size and angular velocity we changed the temperature through the variation of the β . To set the physical scale we used the relation between lattice spacing and inverse lattice coupling ( a √ σ )( β ) fromRef. [31] (with the value of string tension √ σ = 440 MeV). Simulations are performed with the use of Monte Carloalgorithm, each sweep consists of one heatbath update and two steps of the overrelaxation updates. For each set ofparameters, the typical statistics are about 6000–12000 conﬁgurations, separated by 20 sweeps. The statistical errorswere estimated using jackknife method.7

Table II: The list of lattice parameters used in the simulations with OBC, PBC and DBC .OBC PBC DBCLattice Ω I , MeV Lattice Ω I , MeV Lattice Ω I , MeV × ×

0, 30, 45, 60, 75, 90 × ×

0, 15, 30, 45, 50 × ×

0, 15, 30, 45 × ×

0, 30, 45, 60, 75, 90 × ×

0, 10, 20, 30, 38 × ×

0, 10, 20, 30, 36 × ×

0, 30, 45, 60, 75, 90 × ×

0, 12, 20, 24, 30 × ×

0, 12, 24, 30, 36 × ×

0, 15, 30, 45, 60, 75 × ×

0, 10, 15, 20, 24 × ×

0, 12, 24, 30, 36 × ×

0, 15, 30, 45, 50 × ×

0, 15, 30, 45, 50 × ×

0, 10, 20, 24, 30 × ×

0, 12, 24, 30, 36 × ×

0, 15, 30, 45, 50 × ×

0, 10, 15, 20, 24 × ×

0, 12, 24, 30, 36 × ×

0, 15, 30, 45, 50 × ×

0, 15, 30, 45 × ×

0, 10, 20, 24, 30 × ×

0, 10, 20, 30, 38 × ×

0, 12, 24, 30, 36 × ×

0, 10, 15, 20, 24 × ×

0, 12, 20, 24, 30 × ×

0, 30, 45, 60, 75, 90 × ×

0, 15, 30, 45, 50 × ×

0, 30, 45, 60, 90 × ×

0, 10, 20, 30, 38 × ×

0, 15, 30, 45, 50 × ×

0, 12, 24, 30, 36 × ×

0, 30, 45, 60, 75, 90[1] A. L. Watts et al. , Rev. Mod. Phys. , 021001 (2016), arXiv:1602.01081 [astro-ph.HE] .[2] I. A. Grenier and A. K. Harding, Comptes Rendus Physique , 641 (2015), arXiv:1509.08823 [astro-ph.HE] .[3] G. Basar, D. E. Kharzeev, and H.-U. Yee, Phys. Rev. B , 035142 (2014), arXiv:1305.6338 [hep-th] .[4] K. Landsteiner, Phys. Rev. B , 075124 (2014), arXiv:1306.4932 [hep-th] .[5] Y. Jiang, Z.-W. Lin, and J. Liao, Phys. Rev. C , 044910 (2016), [Erratum: Phys.Rev.C 95, 049904 (2017)],arXiv:1602.06580 [hep-ph] .[6] F. Becattini, F. Piccinini, and J. Rizzo, Phys. Rev. C , 024906 (2008), arXiv:0711.1253 [nucl-th] .[7] M. Baznat, K. Gudima, A. Sorin, and O. Teryaev, Phys. Rev. C , 061901 (2013), arXiv:1301.7003 [nucl-th] .[8] L. Adamczyk et al. (STAR), Nature , 62 (2017), arXiv:1701.06657 [nucl-ex] .[9] A. Vilenkin, Phys. Rev. D , 1807 (1979).[10] D. Kharzeev, J. Liao, S. Voloshin, and G. Wang, Prog. Part. Nucl. Phys. , 1 (2016), arXiv:1511.04050 [hep-ph] .[11] G. Prokhorov, O. Teryaev, and V. Zakharov, Phys. Rev. D , 071901 (2018), arXiv:1805.12029 [hep-th] .[12] G. Y. Prokhorov, O. V. Teryaev, and V. I. Zakharov, JHEP , 146 (2019), arXiv:1807.03584 [hep-th] .[13] O. Rogachevsky, A. Sorin, and O. Teryaev, Phys. Rev. C , 054910 (2010), arXiv:1006.1331 [hep-ph] .[14] O. V. Teryaev and V. I. Zakharov, Phys. Rev. D , 096023 (2017).[15] S. Ebihara, K. Fukushima, and K. Mameda, Phys. Lett. B , 94 (2017), arXiv:1608.00336 [hep-ph] .[16] M. Chernodub and S. Gongyo, JHEP , 136 (2017), arXiv:1611.02598 [hep-th] .[17] Y. Jiang and J. Liao, Phys. Rev. Lett. , 192302 (2016), arXiv:1606.03808 [hep-ph] .[18] H. Zhang, D. Hou, and J. Liao, Chin. Phys. C , 111001 (2020), arXiv:1812.11787 [hep-ph] .[19] X. Wang, M. Wei, Z. Li, and M. Huang, Phys. Rev. D , 016018 (2019), arXiv:1808.01931 [hep-ph] .[20] M. Chernodub, (2020), arXiv:2012.04924 [hep-ph] .[21] X. Chen, L. Zhang, D. Li, D. Hou, and M. Huang, (2020), arXiv:2010.14478 [hep-ph] .[22] Y. Fujimoto, K. Fukushima, and Y. Hidaka, (2021), arXiv:2101.09173 [hep-ph] .[23] Y. Nambu and G. Jona-Lasinio, Phys. Rev. , 345 (1961).[24] Y. Nambu and G. Jona-Lasinio, Phys. Rev. , 246 (1961).[25] A. Yamamoto and Y. Hirono, Phys. Rev. Lett. , 081601 (2013), arXiv:1303.6292 [hep-lat] .[26] V. Braguta, A. Kotov, D. Kuznedelev, and A. Roenko, JETP Lett. , 6 (2020).[27] M. Luscher and S. Schaefer, JHEP , 036 (2011), arXiv:1105.4749 [hep-lat] .[28] M. Lüscher, JHEP , 105 (2014), arXiv:1404.5930 [hep-lat] .[29] A. Florio, O. Kaczmarek, and L. Mazur, Eur. Phys. J. C , 1039 (2019), arXiv:1903.02894 [hep-lat] .[30] G. Boyd, J. Engels, F. Karsch, E. Laermann, C. Legeland, M. Lutgemeier, and B. Petersson, Nucl. Phys. B , 419(1996), arXiv:hep-lat/9602007 .[31] R. Edwards, U. M. Heller, and T. Klassen, Nucl. Phys. B517