Rigidly-rotating scalar fields: Between real divergence and imaginary fractalization

The thermodynamics of rigidly rotating systems experience divergences when the system dimensions transverse to the rotation axis exceed the critical size imposed by the causality constraint. The rotation with imaginary angular frequency, suitable for numerical lattice simulations in Euclidean imaginary-time formalism, experiences fractalization of thermodynamics in the thermodynamic limit, when the system’s pressure becomes a fractal function of the rotation frequency. Our work connects these two phenomena by studying how thermodynamics fractalizes as the system size grows. We examine an analytically-accessible system of rotating massless scalar matter on a one-dimensional ring and the numerically-treatable case of rotation in the cylindrical geometry and show how the ninionic deformation of statistics emerges in these systems. We discuss a no-go theorem on analytical continuation between real-and imaginary-rotating theories. Finally, we compute the moment of inertia and shape deformation coefficients caused by the rotation of the relativistic bosonic gas. In all cases, we show that finite-mass effects are quantitative, leaving our conclusions qualitatively unchanged.


I. INTRODUCTION
Effects of rotation on the state of physical bodies have been a subject of passionate interest throughout the decades.In metals, the uniform rotation acts on electrons via a centrifugal force that produces a slight but experimentally perceptible gradient of electric potential measured at ∼ 10 2−3 Hz [1].At the level of electronic spins, one of the numerous examples of rotation-generated phenomena is the Barnett effect [2] which -with its celebrity reciprocal, the Einstein-de Haas effect [3] -relates the mechanical torque and magnetization in ferromagnets.The nuclear analog of the Barnett effect substantially affects the polarization of the protons (ions of hydrogen) in the water rotating with the frequency ∼ 10 4 Hz [4].
However, the fastest rotation of matter has been produced in noncentral collisions of relativistic heavy ions that create quark-gluon plasma in which the vorticity reaches the values ∼ 10 22 Hz [5][6][7].The fast rotation affects the local properties of quark-gluon plasma, leading to various spin polarization phenomena, allowing us to probe experimentally the interior of rapidly rotating plasma in terms of its local vortical structure [8,9].
There are various theoretical shreds of evidence that fast rotation also affects the chiral [10][11][12][13][14][15][16] and (de)confining transitions [17][18][19][20][21][22][23][24][25][26][27] of the quark-gluon plasma.Theoretical methods, however, prevailingly assume a rigid rotation that makes every physical point rotate about a fixed axis with the same angular velocity.While the rigid character of rotation substantially simplifies the analytical treatment of the problem [28,29], the consensus on the thermodynamic properties of quarkgluon plasma, even in this simplest case, is still absent, thus opening a gap between numerical and various analytical calculations.Moreover, the latest first-principle simulation reveals the instability of the rigidly rotating gluon plasma below the "supervortical" critical temperature [27], indicating the complexity of rotation in strongly interacting systems.
First-principle information about the quark-gluon plasma comes from lattice simulation in the Euclidean imaginary time formalism where the real angular momentum Ω brings the sign problem [30] which makes the numerical simulations impossible.This inconvenience can traditionally be overcome by turning the angular momentum into the complex plane and considering the purely imaginary rotation Ω I = −iΩ in full analogy with the baryon chemical potential [17,19,21,23,30,31].
The imaginary rotation differs, however, from the imaginary baryonic chemical potential: the analytical continuation to real rotation used in numerical lattice simulations has some unusual features including the emergence of (stable) ghost-like excitations [19] characterized by "ninionic" deformation of statistics and the appearance of the fractal features of thermodynamics under imaginary rotation.The fractalization imposes a no-go theorem on an analytical continuation for rotating systems in the thermodynamic limit [32].
In our work, we discuss the effect of rotation on the thermodynamics of the simplest possible system represented by the massless scalar fields.This limit is convenient as it allows many of the calculations to be performed analytically.We explore the effect of finite mass in all considered setups and confirm that all features remain qualitatively identical to those in the massless case.
First, we briefly introduce the real and imaginary rotation in Sec.II.Then, in Sec.III, we analyze the interrelation between the fractal features, analytical continuation, and the causality constraint for the model formulated on a one-dimensional ring which can be treated analytically.
Section IV approaches real and imaginary rotation within the scope of the relativistic kinetic theory applied to the three-dimensional rotating gas.It also addresses the mechanical features of the rotating gas, including its moment of inertia and shape-deformation coefficients.This analysis is followed by Sec.V, where we pursue, for simplicity, a "hybrid" quantization approach based on the cylindrical waves with continuous momentum in a spatially unbounded region.We show the advantages of both discussed approaches in Sec.VI, where the rotating gas in the cylindricallybounded region and the discrete quantization of the transverse modes is treated numerically in great detail.Furthermore, we reveal the analytical fractalization of the thermodynamics numerically in the threedimensional rotating gas and explicitly show strong parallels with the fractalization of thermodynamics in the analytically-accessible one-dimensional ring.For simplicity, we implement the boundary using Dirichlet boundary conditions [33], however it is reasonable to expect that similar results regarding fractalization can be expected when von Neumann or, more generally, Robin boundary conditions are considered [34][35][36].
Our last section is devoted to conclusions and a summary of our results.Throughout the article, we work with the conventions ℏ = c = k B = 1.

A. Real rotation and imaginary rotation
Let us consider a quantum-mechanical system of bosonic particles rotating uniformly (rigidly, as a solid body) with the constant angular frequency Ω about the z axis.For simplicity, we can assume that the system of particles is rotating inside a cylinder possessing reflective boundary conditions, such as Dirichlet boundary conditions [33].In the co-rotating frame, the free energy of the system takes the following form: where β = 1/T is the inverse temperature, V is the volume of the system, ω α,m is the energy spectrum of the particles in the laboratory reference frame, and α is a collective notation of quantum numbers other than the projection of angular momentum, m ≡ m z ∈ Z.We work with zero-charge systems so that the chemical potential does not enter the free energy of the system (1).We also ignore the zero-point contribution associated with the vacuum Casimir energy since it does not affect the thermodynamics of the system.In order to determine the thermodynamic characteristics (for example, energy, pressure, entropy, angular momentum, moment of inertia, etc), it is sufficient to evaluate the statistical integral (1).For bosonic particles, the contribution of each quantum level to the thermodynamic quantities is given by the Bose-Einstein distribution where ω is the energy of the quantum level.In a rigidly rotating system, the statistical weight is determined by the energy in the co-rotating reference frame: thus demonstrating explicitly how rotation with Ω ̸ = 0 affects the statistical particle distribution.
It is convenient to calculate the thermodynamic properties of rotating particles using the imaginary-time formalism in which the time coordinate is turned to a complex variable via the Wick transformation, t → τ = it.The imaginary time τ is compactified to a circle of length β = 1/T related to thermal equilibrium temperature T .The compactification imposes the matching conditions on the fields: all scalar fields ϕ are periodic functions along the thermal direction, while all fermionic fields (not considered in this article) obey anti-periodic boundary conditions.The Bose-Einstein statistical distribution (2) for bosonic fields can be recovered automatically from the periodic boundary conditions (4) [37].The imaginary-time approach is intensively used in numerical lattice simulations of quantum field theories where the partition function is formulated in terms of a statistical integral in Euclidean spacetime [38].The lattice simulations are especially useful for obtaining information about non-perturbative effects that cannot be treated with standard perturbative methods [38].
However, the imaginary-time techniques cannot be directly applied to rotating systems because the action of the Euclidean theory becomes a complex quantity at a nonzero angular frequency, Ω ̸ = 0, thus exhibiting the so-called "sign problem" [30].The latter property does not allow us to treat the partition function of a rotating system as a statistical integral, bringing us to an inconvenient similarity with finite-density systems, where the (baryonic) chemical potential also makes the Euclidean action a complex quantity [38].The only practical way to avoid the sign problem for rotation is to consider the angular frequency as a purely imaginary variable: Ω = iΩ I . ( The shift of the angular frequency to the complex plane (5) restores the real-valuedness of the Euclidean action [17,30].Having calculated the desired quantities at a set of imaginary Ω I , one can then apply an analytical continuation to map the results back to the realistic case of real rotation [17,21].This prescription, applied to the angular frequency Ω, follows a standard set of practices invoked to avoid the sign problem in simulations of finite-density systems [39].
In the context of the imaginary-time formalism, there are two methods how one can implement the imaginary rotation (5).The first approach, originally proposed in Ref. [30] and adopted in various numerical Monte Carlo simulations of (quark-) gluon plasmas [21,26,27,30], consists in (i) considering the system in a non-inertial co-rotating reference frame in Minkowski spacetime; (ii) turning the system, via a Wick transformation, to the curved Euclidean spacetime with a complex metric tensor; (iii) implementing the substitution (5) which makes the metric tensor real-valued again; (iv) simulating the thermodynamics at a set of non-zero Ω I with the standard periodic boundary conditions (4); (v) fitting the obtained numerical results by a reasonable analytical function and, finally, (vi) making an analytical continuation of the lattice results to the real-valued frequency by setting The second approach implements the imaginary rotation in the imaginary-time formalism in a more straightforward way using the property that the imaginary frequency Ω I corresponds, after all, to a uniform rotation of a subspace of a timeslice of the Euclidean spacetime about a certain fixed axis [19,23].As the imaginary time variable τ advances for a full period from τ = 0 to τ = β, the system experiences a spatial rotation by the angle: The turn of the space necessitates a modification of the standard bosonic boundary conditions (4) which should now incorporate a translation in imaginary time with the uniform rotation of the Euclidean spacetime.Under the imaginary rotation, the bosonic wavefunction appears to satisfy the rotwisted boundary condition: where the 3 × 3 matrix written in Cartesian coordinates, corresponds to the global rotation of the whole spatial Euclidean subspace, x → x ′ = Rχ x, by the angle (7) around the z axis.In the absence of rotation, the transformation (9) becomes a unit matrix and the boundary condition (8) reduces to the standard periodic condition for bosons (4).The rotwisted boundary conditions, visualized in Fig. 1, have already been discussed in the context of the Euclidean lattice simulations of field theories [23,26].The boundary conditions ( 8) are obviously invariant under 2π shifts of χ, or equivalently, shifts by one unit in ν: and, in the parity-unbroken systems, under the reversal of the rotation angle: The latter condition holds for a system of neutral particles that we consider.The symmetry under clockwise and counterclockwise rotations (11) can be broken, for example, for charged particles subjected to a background magnetic field which leads, in particular, to a rotation diode effect in semiconductors [40].

B. Imaginary rotation and ninionic statistics
The differences between the two implementations of the imaginary rotation are two-fold: one can either consider the curved Euclidean spacetime with corotating coordinates and ordinary boundary conditions (4) implemented along the compactified time (the first approach) or use Cartesian coordinates with the rotwisted boundary conditions (8) following the second approach.In our article, we consider field theories subjected to imaginary rotation introduced via the rotwisted boundary condition.
The boundary conditions imposed on the fields in the imaginary time direction have one-to-one correspondence with the statistical distribution of the particles.For example, the equal-time commutation relations for bosonic fields imply the periodic boundary conditions (4), which lead to the Bose-Einstein distribution for bosons (2).Analogously, anti-commuting fermionic variables possess anti-periodic conditions in the compactified time direction which correspond to the Fermi-Dirac statistics [37].
Therefore, it is appropriate to ask which statistical distribution corresponds to the rotwisted boundary conditions (8)?
It turns out that the imaginary rotation deforms the statistical distribution of fermions and bosons, leading to a "ninionic" deformation which matches neither the bosonic nor the fermionic statistical distributions [32].For example, for bosons, the ninionic deformation takes the following form: where ω ≡ ω α,m is associated with the energy of the quantum state in the laboratory reference frame and ξ = mχ = 2πmν is the deformation parameter associated with the "statistical angle" χ, Eq. ( 7).The latter depends on the angular velocity Ω I of the imaginary rotation in Euclidean spacetime.Notice that at the zero (modulo 2π) statistical angle, the ninionic deformation ( 12) of the bosonic distribution (2) disappears: ω with an integer k ∈ Z.The ninionic deformation ( 12) can be understood as the real part of the bosonic occupation number (2), at an imaginary chemical potential µ = ξ/β.Given the unusual form of the ninionic deformation of the bosonic statistical distribution (12), it is appropriate to ask how this deformation modifies the statistical properties of the thermal state?What are the consequences which are brought to the theory by the introduction of the new dimensionless parameter, the statistical angle ( 7)?The answer to this question, which depends on the volume of the rotating system, is one of the aims of our article.
In respect of the causality, the rigid rotation with real-valued angular velocity is a well-defined notion only for transversely-bounded systems.On the contrary, the imaginary rotation does not impose any bounds on the size of the system due to the absence of the notion of the light cone in the Euclidean space (in other words, there is no causality constraint in the imaginary time formalism because it has no notion of real time).Therefore, the imaginary rotation does not lead to causality problems [30] and can be formulated in the thermodynamic limit in the whole Euclidean space [23].The relation between imaginary and real rotation in terms of the analytical continuation is another aim of our paper.

C. Ninionic statistics and fractal thermodynamics
Sticking to an infinite-volume system, one can show, both in the scope of a classical interacting field theory [31] as well as in a free bosonic quantum field theory [32], that the imaginary rotation characterized by the nonvanishing value of the statistical angle χ modifies the relation between the physical temperature T → T I (β, χ) ≡ T (β, iβΩ I ) and the length of the compactified direction β: where is the Thomae function.In other words, function (15) gives zero for all irrational numbers and equals to a nonzero number 1/q determined by the denominator q of the rational argument x = p/q ∈ Q with two natural coprime numbers p, q ∈ N.
The Thomae function (15), shown in Fig. 2, is known also under other names, such as the raindrop function, the modified Dirichlet function, the popcorn function, etc.This function has the amazing counter-intuitive property stating that the function is discontinuous if its argument x is rational, and it is continuous provided x is irrational.The Thomae function emerges naturally in various physical systems such as one-dimensional disordered crystals, Hubbard model of particles on a ring, phyllotaxis problem, and fractional quantum Hall effect (we refer to Ref. [41] for a comprehensive discussion).The Thomae function (15) possesses a nontrivial fractal structure [41][42][43] which equips the thermodynamics of imaginary rotation with fractal properties.The fractalization (and "defractalization") of thermodynamics of imaginary rotating systems will also be addressed in this paper.
Notice that the behavior of the physical temperature (14) as a function of the statistical angle ( 7) is determined solely by the denominator q of the rational number χ/(2π) and not by its numerator.Irrational (in units of 2π/β) frequencies correspond to zero temperature (14).
In the absence of the imaginary rotation, χ = 0, one gets the standard relation between temperature T and the length of the imaginary time direction β, as expected: The ninionic deformation of bosonic statistics can be readily understood in the imaginary time formalism for free massless bosons in the thermodynamic limit (on an infinite spatial line) where the particles possess the linear energy dispersion, ω k = |k|.In this conformal system, the thermal pressure of bosons P is equal to their energy density, E ≡ P , taking a well-known expression in the absence of imaginary rotation: The temperature of the system is given by the inverse length β of the compactified imaginary-time direction (16).Looking ahead a little, one can also discuss the thermodynamics of the same system compactified into a ring of an infinitely large radius R which is subjected to the imaginary rotation with the angular velocity Ω I .The pressure can be derived via the ninionic statistics (13): where the ninionic parameter ξ m = χm ≡ βΩ I m is expressed via the angular momentum m and the statistical angle (7).The energy spectrum in the statistical sum (18), corresponds to the laboratory frame.In the thermodynamic limit, R → ∞, the energy gaps of the discrete spectrum (19) shrink, the variable m/R becomes the continuum momentum k, and the sum in Eq. ( 18) reduces to an integral thus bridging the gap between the exotic (18) and standard (17) statistical sums in the thermodynamic limit.
However, the presence of the imaginary rotation Ω I makes the system nontrivial even in the thermodynamic limit.Indeed, the pressure of the imaginary rotating system (18) has the same expression as the pressure of the non-rotating one (17) with only one important difference that the temperature of the former ( 14) becomes a fractal function (15) of the imaginary angular frequency Ω I .In the next section, we discuss the particularities of fractalization of thermodynamics by imaginary rotation working with an analytically-solvable example of a free massless particle confined to a one-dimensional ring.

III. REAL AND IMAGINARY ROTATIONS ON THE RING
A. Classical (non-quantum) system Let us consider an ensemble of bosonic massless particles, constrained to move on a ring of fixed radius R with angle coordinates φ.Imposing that the system rotates with angular velocity Ω, its free energy at inverse temperature β reads where J z = −k φ represents the angular momentum of a particle spinning with azimuthal momentum The above integral can be readily evaluated, Other thermodynamic quantities such as total entropy S, thermodynamic pressure P, and total angular momentum M, can be obtained starting from the relation [44] where V = 2πR is the "volume" of the ring.Denoting by an overhead bar volume-averaged quantities, A ≡ A/V , we have while the average energy E can be found from the Euler relation, Taking into account the thermodynamic relations in Eq. ( 23), the above relation can be recast as Starting from Eq. ( 21), it is easy to derive 3. Effect of particle mass M on free energy and shape coefficients K2 and K4, computed in relativistic kinetic theory.The blue line shows the ratio F(0, M )/F0.with F0 = −P0 being the average free energy at vanishing mass [cf.Eq. ( 17)].
while E = P for this one-dimensional system.It can be checked that E agrees with the energy density E computed using the well-known formula Considering now the free energy at small values of the rotation parameter Ω, we expand the free energy density in a power series in the velocity of a corotating particle at the system boundary ρ = R, following Ref.[27]: where K 2n are Ω-independent dimensionless coefficients, F(0) is the free energy density in the absence of rotation, and we took into account that F is an even function of Ω.By construction, one gets K 0 = 1.Comparing Eqs.(29) and (26), it can be seen that Clearly, F(0) and the K 2n coefficients are independent of the system size.
The results for the coefficients K 2n are also valid in a multicomponent non-interacting gas since they reflect the rotational response normalized per degree of freedom.The zero-rotation limit of the moment of inertia, I 0 of a one-component bosonic gas evaluates to which follows from the definition M(Ω) = I(Ω)Ω of the moment of inertia I(Ω) in terms of the average angular momentum M, obtained via the thermodynamic relation (23).
It is interesting to notice that recent first-principle simulations [27] indicate that in the high-temperature limit, the rotating gluon gas possesses the dimensionless moment of inertia K 2 = 2.23 (39) consistent with our estimation (30): This result is not unexpected since at sufficiently high temperatures, the gluon plasma becomes a weaklyinteracting gas of gluons.It is worth mentioning that the value of K 2 receives corrections only if the quanta are massive or if interactions are taken into account.
The next non-zero coefficient in the series (30), corresponds to the correction to the free energy caused by the deformation of the rotating gas due to rotation.This correction also affects the moment of inertia, I(Ω) = I 0 + I 2 v 2 R /2 + . . ., with the universal noninteracting coefficient I 2 /I 0 = 4.
We now consider the effect of the particle mass on the K 2n coefficients.To this end, we evaluate F from Eq. (20) for the case when k 0 = k 2 φ /R 2 + M 2 , with M being the particle mass.Performing a Lorentz transformation to the local rest frame, we get Changing the integration variable to z = k 0 /m and expanding the logarithm using ln(1 − x) = − ∞ j=1 x j /j, we arrive at We can therefore evaluate the coefficients as where we introduced the notation K 2n (M ) ≡ F(0, M )K 2n (M ).As shown in Fig. 3, increasing M decreases the average free energy F(0, M ), while the shape coefficients K 2 (M ) and K 4 (M ) increase with respect to their massless limits.

B. Relativistic rotation, particle spectrum
In this section, we consider a free massless particle on a ring of a fixed radius R with the angle coordinate φ as shown in Fig. 4. For a static ring, the particle wavefunction is described by the Klein-Gordon equation: which is formulated in the inertial, laboratory frame.Let us consider the ring rotating with the constant angular velocity Ω.The coordinates associated with the co-rotating reference frame (denoted by a tilde) are related to the laboratory coordinates as follows: In the co-rotating frame, the Klein-Gordon equation (37) transforms into the following equation: which possesses the energy spectrum in the rotating reference frame: corresponding to the following eigenfunctions: The energy spectrum (40) is bounded from below provided the causality condition is satisfied: The thermodynamics of the system is determined in the rotating reference frame where all statistical distributions are set by the energy in the co-rotating frame ωm rather than by its laboratory-frame counterpart (19).

C. Free energy of rotating scalar field
We consider statistical mechanics of scalar particles in the rotating environment.The corresponding statistical sum, is formulated via the sum over states labeled by the angular momentum m with the occupation number n m of system's levels that possess the total energy Ẽm,nm = ωm n m in the rotating reference frame and the total angular momentum L n,mn = mn m .In Eq. ( 43), F stands for the free energy in the co-rotating reference frame.In the statistical sum, we do not take into account a m = 0 contribution which corresponds to the zero-energy mode and contributes to the zero-point (Casimir) vacuum energy [45].We concentrate on the thermal part of the free energy which possesses interesting fractal properties in the thermodynamic limit.The thermodynamic free energy (43), can be evaluated explicitly: via the Dedekind η function: Notice that the Dedekind function ( 46) is defined only in the upper complex plane Im z > 0, which implies that the free energy ( 45) is well-defined if and only if the causality condition (42) is satisfied.The causality condition is absent for the case of the imaginary rotation (5), when the angular frequency Ω becomes a purely imaginary quantity.Indeed, according to the analytical properties of the Dedekind function (46), the free energy ( 45) is a well-defined analytical function for any real value of the imaginary angular frequency Ω I at any radius of the ring R. Consequently, the rigid imaginary rotation, contrary to rotation in Minkowski spacetime, can be formulated in the thermodynamic limit.
For convenience of our subsequent analysis, we consider all physical quantities in units of the inverse length 1/β of the imaginary time direction.We introduce the dimensionless length of the ring (the one-dimensional volume) L and the frequency of rotation ν, respectively: The normalized frequency ν corresponds to the normalized statistical angle (7).
The free energy density F I = F I /2πR is given by After taking the thermodynamic limit R → ∞ and in the absence of rotation, F I → F 0 = −P 0 = −π/6β 2 , as established by Eq. (17).Note that F 0 agrees with the classical expression F(0) in Eq. (30).The first term in Eq. ( 48) could be erroneously mistaken for the regularized zero-point (Casimir) energy contribution to the free energy.To show that this identification is not correct, let us consider the trivial case ν = 0 which corresponds to a vanishing statistical angle, χ = 0 (nontrivial angles will be considered shortly after).The low-temperature limit, β → ∞, for a ring with a fixed radius R corresponds to a vanishing parameter L. Using the following relation, valid for vanishingly-small positive L, (where the ellipsis denote subleading terms in the limit L → 0), one gets from Eq. ( 48) that the normalized free energy vanishes in the low-temperature limit: For the sake of convenience, we present here the expressions for the normalized Casimir free energy and the Casimir pressure, respectively: which coincides with the first term in Eq. (48).The thermodynamic contribution (48) does not contain the zeropoint energy since the latter should diverge in the L → 0 limit (51), which is not the case (50).Therefore, Eq. ( 48) represents a purely thermodynamic contribution.

D. Fractalization of thermodynamics
Using the property η(−z * ) = [η(z)] * valid for any complex number z from the upper complex semi-plane, Im z > 0, one gets for the thermodynamic part of the free energy density (48), the following expression: where we used the notations in (47).The thermodynamic limit, L → ∞, of the free energy on the ring (52) can be deduced from the beautiful result of Ref. [43], which relates the Dedekind η function (46) with the Thomae f T function (15) as the following limit: Applying (53) to the thermal part of the free energy density (52), we get that the (normalized) thermodynamic energy density "fractalizes" in the thermodynamic limit: lim The non-analyticity of the Thomae function f T has a fractal nature [43] implying the fractalization of thermodynamics under imaginary rotation [32].The result in Eq. (54) implies that close to the thermodynamic limit, L ≫ 1, the non-analytical fractal part of the free energy density dominates over the analytical term in the thermodynamic free energy (52).Interpreting the expression for F in light of Eq. ( 17), we are led to define a rotationdependent temperature T T (ν) via which depends on the statistical parameter ν ≡ βΩ I /2π = p/q via the discontinuous Thomae function f T as given in Eq. (14).Notice that in the thermodynamic limit, the thermodynamics of the system is determined by the ninionic statistics (12).This fact can be seen from the expression for pressure, given in Eq. ( 54), which coincides with Eq. (18).
The fractalization (54), characterized by the nonanalytical behaviour of free energy, is achieved only in the thermodynamic limit when the radius R of the ring becomes infinitely large.At any finite R, all thermodynamic characteristics of the systems are analytical.Thus, it is instructive to see how thermodynamics acquires its fractal properties under imaginary rotation as the radius of the ring increases.48), shown in units of pressure P0 for an infinite ring L → ∞ in the absence of rotation, Eq. ( 17), as a function of the normalized statistical angle ν = χ/(2π) ≡ βΩI /(2π) for various (normalized) lengths L of the ring (47).The plots at finite values of L are given for the analytical behaviour (48) in terms of the Dedekind η function (46).The behaviour in the thermodynamic limit, L → ∞, corresponds to the nonanalytical fractal result (54) expressed via the Thomae function (15).The behaviour near the points ν = 0 and ν = 1 is not shown to preserve a convenient vertical scale.
In Fig. 5 we show the average free energy of bosons (48) as a function of the normalized statistical angle ν, Eq. ( 47), at various (normalized) lengths L of the ring.Pressure, which is a smooth analytical function of ν at small radius of the ring, develops a series of minima and maxima as the length of the ring increases.For large sizes L ∼ 10 3 , pressure of bosonic particles develops selfsimilar features.At L ∼ 5 × 10 3 , the pressure becomes almost indistinguishable from its limiting form (L → ∞) given by Eq. ( 18) and governed by the fractal properties of the Thomae function (15).In this limit, the thermodynamic pressure becomes a fractal dictated by the ninionic statistics (12).
To see how the other thermodynamic quantities fractalize in the imaginary rotation case, it is convenient to rewrite the free energy in Eq. ( 44), as follows.Employing a series expansion of the Bose-Einstein factor, ln we arrive at Switching now to imaginary rotation, Ω = iΩ I , we have The pressure and average angular momentum can be derived from Eq. ( 22): where we introduced M I = −iM = ∂F I /∂Ω I , while S I = β(P I − F I ) + 2πνM I .Using Eq. ( 24), it can be shown that E I = P I .
To reveal the fractal features of the expressions in Eqs. ( 58) and ( 59), we consider ν to be a rational number, represented by the irreducible fraction ν = p/q.Then, the trigonometric function sin 2 (πjν) and sin(πjν) cos(πjν) have a periodicity with respect to j → j + q.To take advantage of this periodicity, we write j = qQ + j ′ , with 0 ≤ Q < ∞ and 1 ≤ j ′ ≤ q, covering the entire 1 ≤ j < ∞ summation range.Out of the new range for j ′ , the value j ′ = q is special, since then the trigonometric functions vanish: In what follows, we split our observables A ∈ {F I , P I , M I , S I , E I } as where A q corresponds to the j ′ = q contribution and δA collects all terms with 1 ≤ j ′ < q.In particular, we have where the right arrow → implies taking the thermodynamic limit L → ∞.The above also imply that The presence of the denominator q of the irreducible fraction p/q = ν in Eqs. ( 62)-( 63) is the hallmark of fractal thermodynamics and hence of loss of analyticity.The remainders δA in Eq. (61) play the role of "defractalizing" the expectation values A for small systems (small L).At large L, it is not difficult to see that δA decays to 0. For example, in the case of the free energy, we have The above sums are readily evaluated when considering the derivative of δF I L with respect to L: Integrating now the series expansion appearing on the second line of the above equation, we find the leadingorder behaviour for L → ∞: Finally, we evaluate the K 2n coefficients defined in Eq. ( 29), which in the present case reduce to where, as before, K 2n ≡ F(0)K 2n .The above sums are performed numerically and the results can be seen in Fig. 6.Strikingly, the shape coefficients K 2n (L) exceed their thermodynamic limit, K 2n (L) > K 2n (∞) = 2n!, for any finite value of L.
For the purpose of evaluating the coefficients F(0, L) and K 2n (L) analytically, we rewrite F from Eq. ( 58) as The coefficients K 2n can be computed by expanding ∆F in powers of ν, where we took into account that F is an even function of ν.Specifically, we have Taking into account that (νL) 2n = (−1) n v 2n R , we can identify We further expand the coefficients ∆F (2n) in powers of L −1 : In the above, we denoted where the derivative is taken with respect to x = 2πj/L, and we took into account that ∆F (2n) ∼ x −2n−1 when x → 0. We arrive at (74) The sum over j can be taken in terms of the polylogarithm, leading to it can be seen that the k = 0 term gives the leading-order, L-independent contribution to F(0)K 2n , while the k = 0 and k = 1 terms contribute to the next-to-leading-order L −1 ln L contribution.The coefficient of the L −1 term receives contribution from all values of k: Alas, the sum over k seems not to converge.This is easily visible in the case of F(0), when where B k represent the Bernoulli numbers, with ] diverges, we refrain from computing the L −1 term and list only the leading and next-to-leading terms, while warning the reader that the L −1 ln L term may also be flawed: where we took into account that ∆F (2,0) = ∆F (2,1) = −2 and ∆F (4,0) = ∆F (4,1) = 24.Thus, we can extract K 2 and K 4 as where we neglected O(L −1 ) contributions.A comparison between the above analytical estimates and the exact numerical evaluation of F(0), K 2 and K 4 is shown in Fig. 6.

E. Analytical continuation: the disk of analyticity
Finally, let us discuss how the fractalization of thermodynamics in the thermodynamic limit leads to the absence of the analytical continuation from the imaginary angular frequencies to the real ones.In other words, we would like to see that the thermodynamic quantities obtained in an infinite volume limit at imaginary rotation cannot be directly connected to the thermodynamics of real rotation.Qualitatively, the validity of this statement can be deduced from both mathematical and physical arguments.
Mathematically, it is clear that a non-analytical function cannot be analytically continued to an analytical domain because the result will depend on the prescription used for the continuation procedure.Moreover, at Analytical FIG. 6.Effect of ring radius R = βL/2π on the average free energy F(0, L) and on the coefficients K2 and K4, computed in the limit of vanishing rotation from Eq. ( 67).The quantity F0 = −P0 is given in Eq. ( 17).The dashed gray curve corresponds to the anayltical estimates in Eqs. ( 80) and ( 81).
the imaginary-rotating side in the thermodynamic limit, the pressure P cannot be expressed as a function of the imaginary velocity squared, Ω 2 I , Eq. ( 18), which renders inapplicable the continuation prescription to the real rotation summarized in Eq. (6).
Physically, the causality condition ( 42) is incompatible with the continuation prescription (6) in the thermodynamic limit, R → ∞, for any finite Ω I .However, outside of the thermodynamic limit at any finite R, the analytical continuation does exist.Let us briefly discuss this point for the example of the ring.
Since the length of the ring is always a positive number, L > 0, the argument of the Dedekind η function in the free energy density (52) always belongs to the upper part of the complex plane, where the Dedekind function is an analytical and well-defined function.Therefore, the imaginary rotation is well-defined at any imaginary angular frequency ν, contrary to its real counterpart (45).
It is convenient to write explicitly the normalized free energies β 2 F for real (F Re ) and imaginary (F Im ) angular frequencies, respectively: Here we defined the following real-valued quantities: which represent the normalized angular frequencies for real and imaginary rotation (x R and x I , respectively), and the inverse size of the ring y > 0. The analytical continuation can be formulated in terms of the relation between ( 82) and ( 83).Despite similarity of Eqs. ( 82) and ( 83), these quantities have different properties.The free energy for real rotation (82) is defined only in the strip −y < x R < y because the Dedekind eta function is defined only in the upper part of the complex plane (excluding the real axis).Physically, the same condition coincides with the causality requirement (42).The pressure of the gas under imaginary rotation ( 83) is defined for any real-valued x I ∈ R.
Equations ( 82) and ( 83) can be written in the following unified form: with For any y > 0, the analyticity properties of the free energy density at real rotation (82) imply that function (85) can be expanded in a series of powers of z around the point z = 0 in the disk |z| < |z 0 | of radius |z 0 | = y > 0. In the original notations, the disk of analyticity can be defined as the generalization of the causality condition (42) in the plane of complex angular frequencies: The free energy F can be written as the following series with the first two coefficients in the explicit form: One can check explicitly that the coefficients of the series (88) diverge in the thermodynamic limit which is consistent with the shrinking radius of convergence (87) as R → ∞.We show the first three non-zero coefficients in Fig. 7.One can also rewrite the free energy (88) in terms of physical variables: where Ω = Ω + iΩ I and and the radius of convergence in the complex Ω plane is determined by Eq. ( 87): Ω c = 1/R.The radius shrinks to zero, Ω c → 0 as R → ∞, thus implying the absence of the direct analytical continuation between real and imaginary angular frequencies in the thermodynamic limit.Thus, thermodynamics of an infinite-volume system subjected to imaginary rotation is not directly connected to the thermodynamics of real rotation.

F. How fractalization emerges as volume grows
Figure 5 shows that at any fixed statistical angle χ = 2πν (or, equivalently, at any imaginary frequency Ω I ) and any finite radius L, the free energy is described by a smooth analytical function of χ.For a rational normalized angle ν = p/q with coprime integer numbers p and q (0 < p < q), the pressure depends both on the numerator p and the denominator q (we remind that in these units, Ω I = (2π/β)p/q).However, as the length L of the ring increases, the pressure turns into a fractal, implying that it loses the sensitivity to the numerator p and keeps only the dependence on the denominator q that defines the imaginary frequency.
At small and moderate ring sizes up to L ≃ 4, the free energy F = F(L, ν) depends on the normalized statistical angle ν monotonically, with −F(L,  17), for a nonrotating ring, in the infinite-volume limit, L → ∞.For rational ν, the free energy (54) in the infinite-volume limit takes fractal values (shown by the arrows) dictated by the Thomae function (15).The inset shows the zoom in on the small-radius region.
for 1/2 > ν a > ν b in the mentioned set of values.In other words, in the analytical region, the thermodynamics of the system behaves analytically, exhibiting a dependence on the actual value of the rational-valued normalized statistical angle and not on its numerator or denominator separately.
As we have seen above, the transition to the fractal regime is associated with the loss of the analytical continuation from imaginary to real rotation.For purely imaginary rotation, the convergence segment (87) for the variable ν is defined by the condition: implying that for the largest value, ν = 1/2, the nonanalytical regime should come into play at the ring length L = 2, while for the smallest nonzero value, ν = 1/10, the critical length is larger, L = 10.This region of the lengths -shown in the inset of Fig. 8 -is characterized by the breaking of the monotonic behavior of pressure on the statistical angle, which is a precursor of the fractal features observed at larger lengths of the ring.At higher values of L, the behavior of pressure on ν becomes more peculiar.To see this in detail, it is convenient to start from the non-rotating case, ν = 0/10 = 0, and associate it to the pair (p, q) = (1, 1) since ν = 0/10 ≡ 0 and ν = 1/1 ≡ 1 correspond to the same static case related to each other by the translation symmetry, ν → ν + 1.The ν = 0 pressure, characterized by the denominator q = 1, is shared both by real, Ω = 0, and imaginary, Ω I = 0, static cases.In Fig. 5, it provides us with a benchmark value for the pressure in the large-L limit.
The values of the statistical angle ν = 1/10 and ν = 3/10 correspond to rotations with different imaginary angular frequencies Ω I = π/(5β) and 3π/(5β), respectively, but they share the same denominator q = 10.According to the fractalization property (14), both these cases -which are characterized by the pairs of coprimes (p, q) = (1, 10) and (p, q) = (3, 10), respectively -should correspond, in the thermodynamic limit, to the pressure of free bosonic gas at the same temperature T = 1/(10β) which is ten times smaller than the temperature in the non-rotating Ω I = 0 (ν = 0) case.The pressure for ν = 1/10 and ν = 3/10 is, consequently, 1/q 2 ≡ 1/100 of the gas pressure in the absence of imaginary rotation.The described features are clearly seen in Fig. 5: the ν = 1/10 and ν = 3/10 pressures, very different at low L ∼ 1, start to approach each other at L ∼ 10, converging into a single curve already at L ∼ 50.This asymptotic behaviour has fractal features as the thermodynamics of the gas is sensitive only to the denominator of the rational (properly normalized) angular frequency.
The cases ν = 2/10 ≡ 1/5 and ν = 4/10 ≡ 2/5 correspond to the coprime pairs (p, q) = (1, 5) and (p, q) = (2, 5) that share the same denominator q = 5.The pressure for these imaginary angular frequencies collapse to a single line even earlier, at L ∼ 7, as it can be seen from the inset of Fig. 8.In both cases, pressure approaches the result for a free bosonic gas with temperature T = 1/(5β) which is q 2 = 25 times smaller than the pressure of the non-rotating gas.
Finally, the normalized statistical angle ν = 5/10 ≡ 1/2 gives the denominator q = 2, temperature T = 1/(2β) and a gas pressure which is q 2 = 4 times smaller than the one of the non-rotating gas.
The monotonic analytical behaviour of the free energy F I (ν) ≡ F I (ν, L), seen at small lengths of the ring L, [small L (analytical)] : (94) (95) ), as it is clearly seen in Fig. 5.

G. Negative thermodynamic pressure of ninions
Apart from the fractal features of the thermodynamic limit -already anticipated from the analytical approach discussed earlier -the pressure at finite volumes L ∼ 1 . . . 10 appears to possess an unexpected feature.Namely, there are regions of the statistical angle χ where the thermal contribution to pressure is negative, as it is clearly seen in Fig. 5.In this sense, the "ninions" -the auxiliary particles which are associated with the ninionic deformation of the standard statistical distribution ( 12) -provide us with the similar phenomenon as the Casimir effect with one important difference: the negative "ninionic" pressure is produced by thermal, and not quantum, fluctuations.As temperature rises, the negative pressure rises as well.
The effect of the negative pressure appears in the analytical region (93) as it is seen in Fig. 8 and especially in the inset of this figure.This unusual behavior is an exotic property of ninions which is not associated with the fractal statistics.

H. Effect of finite particle mass
We now consider the effect of finite mass M on the fractalization properties discussed in Subsec.III D. Due to the relativistic dispersion relation ω m = M 2 + m 2 /R 2 , the resummation of the free energy as in Eq. ( 57) is no longer possible.Writing F I for imaginary rotation parameter Ω = iΩ I directly from Eq. (44) gives Due to its slow convergence, the above sum is not amenable to asymptotic analysis.Instead, we rely on brute force to evaluate Eq. ( 96).It can be seen in Fig. 9 that the plateau (fractal) value is reached for any values of βM .For convenience, the plot reports the ratio F I (L)/F 0 , where F 0 = −P 0 = −π/β 2 represents the free energy of a boson gas [see Eq. ( 17 Finally, we analyze the effect of the particle mass M onto the free energy in the absence of rotation, F(0, M ), and on the coefficients K 2n .The relevant expressions are (97) Figure 10 indicates a generally decreasing trend of the free energy density F(0, M ), computed for a static system (ν = 0) at finite mass, with respect to its massless, thermodynamic limit value, F 0 = −π/6β 2 .The shape coefficients K 2n (M ) increase with M with respect to their classical, massless values, K 2n = (2n)!.The results are consistent with those observed in the classical analysis shown in Fig. 3 and one can expect to reach agreement in the thermodynamic limit, when L → ∞.

IV. RIGIDLY-ROTATING BOSE-EINSTEIN DISTRIBUTION
We now move on and consider a (3+1)d rigidlyrotating system comprised of uncharged, massless boson particles.In this section, we consider such a system  from the perspective of relativistic kinetic theory, which is introduced briefly in Subsect.IV A. In Subsect.IV B, we discuss the thermodynamic properties of a rigidlyrotating system with real rotation parameter Ω, whose properties for slow rotation are discussed in Sec.IV C. Finally, Subsec.IV D is dedicated to the case of imaginary rotation.

A. Relativistic kinetic theory
Although throughout this article we consider noninteracting bosonic systems, it is instructive to discuss, for a brief moment, an interacting model.This approach will allow us to elucidate thermal distributions in thermodynamic equilibrium and shed some light on the physical nature of imaginary rotating systems.
In relativistic kinetic theory, the system dynamics are described using the relativistic Boltzmann equation [46][47][48]: where f k ≡ f k (x) is the one-particle distribution function and k µ = (k 0 , k) is the on-shell momentum satisfying k 2 = 0.The macroscopic properties of the system can be described using the energy-momentum tensor, where dK = gd 3 k/[(2π) 3 k 0 ] is the Lorentz-invariant integration measure and the degeneracy factor of a single neutral scalar field considered in this paper is g = 1.
The conservation law ∂ µ T µν = 0 demands that k µ be a collision invariant, i.e.
The prototypical collision term is that corresponding to 2-to-2 scattering processes [46,49], where fk = 1 + f k is the Bose enhancement factor and the Lorentz-invariant transition rate W kk ′ →pp ′ can be written in terms of the quantum-mechanical differential cross section dσ/dΩ as with s = (k + k ′ ) 2 and Θ s being the center of mass squared energy and emission angle, respectively, with According to the H theorem, the collision term C[f ] drives the system towards local thermal equilibrium, described for the case of a free bosonic gas by the Bose-Einstein distribution: where T (x) is the local temperature and u µ (x) is the local four-velocity.It is easy to check that C[f ] = 0 when the gas is in thermal equilibrium, i.e. f * = f BE * for * ∈ {k, k ′ , p, p ′ }.In global thermal equilibrium, f k = f BE k at each space-time point and Eq.(98) becomes: where β µ (x) = u µ (x)/T (x) is the temperature fourvector.Thus, in global equilibrium, β µ satisfies the Killing equation, ∂ µ β ν + ∂ ν β µ = 0.In this paper, we seek the solution that corresponds to rigid rotation: where Ω is the angular velocity and β = 1/T 0 is a constant corresponding to the inverse temperature on the rotation axis (where x = y = 0).The equilibrium distribution (104) thus reads: where with ρ = x 2 + y 2 being the distance from the point at (x, y, z) to the rotation axis.

B. Thermodynamics of rigid rotation
We now discuss the properties of the rigidly-rotating global equilibrium state.From the relation β µ β µ = 1/T 2 (x), we can identify the local temperature as [47] T where γ(ρ) is the Lorentz factor of a co-rotating particle at a distance ρ from the rotation axis.Relation (109) corresponds to the Tolman-Ehrenfest law [50,51], which relates local temperature to the metric in a static gravitational field, in the curvilinear background of co-rotating reference frame.Similarly, the local four-velocity reads Both the Lorentz factor and the local temperature T (ρ) = γ(ρ)/β diverge on the light cylinder, as ρΩ → 1.
The energy-momentum tensor T µν (x) can be obtained via In the case of massless particles, the energy E = 3P is expressed via the local pressure P , which reads We now consider the thermodynamic limit of our rigidly-rotating system.Identifying F (ρ) = −P as the local free-energy density, we consider the average free energy F = F/V in a cylindrical volume V = πR 2 L z of height L z and radius R, centered on the rotation axis.The mean free-energy density reads The same result can be obtained starting from the expression for the grand potential F of relativistic bosons in rotation [44], cf.Eq. ( 20), where J z = −k φ is the z component of the particle's angular momentum.Other thermodynamic quantities can be obtained starting from Eq. ( 22).Taking into account that V = πR 2 L z , it can be seen that the radial and vertical directions are not equivalent.Therefore, we replace the term PdV by with the hydrostatic pressure obtained as the weighted average P = (2P R + P z )/3.Thus, the thermodynamic pressures are given by Similarly, the average entropy S, angular momentum M, and energy density E are given by Eqs. ( 23) and (24), respectively.We now evaluate the above quantities using the classical expression in Eq. ( 113): It can be checked that E and F satisfy Eq. ( 25).Furthermore, E is compatible with the classical relation [cf.(27)]: with where we took into account the explicit expressions for the local hydrostatic pressure (112) and the Lorentz factor (109).
C. Slow rotation: moment of inertia and shape We now consider the coefficients F(0) and K 2n introduced in Eq. (29).Comparing F in Eq. (113) for the (3 + 1)d system and in Eq. ( 26) for the (1 + 1)d ring, it can be seen that the Ω dependence is through the same γ 2 (R) factor.Therefore, the coefficients K 2n are identical to those in Eq. (30), The average free energy in the absence of rotation evaluates to The moment of inertai in the zero-rotation limit I 0 evaluates to The coefficient I 2 in the series The positiveness of I 2 > 0 implies that the rotating matter tends to increase its angular momentum with an increase in angular frequency.This property signals the change in the shape of rotating system leading to a spatial redistribution of energy as a result of the rotation, which can already be seen from Eq. (119): rotation tends to increase the contributions to the energy density (118) coming from the outer regions as compared to those coming from the inner ones.The physical situation is somewhat similar -neglecting viscosity effects -to water rotating in a glass: its moment of inertia increases with rotation because the distribution of mass within the glass changes, with the water particles moving away from the axis of rotation, increasing the distance of each mass element from the axis, and, hence, the moment of inertia becomes larger.Finite-size corrections, related to the finite transverse size of the system and, consequently, to quantization of the transverse modes, will be discussed below in Subsects.V H and VI C.

D. Imaginary rotation
We now turn to the case of imaginary rotation.Setting Ω = iΩ I with real Ω I is not possible directly in f BE , because that would lead to a complex-valued distribution function.Instead, we can consider the properties of the system described by the distribution which is nothing but a form of the ninionic deformation of the Bose-Einstein statistics (12).Since β µ (iΩ I )∂ µ = β(∂ t + iΩ I ∂ φ ) still satisfies the Killing equation, the left-hand side of the Boltzmann equation (98) vanishes.Somewhat unsurprisingly, the collision term on the righthand side of the same equation does not vanish.This can be seen by considering the small-Ω expansion of f BE k (Ω) introduced in Eq. (107): where Considering now Ω → ±iΩ I and taking the average as described in Eq. (124) gives Taking this substitution back into the collision term (101) shows that It can be seen that in general, C[f ] does not vanish when f k = f I k , hinting that thermal equilibration will generically reduce the magnitude of Ω 2 I .Keeping in mind that imaginary-rotation states are not in actual thermal equilibrium -in a sense that their deformed distribution (124) does not have the equilibrium Bose-Einstein form (2) and that the collision integral (127) does not vanish -we can still derive the macroscopic energy-momentum tensor, which becomes now diagonal: with where γ I (ρ) is reminiscent of the Lorentz factor of corotating particles, Notice that the Euclidean version of the quantum Tolman-Ehrenfest effect gives a different Lorentz factor [26]: where [x] 2π = x + 2πk ∈ (−π, π], with k ∈ Z, is invariant under the 2π symmetry enforced by the natural periodicity of the imaginary rotation (10).The apparent noncompliance of the kinetic Euclidean Lorentz factor (129) with the periodicity requirement (10) can be traced back to the continuous nature of the angular component k φ of the momentum (108).
From a thermodynamic point of view, the structure of the energy-momentum tensor reveals an underlying equilibrium (perfect fluid) contribution, T µν pf;I = diag(E I , P I , ρ −2 P I , P I ), with hydrostatic pressure P I = E I /3, and a shear-stress tensor π µν I = T µν I − T µν pf;I with components It is instructive to note that, on the rotation axis, ) is independent of Ω I , while at ρ = √ 3/|Ω I |, the energy density reaches 0. At larger distances, E I decreases to a minimum (negative) value −π 2 /(9720β 4 ) (reached at ρ = √ 5/|Ω I |) and afterwards increases asymptotically towards its limit 0. In this largeρ limit, Eq. (128) shows that T µν I behaves as follows: with γ I ∼ (ρ|Ω I |) −1 .Thus, far away from the rotation axis, the azimuthal pressure becomes negative and three times larger in magnitude than the energy density, while the radial and vertical pressures remain positive, each being equal in magnitude to the energy density.We now consider the large-volume limit of our system.The average energy inside a cylinder of radius R is simply which agrees with the expression in Eq. (117) under the substitution γ(R) → γ I (R).Substituting P I = E I /3 in Eq. ( 113) will clearly give a different result for the average free energy than F in Eq. (113).To achieve agreement up to the substitution γ(R) → γ I (R), we must replace P (ρ) by P ρ;I (ρ).This choice is supported also by the more fundamental expression for the free energy obtained by setting Ω → iΩ I in the third line of Eq. (114), i.e.
which is consistent with the expression in Eq. ( 25).Applying the thermodynamic relations in Eqs. ( 23) and (116) gives expressions for quantities analogue to the system pressure, entropy, and angular momentum: with M I = −iM = ∂F I /∂Ω I [see Eq. ( 59)].The above quantities are compatible with the Euler relation ( 24), formulated now for a system under imaginary rotation.

E. Effect of finite mass
It is interesting to check how a finite particle mass affects the properties of the rigidly-rotating boson gas.Starting from Eq. (114), one may perform a Lorentz transformation in the d 3 k integral to the local rest frame, such that Writing d 3 x = ρdρdφdz and d 3 k = k 2 dkdΩ k , the φ, z and Ω k integrals can be performed automatically, leading to (138) Expanding the logarithm as in Eq. ( 35), and introducing y = 1/γ(ρ), with dy = −ργ(ρ)Ω 2 dρ, leads to where we changed the momentum integration variable using kdk = k 0 dk 0 .Changing again the integration variable to z = k 0 /M and using ∞ 1 dz e −az = K 1 (a)/a, we arrive at 11. Effect of particle mass M on free energy and shape coefficients K2 and K4, respectively, computed in relativistic kinetic theory for the (3+1)d system (see Fig. 3 for explanation about notation).The inset magnifies the 0 ≤ βM ≤ 0.2 region to demonstrate the domain of applicability of the asymptotic formulas in Eqs. ( 142) and ( 143) .
We are now in a position to evaluate the average free energy in the absence of rotation, F(0, M ), as well as the shape coefficients K 2n (M ) = d 2n F(Ω, M )/dΩ 2n ⌋ Ω=0 : At small M , it is possible to derive the asymptotic expansions where K 2n ≡ F(0, M )K 2n and F 0 = −π 2 /90β 4 is the free energy density for a massless bosonic gas.The coefficients K 2 and K 4 have therefore the following asymptotic behaviour: The above estimates are validated in Fig. 11 by comparison to the exact formulas in Eq. ( 141).The results obtained in this section for the (3+1)d system are in good qualitative agreement with those obtained for the (1+1)d ring (cf.Figs. 3 and 10).

V. IMAGINARY ROTATION OF THE SCALAR FIELD VS. REAL ROTATION A. Mode solutions
We consider a real, massless scalar field φ.The decomposition of the field operator φ(x) reads as follows: where f j (x) represent a complete basis of orthonormal mode solutions of the Klein Gordon equation, These modes are taken as eigenfunctions of the Hamiltonian H = i∂ t , momentum component P z = −i∂ z , and angular momentum component J z = −i∂ φ : with q j = ω 2 j − k 2 j .The one-particle operators âj satisfy the canonical commutation relations, where δ(j, j ′ ) = ω −1 j δ(ω j − ω j ′ )δ(k j − k j ′ )δ mj ,m j ′ .The sum over the quantum numbers is abbreviated as In this and subsequent subsections, we pursue, for simplicity, a "hybrid" quantization approach: we work in the basis of the cylindrical waves (146) with continuous transverse momentum q j which is typical for the unbounded system while restricting the integral (148) over the longitudinal momentum, |k j | ⩽ ω j , to preserve the hermiticity of the Hamiltonian.This set of modes is not suitable to describe rigid real rotation, since in that case, the system must be enclosed inside a boundary in order to preserve causality [33], as will be discussed in Sec.VI.In this section, we will focus primarily on the study of rigid rotation with imaginary angular velocity, for which no causality issues arise and the set of eigenmodes presented above is perfectly applicable.

B. Thermal states
The statistical operator for a thermal state at inverse temperature β which rotates with angular velocity Ω is ρ(β, Ω) = e −β(: H:−Ω: Ĵz :) , where we took the normal-ordered operators Using the commutation relations it is not difficult to establish that where ωj ≡ ωj (Ω) = ω j − Ωm j represents the co-rotating energy (3).The thermal expectation value (t.e.v.) of an arbitrary operator A(x) is where Z = Tr(ρ) is the partition function.Using Eq. ( 152), the t.e.v. of the product of two one-particle operators can be seen to satisfy Using the commutation relations in Eq. ( 151), we establish Introducing the functions the scalar condensate2 becomes Considering now the conformal energy-momentum tensor, defined classically as [52] T its expectation value T µν = ⟨: T µν :⟩ can be expressed as 000 , (159b) where we introduced the notations: while all other components vanish.
Turning back to the definition of the functions G abc in Eq. ( 156), we immediately notice divergences associated with the Bose-Einstein factor [e β ω − 1] −1 .For each value of m such that Ωm > 0, there will be a value of ω where this factor diverges.The only notable exception is the rotation axis, where J 2 m (qρ) = δ m0 and [e βω − 1] −1 has the usual Bose-Einstein infrared divergence when ω = 0. Thus, we are led to conclude that thermal rigidly-rotating states of the scalar field are illdefined at each point outside the rotation axis due to long wavelength (super-horizon) modes, for which ω ≤ Ωm [53].We will come back to this issue in Sec.VI when we will discuss the Klein-Gordon field enclosed inside a cylindrical boundary.

C. Evaluation for imaginary rotation
We now seek to construct states which undergo imaginary rotation, Ω = iΩ I , where Ω I ∈ R. As also noted in Sec.IV D, a state under imaginary rotation leads to complex expectation values of physical observables.This problem can be alleviated by considering the hermiticized version of ρ, namely which is equivalent to averaging over the results obtained for positive and negative values of Ω I .Under the above hermitization, the t.e.v. in Eq. (155) becomes which is similar to the relativistic kinetic theory distribution f im k in Eq. ( 124) under the substitution k φ → −m j .The thermodynamic state corresponding to the t.e.v. ( 162) is characterized by ninionic statistics (12).In what follows, we perform the calculations considering averages using the statistical operator ρ(β, iΩ I ), keeping in mind that the final result is obtained by taking the real part.
In order to analyse the functions G I abc , obtained by replacing Ω → iΩ I in Eq. ( 156), the Bose-Einstein factor can be expanded in a power series, as follows: where ω = ω − iΩ I m has a positive real part, ω > 0. Writing it can be seen that the power of m can be accounted for by taking derivatives with respect to the rotation parameter: On the other hand, the sum over m can be performed in G j;I ab0 using the summation theorem for Bessel functions [see Eq. ( 8.531.1) in Ref. [54]; and Eq. ( 10.23.7) in Ref. [55]]: leading to where x = jβω and θ is defined by (k, q) = ω(sin θ, cos θ), while α j is given by where, for convenience, we reproduced Eq. ( 47) and introduced other notations to be used later (notice that 0 ⩽ l ⩽ L).
In order to perform the integral with respect to θ in Eq. ( 167), we replace the Bessel function J 0 (x) by its series expansion, The integral with respect to θ can be performed now term by term using the relation (valid for Re γ > −1) Using the following identities for the gamma functions, we arrive at corresponding to the cases b = 0 and 2 in Eq. ( 167).The summation can be trivially performed, π/2 finally arriving at or equivalently, with Re(z) = (z + z * )/2 and Im(z) = (z − z * )/2i being the real and imaginary parts of a complex number z, we obtain Specifically, for ϕ 2 and T µν , we require: Furthermore, G j;I 001 , G j;I 101 , and G j;I 002 can be obtained by means of Eq. ( 165): Finally, the functions G (1);j;I 000 and G (2);j;I 000 corresponding to the expressions in Eq. (160) are G (1);j;I 000

G
(2);j;I 000 2α 2 j dG (1);j 000 Substituting Eq. (180c) into Eq.( 157) allows ϕ 2 to be expressed as while the non-vanishing components of T µν , given in Eq. ( 159), reduce to The above results show that the diagonal components of T µν are real-valued and even with respect to ν = βΩ I /2π, while T tφ I is imaginary and odd with respect to ν → −ν.Under the hermitization (161), it is clear that T tφ I vanishes and T µν I remains diagonal.As in the case of the classical relativistic kinetic theory (RKT) analysis in Sec.IV D, the resulting state is not isotropic.Identifying as in the classical case E I = T tt I and the perfect-fluid contribution T µν pf;I = diag(E I , P I , ρ −2 P I , P I ) with P I = E I /3, the quantum shear-stress tensor π µν I = T µν I − T µν pf;I = diag(0, π ρρ I , π φφ I , π zz I ) can be obtained as At large distances from the rotation axis, α j → ∞, implying that The structure of the above result is similar to that obtained in Eq. ( 132), with the important difference that the quantum field-theoretical (QFT) T µν I depends on ν through the harmonic function sin(πjν).This property implies that T µν I obtained in QFT is periodic with respect to ν, with period ∆ν = 1, in agreement with the symmetry (10) expected on general grounds.This periodic behaviour is fundamentally different from that observed in RKT (Subsect.IV D), where no periodicity in ν can be seen.

D. Values on the rotation axis: no analytical connection between real and imaginary rotations
On the rotation axis, we have α j = 0. Using the relations ∞ j=1 j −2 = π 2 /6 and ∞ j=1 j −4 = π 4 /90, we find that the expectation value of ϕ 2 is not affected by the imaginary rotation while the energy-momentum tensor acquires a nontrivial dependence on the imaginary angular frequency: 12. The condensate ϕ 2 and the components of the energy-momentum tensor T µν on the axis of rotation of the cylinder, ρ = 0, under (a) imaginary rotation with ν ≡ νI = βΩI /2π and (b) real rotation with νR = βΩ/2π, normalized with respect to their values in the absence of rotation.All quantities under imaginary rotation are periodic (ν → ν + 1) in agreement with Eq. (10).The dashed lines extending in the region |νR| > 1 indicate the expected behaviour if the components of Tµν were periodic with respect to νR.
with Li s (x) = ∞ k=1 x k /k s being the polylogarithm and B n ≡ B n (0) being the Bernoulli numbers: In terms of the Bernoulli polynomials, the components of the energy-momentum read as well as ρ 2 T φφ I ρ=0 = T ρρ I | ρ=0 .These are the results on the axis of rotation for an infinite-volume system subjected to the imaginary rotation.
We now compare the results in Eq. ( 184) to those derived on the basis of real rotation in Refs.[56][57][58][59] using a perturbative approach for slow rotation, reproduced below for definiteness: and ρ 2 T φφ ρ=0 = T ρρ | ρ=0 , with ν R defined in Eq. ( 169).The above results can be derived from Eq. (184) using the following replacements: It is remarkable to observe that the diagonal components of T µν satisfy however contrary to the same quantities evaluated under imaginary rotations, they do not exhibit periodicity with respect to ν R .Figure 12(b) shows that when |ν R | > 1, T zz increases dramatically, while T tt and T ρρ = ρ 2 T φφ eventually become negative.The dashed lines extending in the region |ν R | > 1 indicate the expected behaviour if T µν were periodic with respect to ν R .Before ending this subsection, we remark that the presence of odd powers of ν in the expressions for T µν I is unexpected and seemingly unsupported by the formulas in Eq. (181).For example, in the case of T tt I given by Eq. (181b), a Taylor expansion of the summand around ν = 0 fails to capture the ν 3 term revealed in Eq. (184b).Moreover, since ν always appears multiplied by the summation variable j, such an approach can reliably produce only the first two terms, proportional to j −4 ν 0 and j −2 ν 2 .The third term proportional to j 0 ν 4 cannot be computed due to the divergence of the sum over j.We are thus led to believe that the ν 3 term appearing in Eq. ( 184) is related to an inherent non-analytic behavior of T µν I with respect to the rotation parameter ν.We remark that the results quoted in Eqs.(189) for the case of real rotation were obtained also using a Taylor series approach and may therefore omit similar non-analytical ν 3 R -like terms.

E. High temperature expansion
Let us now consider the large temperature expansion, when β → 0. Since β comes multiplied by j under the summation sign in Eqs.(181), higher-order terms with respect to β come with higher powers of j.Since the summation over j and the power series with respect to jβ in general do not commute, this procedure allows only the coefficients of the β −4 and β −2 terms to be extracted.The results are As discussed in the previous subsection, the ν 3 terms revealed in Eq. ( 184) are not captured by the perturbative series expansion approach.Nevertheless, the results reported in Eq. ( 192 We now consider the case when ν = p/q is a rational number, where p/q is an irreducible fraction, as discussed in Sec.III D. Writing j = Qq + j ′ , with 0 ≤ Q < ∞ and 1 ≤ j ′ ≤ q, the trigonometric functions taking as argument jπν = πpQ+j ′ πp/q depend only on j ′ .Specifically, Eq. (181) becomes where we introduced the notation while j ′ was relabeled as j for convenience.The sum over Q introduced by the procedure shown in Eq. ( 193) can be performed using where /ϕ 2 0 l = 2πρ/β 0/10 (q = 1) 5/10 (q = 2) 2/10 (q = 5) 4/10 (q = 5) 1/10 (q = 10) 3/10 (q = 10) 10 −5 |T tt |/T tt 0 l = 2πρ/β 0/10 (q = 1) 5/10 (q = 2) 2/10 (q = 5) 4/10 (q = 5) 1/10 (q = 10) 3/10 (q = 10) FIG. 13.Fractalization of thermodynamics with increasing volume: The thermal expectation values of (a) ϕ 2 and (b) T tt under imaginary rotation normalized with respect to their values in the absence of rotation (ϕ 2 0 = 1/12β 2 and T tt 0 = π 2 /30β 4 ) as functions of dimensionless distance l = ρ/(2πβ) from the rotation axis of the cylinder.The lines and points show results for rational values of ν = βΩI /2π of the form r/10, with 0 ≤ r ≤ 5, corresponding to irreducible fractions p/q with q = 1, 2, 5 and 10, which are identical to the imaginary frequencies used for the rotating ring in Fig. 8.The horizontal dashed black lines represent the expected large-distance plateau given by (a) 1/q 2 and (b) 1/q 4 , which signal the fractal behaviour of thermodynamics.The gray dotted lines represent the relativistic kinetic theory prediction in Eq. (128b).A small segment of the result for T tt when p/q = 1/10 corresponds to negative values and is represented with dashed lines.The values are obtained using Eq.(200).ψ(x) = Γ ′ (x)/Γ(x) is the digamma function and the primes denote differentiation with respect to the argument, e.g.ψ ′′ (x) = d 2 ψ(x)/dx 2 .Also, Im and Re denote the real and imaginary parts of their arguments, respectively: Im When considering the summation over j appearing in Eq. (193), a special case corresponds to j = q, when the sine term s j = sin(πjp/q) cancels.In this case, we employ Since s j = 0 implies also x j = 0, Eqs.(193) show that the j = q contribution becomes l-independent.This allows all expectation values to be split as where the first terms are coordinate-independent and correspond to a bosonic gas at rest with inverse temperature qβ: The factor q represents the denominator of the irreducible fraction ν = p/q.More importantly, because these terms are independent of the transverse distance to the rotation axis given by l, they become dominant at large distances from the rotation axis, giving rise to fractal structures in the thermodynamic (infinite volume) limit.It is noteworthy that the fractal terms are completely absent in the relativistic kinetic theory analysis in Sec.IV D and thus represent a purely quantum effect.The second terms in Eq. (198) "defractalize" the result close to the rotation axis and are computed via: where ψ j , x j and s j were introduced in Eqs. ( 194) and (196).
For ν = 1/2, we have: The behaviour of the scalar condensate ϕ 2 I and energy density T tt I as functions of l is illustrated in panels (a) and (b) of Fig. 13, respectively, where we consider the cases ν = p ′ /10 with 0 ≤ p ′ ≤ 5, corresponding to irreducible fractions p/q with q = 1, 2, 5 and 10.As l > 1, visible differences between the curves corresponding to various values of ν can be seen.Contrary to the classical case shown in Eq. ( 132), the far-field behavior of ϕ 2 I and T µν I is dominated by quantum effects.An estimate of how these observables approach their asymptotic values ϕ 2 q and T µν q can be obtained by considering the decay of the "classical" part from Eq. ( 132) to values of the same order of magnitude as ϕ 2 q and T µν q , which occurs at values l ≳ l q , where The q 2 dependence of l q is confirmed for both ϕ 2 I and T tt I , however the p dependence appears to be negligible.The emerging fractal behaviour exhibits a stark contrast to the classical result in Eq. (128b) derived within relativistic kinetic theory, which is also shown in Fig. 13(b) using dashed gray lines.Sizeable deviations can be seen for the curves with smaller value of q, which reach the fractalized plateau at smaller values of L. In the p/q = 1/10 case, the RKT curve follows closely the QFT one, providing a good approximation also in the region where T tt I becomes negative.Noting that the RKT result for p/q = 3/10 falls off too rapidly compared to the QFT curve leads us to conclude that the classical RKT description becomes valid only in the limit ν → 0.
As discussed above, the fractal behaviour manifests itself at large distances from the rotation axis, i.e., as l → ∞. Figure 14 illustrates the expectation values of ϕ 2 I /ϕ 2 0 (left) and T tt I /T tt 0 (right) with respect to ν for various values of l.We considered ν = p ′ /q ′ with 1 ≤ q ′ ≤ 20 and 0 ≤ p ′ ≤ q ′ , covering all irreducible fractions p/q with 1 ≤ q ≤ 20.These results are represented with purple circles.We also considered a set ν j (0 ≤ j ≤ n = 20) of "irrational" values of ν, represented using green squares, obtained as: where −0.01 < δν j < 0.01 is a random number. 3In order to employ a logarithmic scale on the vertical axis, we represented the absolute values of our observables, with the convention that filled and empty symbols are used when the observables are positive and negative, respectively.Since ϕ 2 I > 0 for all values of ν and l, this discussion applies only to T tt I (see, e.g., the (p, q) = (1, 10) curve in Fig. 13).
For l ≲ 1, both ϕ 2 I and T tt I exhibit a smooth dependence on ν.As l is increased, the expectation values for the case when ν = p/q is an irreducible fraction become frozen on their corresponding asymptotic values (1/q 2 for ϕ 2 I /ϕ 2 0 and 1/q 4 for T tt I /T tt 0 ), earlier for smaller values of q than for larger values of q.In contrast, the expectation values corresponding to the irrational values of ν continue their decreasing trend towards 0.
Strikingly, the thermodynamics of the scalar field in the (3+1)d cylindrically symmetric space, obtained numerically and shown in Fig. 14, resembles drastically the one at the ring, obtained analytically and represented in Fig. 5, with the fractalization features becoming more pronounced as the distance l to the rotation axis is increased.

G. Thermodynamic limit
For the purpose of analyzing the large-volume limit of our system, we consider a fictitious cylinder of radius R ≡ βL/2π and of large vertical extent L z , centered on the rotation axis.The volume-averaged scalar condensate and energy density can be computed by integrating Eqs.(198), (200a) and (200b) over this cylinder and dividing by the total volume V = πR 2 L z : where ϕ 2 0 and T tt 0 were introduced in Eq. (199), s j was defined in Eq. (194), while X j corresponds to the old x j evaluated at the volume boundary: 3 For j = 0 and j = 20, we employed ν 0 = |δν 0 | and ν 20 = 1 − |δν 20 |, respectively, in order to ensure that 0 ≤ ν j ≤ 1.
Furthermore, we keep the notation Γ j and ψ j introduced in Eq. ( 196), but now we understand that these functions take the argument j q + iX j .We now compute the average free energy F I from E I starting from Eq. ( 25): where F 0 (β) = −π 2 /90β 4 is the free energy of a bosonic gas in the absence of rotation.
The entropy and angular momentum given by Eq. ( 23) require taking derivatives of F I with respect to β and Ω I at constant Ω I and β, respectively.This is not possible at the level of the fractalized form in Eq. ( 206).Thus, we seek to obtain the free energy after rewriting E I for general (not necessarily rational) values of ν: It can be checked that writing ν = p/q and j = qQ + j ′ gives Eq. (204b).Applying now Eq. ( 25) leads to where the term π/2 appearing on the second line represents an integration constant such that lim Ω I →0 F I = − ∞ j=1 1/(π 2 β 4 j 4 ) = F 0 .Using Eq. ( 23), the average entropy S I and angular momentum M I are Considering as before that ν = p/q and writing j = qQ + j ′ , with 0 ≤ Q < ∞ and 1 ≤ j ′ ≤ q, we get where s j and c j were introduced in Eq. ( 194), ψ j in Eq. ( 196) with x j replaced by X j , while X j is defined in Eq. ( 205).Furthermore, the following notation was introduced: Comparing Eqs.(210a) and ( 206), it can be seen that The above identity is easily established by noting that Integrating the above expression with respect to X j and demanding S Q (X j = 0) = 0 gives Eq. (212).A similar expression can be obtained for S ′ Q .Taking the derivative with respect to X j eliminates the arctangent, such that Integrating the above relation with respect to X j gives For the case when ν = p/q = 1/2, we find while M I = 0.
H. Slow rotation: moment of inertia and shape In the case of slow rotation, we take the free energy in the absence of rotation, F(0), as well as the first two coefficients, K 2 and K 4 , by performing a series expansion with respect to ν in Eq. ( 208): where we took into account that Ω 2 I → −Ω 2 for the purpose of extracting K 2n .We also employed the notation K 2n ≡ F(0)K 2n .
The leading term F(0) and the first coefficient K 2 , which corresponds to the dimensionless moment of inertia, can be computed exactly: In the case of the rotational shape coefficient K 4 , the summation of the coefficient of L −4 cannot be performed.However, its leading-order behavior and the first Ldependent correction can be extracted as follows: Comparing Eqs. ( 218)-( 219) to (30) [see also Eq. ( 120)], we see that the classical result K 2n = (2n)!receives Ldependent corrections that vanish as the transverse size of the cylinder becomes infinite, L → ∞.As in the case of the (1 + 1)d ring, shown in Fig. 6 [see also Eq. ( 81)], both K 2 (L) and K 4 (L) exceed at finite L their thermodynamic limits, K 2n (∞) = (2n)!.In Subsect.VI C, we discuss the dimensionless moment of inertia K 2 and the rotational shape change coefficient K 4 in a cylinder of a finite radius, taking into account the proper quantization of the radial modes.We will see that the behaviour of K 2 qualitatively agrees with Eq. (218) while the coefficient K 4 becomes finite and converges to the classical result (120) in the infinite volume limit.

I. Effect of finite mass
In this subsection, we discuss the effect of a finite mass on the fractalization of thermodynamics (seen at the level of the local observables), as well as on the slow rotation coefficients, K 2n .For this purpose, we start the analysis with the component T tt I , which is still given by Eq. (159a).The only difference is that the ω integration in the functions G I abc introduced in Eq. ( 156) now starts from ω = M , Focusing on the functions G I ab0 = ∞ j=1 (−1) j G j;I ab0 , as pointed out in Eq. ( 164), and taking the same steps as described in Subsec.V C, we arrive at Eq. ( 167), which is modified to where we introduced y = [x 2 −(jβM ) 2 ] 1/2 corresponding to the momentum magnitude, p = √ ω 2 − M 2 .In the case b = 0, which is relevant to the computation of T tt I , the θ integral can be performed using Eq.(174a): The integral in (222) proves difficult for general values of the parameters.We now focus on T tt I , requiring G I 200 and the function G (2);I 000 = ρ∂ ρ ρ∂ ρ G I 000 , for which we find: This allows T tt j;I to be obtained using Eq.(159a) as x 2 sin(α j y) 0/10 (q = 1) 5/10 (q = 2) 2/10 (q = 5) 4/10 (q = 5) 1/10 (q = 10) 3/10 (q = 10) FIG. 15.Ratio T tt I )/T tt 0 between the energy density at finite mass, given in Eq. ( 224), and T tt 0 = π 2 /30β 4 .The rotation parameter takes the values ν = p ′ /10, with 0 ≤ p ′ ≤ 5, leading to irreducible fractions p/q with q ∈ {2, 5, 10}.The dotted colored lines and points correspond to regions where T tt I (M ) < 0, while the dashed black lines represent the thermodynamic limit derived in Eq. ( 225).
To study the fractal properties of T tt I , we consider a rational rotation parameter, ν = p/q, with p < q being irreducible coprime numbers, and write j = qQ + j ′ , with 1 ≤ j ′ ≤ q and 0 ≤ Q < ∞.Since α j = l πj sin(πjν) vanishes when he j ′ = q, the corresponding contribution T tt q becomes coordinate-independent and evaluates to The results of the numerical evaluation of T tt I (M ) for ν = p/q with q = 10 and 0 ≤ p ≤ 5 are shown in Fig. 15.The fractal pattern can easily be seen as L is increased.The large-L value coincides with the expression in Eq. ( 225), shown with horizontal dotted black lines.
We now attempt to evaluate the average free energy F(0, M ) in the absence of rotation, as well as the shape coefficients K 2 and K 4 .For this purpose, we evaluate the average energy E j;I = 2 R 2 R 0 dρ ρT tt j;I starting from Eq. ( 224): where A j = (L/πj) sin(πjν).Computing the free energy via F I = β −1 dβ E I is cumbersome due to the the βM dependence.Since we are interested only in the slow rotation limit, we can perform an expansion with respect to v R = iνL, where The integral with respect to x can be performed exactly for all terms displayed above.However, the x −1 factors appearing in the L −2 and L −4 terms (which vanish in the thermodynamic limit) lead to results in terms of the exponential integral function, Ei(−jβM ), which makes analytical manipulations difficult.Instead, we focus on the leading-order contributions, for which we get where the arrow → indicates that the thermodynamic limit L → ∞ was taken.The summation over j can be performed in Eq. ( 229) in terms of the polylogarithm functions.However, at this stage we can already obtain the coefficients of the free energy density, where K 2n ≡ F(0, M )K 2n , by using the relation Performing now the summation over j via F I = ∞ j=1 F j;I , we arrive at Evaluating the above for small values of the mass, we get the following asymptotic behaviour:

B. Scalar condensate, energy-momentum expectation values and fractalization
Figure 16 shows the main features of ϕ 2 I (top panel) and T tt I (lower panel) as functions of l = 2πρ/β for several different radii R chosen such that the quantity takes the values L = 10, 100, and 1000.Only the case of rational (imaginary) rotation parameter is considered, with ν = βΩ I /2π = p ′ /10 and 0 ≤ p ′ ≤ 5, giving rise to all irreducible fractions p/q ≤ 1/2 with q = 1, 2, 5, and 10.The dashed black lines represent the results obtained in the unbounded case, computed based on Eqs.(198), (200a), and (200b).As expected, the Dirichlet boundary conditions considered in this section affect the behaviour of the observables close to the boundary.Specifically, ϕ 2 I = 0 when ρ = R, since G I 000 (ρ = R) vanishes identically by virtue of the quantization condition J m (q mn R) = 0; while T tt I is decreased by about a factor of 10 compared to its bulk value.In both panels, the bounded and unbounded results stay in good agreement throughout most of the cylinder if L is sufficiently large.In particular, ϕ 2 I exhibits a notably smaller value on the rotation axis when L = 10 compared to the unbounded case, while for the L = 100 and 1000 cases, good agreement can be seen.
As mentioned in Sec.V B, the boundary permits the study of a system undergoing real rigid rotation, as long as ΩR = νL ≤ 1 and the light cylinder is excluded from the system.It is thus interesting to compare expectation values computed for imaginary and real rotation, ν I and ν R .To keep the comparison meaningful, both ν I and ν R are restricted to be lower than or equal to 1/L.Fig. 17 shows the radial profile T tt (ρ) for three cylinders, with L = 1 (a), 10 (b) and 100 (c), in the case of slow (νL = 0.1, blue), medium (νL = 0.5, red) and fast (νL = 1, green) rotation.At small L = 2πR/β, the boundary effects dominate over thermal ones and T tt decreases monotonically from the rotation axis towards the boundary.Furthermore, T tt (ρ = 0) is strongly suppressed (by four orders of magnitude at L = 1) compared to its value for a boson gas at rest, T tt 0 = π 2 /30β 4 .The effect of imaginary rotation is negligible, while in the case of real rotation, T tt (ρ) increases slightly at νL = 1.At L ≳ 10, the bulk of the system is dominated by thermal effects.Panels (b) and (c) also show the RKT result for T tt : where γ 2 = 1/(1 − ν 2 R l 2 ) and γ 2 I = 1/(1 + ν 2 I l 2 ).In the case when L = 100, the QFT results deviate from the RKT ones only in a small vicinity of the boundary.Such good agreement is also a consequence of the fact that at large L, ν is constrained to be small.As discussed in Sec.V F, RKT is expected to agree with QFT for small values of ν and sufficiently far from the boundary (see also Fig. 13. Next, the value of T tt (ρ = 0; ν) on the rotation axis for both real and imaginary rotation is shown in Fig. 18 for (a) L = 1, (b) L = 10 and (c) L = 100.As before, in the case L = 1, the value of T tt is suppressed by over three orders of magnitude.Here, the rotation parameter covers an entire period of the system undergoing imaginary rotation.In the case of real rotation, no such periodicity arises, contrary to the expectation based on the result in Eq. (189b) for the unbounded case.The inset in panel (a) shows the effect of rotation on T tt (ρ = 0; ν) for smaller values of ν.It can be seen that for ν ≲ 0.1, the quantity |T tt (ν) − T tt (0)| has the same behavior for both real and imaginary rotation.At L = 10 and 100, T tt (ρ = 0; ν = 0) approaches T tt 0 = π 2 /30β 4 .The maximum value of ν is greatly reduced.It can be seen that for such a small interval of ν, the result corresponding to imaginary rotation is approximately equal to that obtained for real rotation, mirrored with respect to the value T tt (ρ = 0; ν = 0) obtained in the absence of rotation, as shown with black dotted lines.Further-more, the red dotted lines indicate the analytical predictions in Eqs.(184b) and (189b), scaled by the value T tt (ρ = 0, ν = 0)/T tt 0 on the rotation axis in the absence of rotation, corresponding to the given value of L. The agreement with the analytical predictions is better at larger values of L, which may also be due to the smaller range allowed for ν.
We now consider the thermodynamic system as a whole and discuss the volume-averaged free energy F = F/V , computed via the equivalent of Eq. (114): Applying Eqs. ( 116) and ( 23) leads to the following expressions for the radial pressure P R , average entropy S, and average angular momentum M: while P z = −F.The average energy E and scalar condensate Φ 2 can be obtained by taking the volume average of the expressions in Eq. ( 236): In deriving the above expressions, we employed the integration formula together with the Dirichlet boundary conditions J m (q mn R) = 0. Comparing Eqs.(240e), (240a) and (240b), it is easy to see that E = P/3 with P = 2 3 P R + 1 3 P z being the isotropic pressure.Moreover, the relations in Eq. ( 240) are compatible with the Euler relation (24).Φ 2 I /ϕ 2 0 L = 2πR/β 0/10 (q = 1) 5/10 (q = 2) 2/10 (q = 5) 4/10 (q = 5) 1/10 (q = 10) 3/10 (q = 10) 0/10 (q = 1) 5/10 (q = 2) 2/10 (q = 5) 4/10 (q = 5) 1/10 (q = 10) 3/10 (q = 10) FIG.19.Integrated quantities Φ 2 I 2 0 and EI /T tt 0 for the case when the system is enclosed within a cylindrical boundary located at R = βL/2π, with L shown on the horizontal axis.The colored lines with points show the results for rational (imaginary) rotation parameter ν = p ′ /q ′ shown in the legend (q ′ = 10 and 0 ≤ p ′ ≤ 5), corresponding to the irreducible fraction p/q, with q shown between the parentheses.The black dotted lines represent the results obtained in the unbounded case, shown in Fig. 13.ary effects become localized around a small vicinity of the boundary and the bounded and unbounded results approach each other.While for Φ 2 I , visible discrepancies remain even for L ≳ 10 3 , in the case of E I , F I , S I , and M I , the results obtained in the bounded case start following the ones corresponding to the unbounded case already when L ≳ 10.As in the previous sections, the fractal structure reveals itself at large values of L.
C. Slow rotation: moment of inertia and shape Finally, we discuss the expansion in Eq. ( 29) of the free energy density F in the case of slow rotation.In particular, we focus on the free energy in the absence of rotation, F(0), as well as the first two coefficients, K 2 (related to the moment of inertia) and K 4 (related to rotation-induced deformability), which we evaluate us- 0/10 (q = 1) 5/10 (q = 2) 2/10 (q = 5) 4/10 (q = 5) 1/10 (q = 10) 3/10 (q = 10) 10 −5 ing: where K 2n = K 2n F(0).The coefficients F(0), K 2 and K 4 are studied as functions of the transverse size of the system L with respect to their unbounded counterparts.The ratios F(0)/F 0 , K 2 /2! and K 4 /4! are represented in Fig. 21, where the denominators of these expressions are the classical expectations given in Eq. (120).As we already noted in Subsect.IV C, a strong quenching due to the boundary can be seen at small values of L, which is however less pronounced for the shape coefficient K 4 compared to the moment of inertia coefficient K 2 and the average free energy F 0 .At L = 100, these coefficients already reach their expected asymptotic values K 2n = (2n)!, Eq. (120).Surprisingly, both K 2 (L) and K 4 (L) evaluated at finite L are smaller than their asymptotic values, K 2n (∞) = (2n)!.This property is in contrast to the behaviour seen for the (1 + 1)d ring [see Fig.In a cylinder of a finite radius R, a large-size L → ∞ limit corresponds also to the high-temperature limit, T → ∞ since L = 2πR/β ≡ 2πRT .Figure 21 shows that the dimensionless moment of inertia K 2 approaches the asymptotic value K 2 = 2 from below, indicating that the moment of inertia should decrease as temperature decreases.This effect is related to the presence of an effec-tive energy gap between the states with zero, m = 0, and non-zero, m ̸ = 0, orbital momenta due to the finite size of the system.Therefore, at lower temperatures, the system mostly resides in the m = 0 state and the rotational modes, which contribute to the moment of inertia (244b), are not excited.Since the latter modes do not participate in rotation at low T , the moment of inertia of the system decreases as the system gets colder.This effect should evidently also occur for K 4 and higher coefficients.Interestingly, the same qualitative behaviour for the moment of inertia K 2 is also observed in the first-principle simulations of gluon plasma in the high-temperature phase of Yang-Mills theory: as the temperature increases, the moment of inertia approaches the high-temperature value K 2 = 2, Eq.(120), from below [27].
The above discussion indicates that the restriction in size of the quantum system reduces the K 2n coefficients.We now consider the behavior of the K 2n coefficients for particles of finite mass M .As indicated in Sections III A, III H, IV E, and V I, the K 2n coefficients increase as M increases.In Fig. 22, we demonstrate that the same qualitative behavior is achieved for the spatially bounded system.

VII. CONCLUSIONS
In the present work, we studied the thermodynamic properties of a massless scalar field subjected to rigid rotation and an inter-relation of the real rotation with its imaginary analogue.The latter concept -a rigid rotation with an imaginary frequency [19,23] -has a practical interest since rotating systems cannot be implemented in Euclidean path-integral formalism, suitable, for example, for numerical first-principle calculations on the lattice [17,21,26,27,30].In this sense, rotation shares the deficiency suffered by finite-density systems, namely the sign problem [30], and needs to be implemented in Euclidean spacetime via its imaginary version supplemented with the subsequent analytical continuation to real angular frequencies [17,21,27].
Using the 1 + 1-dimensional toy model of a scalar field under rigid rotation on a ring, we explicitly demonstrated that the analytical no-go theorem [32], which describes the impossibility of continuation of the imaginaryangular-frequency thermodynamics to the real angular frequencies is related to the development of the fractality of thermodynamics for the former.The result applies in the thermodynamic limit.Within this model, thermodynamic functions such as pressure and energy density can be expressed analytically via the Dedekind η function (46).The latter function tends to the fractal, non-analytical Thomae function (15) as its argument approaches the real axis [43], which corresponds to the infinite-volume (thermodynamic) limit.
In the case of the (3+1)-dimensional Minkowski space, we first considered a classical description of scalar particles under rotation using relativistic kinetic theory.In the absence of boundaries, rigid rotation with a real rotation parameter leads to a violation of causality and subsequent divergence of all observables on the light cylinder.Imaginary rotation can be described by a real distribution function only as an average of clockwise and counterclockwise rotations, leading to a seemingly nonequilibrium state.As expected, observables such as the energy-momentum tensor T µν decrease as the distance to the rotation axis is increased, at a faster rate for faster rotation.
Under the quantum field theoretical treatment, rigid rotation with a real rotation parameter leads to the divergence of T µν at each point of the space-time, also inside the light cylinder.Under imaginary rotation, T µν evaluated on the rotation axis can be expressed via Bernoulli polynomials (181).Away from the rotation axis, we were able to demonstrate the fractalization of both the field fluctuations ϕ 2 and of T µν in the case of imaginary rotation, in complete analogy to the 1 + 1-dimensional toy model discussed above.In all cases, our results exhibit a periodicity with respect to the imaginary rotation parameter which is not present in the classical, kinetic analysis.For this reason, we found agreement with kinetic theory only in the limit of slow rotation and only in the vicinity of the rotation axis (before fractalization sets in).
For comparison with the case of real rotation, we took results obtained using a perturbative calculation for slow rotation (with respect to a stationary background) and found on the rotation axis an analytical result (189) which can be related to the one obtained for imaginary rotation in an essentially non-analytical way (190).We conclude that even on the axis of rotation -which appears to be static in the rotating system -the nonanalyticity is strong and unavoidable.
We demonstrated the same analytic-non-analytical transition using numerical calculations for the thermal scalar field system enclosed in a cylinder, undergoing rotation with imaginary angular frequencies.As the radius of the cylinder grows, the pressure becomes a nonanalytical function of temperature expressed, again, via the Thomae function (demonstrated in Figs.[13][14][15][16][17][18][19][20].In this limit, the boundary effects become less important and the results obtained in the unbounded case provide a good approximation for our observables inside the cylindrical boundary.For values of the rotation parameter respecting the causality constraints (i.e., when the light cylinder is outside the boundary), the fractalization features do not appear in the case of imaginary rotation.For sufficiently slow rotation and high temperature, our numerically-obtained results are compatible with the flip Ω 2 I → −Ω 2 from imaginary to real rotation, signaling the restoration of analytical continuation in this limit.
The exotic fractalization properties discussed above are related to a ninionic deformation (12) of statistical distributions at imaginary angular frequencies [32].For this reason, we conclude that the results obtained in the infinite-volume system subjected to imaginary rotation cannot be analytically related to the properties of the physically rotating system (with a real angular frequency).However, we explicitly demonstrated both for the analytically-treatable case on the ring and numerically-accessible case of rotating cylinder that the imaginary rotation in a spatially bounded system in Euclidean space can be continued to the real-frequency domain in Minkowski spacetime provided that the causality is respected for the latter spacetime.
We have also shown that for real-frequency rotation, the dimensionless moment of inertia K 2 , normalized per one degree of freedom, is equal to two, K 2 = 2, in the thermodynamic limit of large radius R of the cylinder (32).The quantity K 2 determines the correction to the free energy (29), due to small nonzero angular frequency, Ω → 0. This result matches well the first-principle result of Ref. [27] on the behavior of gluons in high-temperature limit of Yang-Mills theory.We have also shown that in the thermodynamic limit, the generic expansion of the pressure (free energy) in higher orders in angular frequency ( 29) is characterized by exactly-calculable coefficients K 2n with n = 1, 2, . . ., given by identical expressions (30) for the ring and (120) for the cylinder, which also includes the shape coefficient K 4 = 4! responsible for an Ω 4 correction to the pressure.Below we summarize our main findings: 1. We demonstrated the no-go theorem [32] regarding the impossibility of continuation of the imaginaryangular-frequency thermodynamics to the real angular frequencies using an analytical 1 + 1dimensional toy model, revealing the development of the fractality of thermodynamics for the former.
2. Since fractalization does not show up in the classical, kinetic theory treatment of imaginary rotation, we conclude that this is a purely quantum effect.
3. We found similar fractalization on the unbounded Minkowski space for imaginary rotation, as well as in the thermodynamic (infinite volume) limit of the bounded system.
4. For the case of a causal boundary that excludes the light cylinder, we were able to restore the analytical continuation from imaginary to real rotation in the limit of slow rotation and large temperatures.
5. We attributed the exotic fractalization properties to a ninionic deformation (12) of statistical distributions at imaginary angular frequencies [32].
6.These statements are applicable both for massless and massive fields, with the fractal features become pronounced in the thermodynamic limit for massive particles.While an increasing mass of bosons enhances the moment of inertia K 2 and the shape coefficient K 4 of the rotating gas, this effect vanishes in the high-temperature limit.
The results obtained in this paper shed light on the implications of the effects of imaginary rotation obtained in the context of first-principle lattice simulations and constitute the basis for the analysis of more complicated systems, e.g., free fermions or the chiral phase transition in the effective QCD models such as the Nambu-Jona-Lasinio model or nonlinear sigma models.Extending the present analysis to the case of Dirac fermions is a logical avenue of future research [60].

FIG. 4 .
FIG. 4. Illustration of a particle on a ring of the radius R and the angular coordinate φ.The ring rotates with the angular velocity Ω counterclockwise.

FIG. 5 .
FIG.5.Fractalization of thermodynamics of scalar particles on the ring under the imaginary rotation ΩI : free energy density −FI Eq.(48), shown in units of pressure P0 for an infinite ring L → ∞ in the absence of rotation, Eq. (17), as a function of the normalized statistical angle ν = χ/(2π) ≡ βΩI /(2π) for various (normalized) lengths L of the ring(47).The plots at finite values of L are given for the analytical behaviour(48) in terms of the Dedekind η function(46).The behaviour in the thermodynamic limit, L → ∞, corresponds to the nonanalytical fractal result (54) expressed via the Thomae function(15).The behaviour near the points ν = 0 and ν = 1 is not shown to preserve a convenient vertical scale.

FIG. 7 .
FIG. 7. The first three nonzero coefficients in the series (88) of the free energy of the ring as a function of the (inverse) normalized radius 1/L.The direction of the thermodynamic limit is shown by the arrow.

FIG. 8 .
Figure5shows that at any fixed statistical angle χ = 2πν (or, equivalently, at any imaginary frequency Ω I ) and any finite radius L, the free energy is described by a smooth analytical function of χ.For a rational normalized angle ν = p/q with coprime integer numbers p and q (0 < p < q), the pressure depends both on the numerator p and the denominator q (we remind that in these units, Ω I = (2π/β)p/q).However, as the length L of the ring increases, the pressure turns into a fractal, implying that it loses the sensitivity to the numerator p and keeps only the dependence on the denominator q that defines the imaginary frequency.This curious fractalization transition is shown in Fig.8for the particular set of imaginary frequencies Ω I ≡ 2πν/β = pπ/(5β) with p = 0, 1, . . ., 9. Given the periodicity (10) under ν → ν + 1, as well the reflection symmetry(11) with respect to ν → −ν, applied to the pressure, this particular choice leaves us with six distinct values of the normalized statistical angle: ν = 0/10, 1/10, . . ., 5/10.At small and moderate ring sizes up to L ≃ 4, the free energy F = F(L, ν) depends on the normalized statistical angle ν monotonically, with −F(L, ν a ) < −F(L, ν b ) 10 )<F I ( 3 10 )<F I ( 4 10 )<F I ( 5 10 ), is completely lost for large L giving us the fractal nonanalytical hierarchy: [large L (fractal)] :

FIG. 10 .
FIG.10.Mass dependence of the free energy density F(0, M ), normalized with respect to the classical value F0 = −π/6β 2 , as well as of the shape coefficients K2 and K4.The system size is set to (a) L = 10 and (b) L = 100.

Figure 12 (
a) confirms that T µν I is periodic with respect to the imaginary rotation parameter ν, as implied by the presence of {ν}.The energy density T tt I and the radial pressure T ρρ I are decreased by the imaginary rotation, while the azimuthal and vertical pressures T φφ I = T zz I are increased.An alternative way of characterizing T µν I on the rotation axis is by using the Bernoulli polynomial sin(α j y) − α j y cos(α j y)] .(224)