Building up the elliptic flow: analytical insights

In this paper, we present a fully analytical description of the early-stage formation of elliptic flow in relativistic viscous hydrodynamics. We first construct an elliptic deformation of Gubser flow which is a boost invariant solution of the Navier-Stokes equation with a nontrivial transverse profile. We then analytically calculate the momentum anisotropy of the flow as a function of time and discuss the connection with the empirical formula by Bhalerao {\it et al.} regarding the viscosity dependence of elliptic flow.


I. INTRODUCTION
One of the most intriguing results of ultrarelativistic heavy ion collisions at RHIC and the LHC is the strong collectivity of the created hot and dense matter, especially the considerable elliptic flow [1][2][3][4][5][6]. In non-central collisions, the overlapping region of the colliding nuclei approximately has the shape of an ellipse in the transverse plane. This region expands hydrodynamically, and the initial anisotropy in the pressure gradient gets converted into momentum space anisotropy, resulting in the modulation of the azimuthal angle distribution of charged particles in the final state [7] dN dφ ∝ 1 + 2v 2 cos 2φ .
The coefficient v 2 is called the elliptic flow parameter and is one of the central objects of experimental and theoretical study in heavy-ion physics because it is a sensitive probe of the equilibrium/nonequilibrium properties of the created matter.
While the mechanism to generate v 2 is well understood, little is known about its analytical details. Theoretically, the extraction of v 2 relies heavily on numerical (viscous) hydrodynamic simulations supplemented with some initial condition and the equation of state (see, e.g., [8][9][10][11][12][13][14][15][16][17]). In this paper, we provide a fully analytical description of the early-stage formation of v 2 by deriving and utilizing an approximate elliptic solution of the relativistic Navier-Stokes equation. Such an analysis has long been infeasible due to the difficulty of constructing realistic elliptically-shaped solutions of hydrodynamic equations for which only a few attempts have been made in the literature [18][19][20][21]. The new solution we are able to present here has been achieved along the line of the recent progress in constructing exact solutions of viscous hydrodynamic equations using conformal symmetry [22][23][24][25][26]. Within the region of validity of the solution, one can study explicitly how the elliptic flow develops as a function of proper time and how the shear viscosity affects this evolution.
In Section II, we review the Gubser flow [22,23] which is an exact boost-invariant solution of the relativistic Navier-Stokes equation with nontrivial radial flow velocities. In Section III, we use the so-called Zhukovsky transform to elliptically deform the Gubser flow in the transverse plane and construct an approximate solution. We then calculate in Section IV the momentum space anisotropy of this flow and discuss the connection to the empirical formula proposed by Bhalerao et al. [27] regarding the shear viscosity dependence of the elliptic flow. In the end, we summarise in Section V.

II. GUBSER FLOW
In this section, we briefly review the exact boost-invariant solution of the relativistic Navier-Stokes equation found by Gubser [22,23]. The solution is naturally explained by rewriting the Minkowski metric as whereτ ≡ t 2 −ẑ 2 is the proper time andζ ≡ tanh −1ẑ /t is the space-time rapidity. In this coordinate system, the four-velocityû µ of the fluid (normalized asû µû µ = −1) readŝ withû φ =û ζ = 0 and L is roughly the initial transverse size of the fluid. Theζ-independence of Eq. (3) andû ζ = 0 mean that the flow expands in the longitudinal (ẑ) direction in a boost-invariant way. It also expands in the transverse direction with the transverse velocityv Plugging this velocity profile into the relativistic Navier-Stokes equation and assuming conformal symmetry, Gubser obtained the following exact solution for the energy densitŷ where C > 0 is a constant and the shear viscosityη has been made dimensionless by factoring out the corresponding power of the energy density η 0 ≡η/Ê 3/4 . In Eq. (5), we have defined The following components of the shear tensor will be needed for a later calculation.
As already pointed out by Gubser, the solution has a pathological behavior at large negative values of ρ corresponding to largex ⊥ and/or smallτ . Indeed, when | sinh ρ| ≫ 1 one can approximate so thatÊ The quantity inside the square brackets is proportional to the temperatureT ∝Ê 1/4 and this becomes negative for sufficiently negative values of ρ. This is actually not surprising since the relativistic (first-order) Navier-Stokes equation is known to have unphysical features. 1 As demonstrated in Ref. [24], the problem can be cured by switching to the (second-order) Israel-Stewart equation. For the present purpose, we are not concerned about this issue since one can consider η 0 to be arbitrarily small (or C arbitrarily large) so that the temperature remains positive in a parametrically large region ofx ⊥ andτ .τ In this region, the solution is well-behaved and offers an attractive model for the studies of strongly interacting matter created in heavy-ion collisions as discussed in Refs. [22,23].

III. ELLIPTIC SOLUTION
The Gubser solution described above is cylindrically symmetric around theẑ-axis. Here we relax this restriction and construct an approximate solution which has the shape of an ellipse in the transverse plane. This can be achieved by employing the so-called Zhukovsky (Joukowski) transform 2 which maps a circle onto an ellipse as followŝ where a is a constant and we only consider the region x ⊥ > a. As is manifest in its complex representationω = x + iŷ = ω + a 2 ω , Eq. (11) is a conformal transformation in two-dimensions and therefore the metric is preserved up to a Weyl factor Embedding this in four-dimensions, we obtain where we have relabeledτ = Aτ andζ = ζ. If the last two terms were absent, the metric inside the square brackets would be exactly Minkowskian. Let us find the conditions under which these terms can indeed be neglected. More explicitly, we find Since the elliptic deformation ∼ cos 2φ is an O(a 2 /x 2 ⊥ ) effect, we assume a 2 /x 2 ⊥ ≪ 1 and neglect terms of order O(a 4 /x 4 ⊥ ). Then the (dA) 2 term in (13) can be safely dropped. In order to drop the cross term O(dτ dA) as well, we must additionally assume that x ⊥ ≫ τ and neglect terms of order O(τ a 2 /x 3 ⊥ ) relative to the leading term. Under these conditions, the coordinate systems (τ ,x,ŷ,ζ) and (τ, x, y, ζ) are conformally related, and one can map the Gubser solution expressed in the former coordinates into the latter. 3 The transformation rule is given by To the order of interest, we can approximatê We thus find the energy density Since δE < 0, the curves of constant energy are elliptic, with the major axis in the y-direction. The flow velocity in turn is given by and the shear tensor is Note that δu ⊥ is negative, meaning that both E and u ⊥ are stretched in the y-direction. This may seem contradictory to the standard picture that the stronger pressure (or energy density) gradient in the x-direction develops a stronger flow in the same direction. In fact, there is no contradiction. The negativeness of δu ⊥ (at large x ⊥ ) is a direct consequence of conformal symmetry 4 which dictates a power-law decay of E and the hydrodynamic equation which, in the present accuracy, boils down to (For simplicity here we consider the ideal hydrodynamic equation.) On the other hand, we find where we have neglected u 2 φ ∼ O(a 4 ). Eq. (23) shows that, due to the φ-component of the flow velocity which tends to point to the x-direction (δu φ < 0), the flow is on average stronger in the x-direction, as expected.
Eqs. (19)-(21) define an approximate elliptic solution of the Navier-Stokes equation with the shear viscosity η = η 0 E 3/4 . They satisfy the equation with the accuracy of order O(a 2 /x 2 ⊥ ) relative to the leading (isotropic) solution, 4 The situation is different in the confining (non-conformal) case. Suppose that the energy density decays as a Gaussian instead of a power-law where σy > σx. The Euler equation is then for the φ-dependent part. Thus δu ⊥ is positive in this case.
but break down at order O(a 2 τ 2 /x 4 ⊥ ) which can be checked explicitly. 5 One way to understand this is to notice that there is in fact an O(τ 2 /x 2 ⊥ ) uncertainty in the definition of a 2 . Indeed, we could have employed a nonconstant a 2 (τ ) already in Eq. (11) as long as its τ -dependence is sufficiently weak If Eq. (24) is satisfied, the extra terms that would appear in the transformation as in Eq. (13) are of the same order as, or less than those already neglected. The maximal uncertainty incurred by Eq. (24) is and this affects the Navier-Stokes equation at order O(a 2 τ 2 /x 4 ⊥ ). A similar caveat applies to the eccentricity of the flow which we define as For the ideal solution at early times, we find On general grounds, one expects that ǫ(τ ) decreases with time. Unfortunately, this τ -dependence is not reliably calculable in the present approach due to the above uncertainty in a 2 . We however note that Eq. (25) suggests that ǫ(τ ) decreases quadratically in τ where c > 0 remains undetermined. Eq. (28) appears to be consistent with the result of numerical simulations [9,14] (see also, [28]). Before proceeding to the calculation of elliptic flow with viscous corrections, we make one more simplification. In order to clearly see the role of the viscosity, we shall neglect the terms of order O(τ 2 /x 2 ⊥ ) and O(τ 2 /L 2 ) in the energy density Eq. (19) while keeping powers of η 0 (L/τ ) 2/3 up to quadratic order. Actually, once we go beyond linear order in η 0 , we must be careful with the O(1) terms inside the brackets of Eq. (9). The condition that these terms are parametrically smaller than the O(η 2 0 (L/τ ) 4/3 ) term is where another constraint in Eq. (10) has been included. Thus, strictly speaking, the following analysis is valid in the window as shown in Eq. (29) in the viscous case with η 0 = 0.

IV. VISCOUS EFFECT ON ELLIPTIC FLOW
We now turn to the discussion of the elliptic flow parameter v 2 which characterizes the momentum space anisotropy in the final state. Besides the leading contribution from the ideal flow, we also include the shear viscous correction which plays an important role in the phenomenological study of quark gluon plasma [11,[13][14][15]17].
With an analytic solution of the Navier-Stokes equation at hand, we can get information about v 2 via a closely related quantity [8,9] ǫ p (τ ) ≡ dxdy(T xx − T yy ) dxdy(T xx + T yy ) . Explicitly, we find where π µν is the shear-stress tensor. Unlike v 2 , ǫ p (τ ) in Eq. (30) is defined at all times, and its value at the 'build-up time' where L represents the characteristic transverse size and c s is the speed of sound, is considered to be a measure of v 2 ∝ ǫ p (τ f ). (c s = 1/ √ 3 in a conformal theory.) However, as already mentioned, τ f ∼ O(L) is outside the region of validity for our approximate solution. Nevertheless, we show that Eq. (30) evaluated at τ ≪ τ f provides us with useful analytical insights into the property of ǫ p (τ ) as we extrapolate τ → τ f .
First, let us consider the denominator in Eq. (31). Since ǫ p ∝ a 2 and the numerator is O(a 2 ), we can set a = 0 in the denominator from the beginning. We then recall that in the Navier-Stokes approximation, where σ µν is given by Eq. (21) and in the present approximation. Moreover, u 2 ⊥ ∼ O(τ 2 ) in the factor 1+2u 2 ⊥ can be neglected. It is then straightforward to show that Note that the viscous correction shows up only at O(η 2 0 ). 6 We now turn to the numerator which requires some extra care. After the φ integral, we get where the potentially dangerous O(1/x 2 ⊥ ) terms in Eqs. (20) and (21) which could cause trouble in the remaining dx 2 ⊥ integral have canceled out. In Eq. (37), we recognize two types of contributions with clear but distinct physical 6 In fact, the O(η 3 0 ) term of the denominator (but not the numerator) is divergent because the x ⊥ -integral does not converge as x ⊥ → ∞. This is an artifact of the solution Eq. (9) which becomes unphysical at large x ⊥ , and limits our calculation of ǫp to O(η 2 0 ).
interpretations. The terms proportional to δE and δE 3/4 are due to the spatial anisotropy of the source (energy density). Since the bulk of matter is initially stretched in the y-direction (δE < 0), these terms give a negative contribution to ǫ p . On the other hand, the terms proportional to E 0 and E 3/4 0 are due to the anisotropy of the flow velocity δu. As already discussed in Eq. (23), they give a positive contribution to ǫ p .
We thus evaluate the two contributions separately and find, after dividing by Eq. (36), and Summing the two contributions, we arrive at Surprisingly, ǫ p is negative as a result of the slightly larger contribution from the source anisotropy. This is at odds with the observed behavior in hydrodynamic simulations [8,12,14] which is qualitatively quite consistent with the contribution from the flow anisotropy (39) alone, namely, ǫ p (τ ) is positive, grows quadratically in time ǫ p (τ ) ∼ τ 2 , and the viscosity tends to suppress it. 7 While this discrepancy may seem worrisome, one should notice the large cancelation which resulted in a barely negative value found in (40). This suggests that whether the source anisotropy contribution is large enough to flip the sign of ǫ p is subtle and model-dependent. It is then interesting that our conformal solution reveals and exemplifies the logical possibility that even if can become negative. When this occurs, the simple proportionality between v 2 and ǫ p is far from obvious and may be subject to large 'non-flow' effects (see e.g., Refs. [16,29]).
What are, then, the main characteristics of a given model which determine the sign of ǫ p ? In the model at hand, the transverse flow u 2 x − u 2 y generated by the pressure gradient is weak and insufficient to make ǫ p positive because δu ⊥ is negative (cf., (23)). This is an artifact of conformal symmetry which dictates that the energy density should decay as a power-law. In confining theories including QCD, δu ⊥ will be positive (see Footnote 4), and one therefore expects that the flow anisotropy contribution (39) dominates over the source anisotropy contribution (38). After all, elliptic flow is the anisotropy in the flow velocity, and this is faithfully incorporated in (39), but not in (38). 8 Let us therefore take a closer look at Eq. (39). Parametrically, and to linear order in η 0 , this can be rewritten as where we have used E ∼ T 4 ∼ C 4 /(τ 4/3 L 8/3 ) and the kinetic theory relation η/E 1/4 ∼ η/T ∼ 1/σ (with σ being the cross section of microscopic degrees of freedom) together with an estimate C 3 ∼ dN/dY of the rapidity distribution of particle multiplicity in the final state (see, Eq. (44) of Ref. [22]). Now let us compare Eq. (41) with the following empirical formula proposed in heavy-ion collisions [27] (and also in high-multiplicity proton-proton collisions [33]) where K 0 > 0 is a number of order unity and S ⊥ is the initial transverse overlap area of the two colliding nuclei. Despite its simplicity, the formula in Eq. (42) is quite successful in fitting the RHIC/LHC data and the results of , as well as their uncertainty have been canceled. It is also straightforward to see that the genuine O(a 2 τ 2 /L 4 ) corrections of E and u µ can only enter at even higher order in the results in Eqs. (38)-(40). 8 Another difference from hydrodynamic simulations is the treatment of the initial velocity. In typical numerical simulations, the transverse velocity at the initial time τ = τ 0 is set to be zero, whereas in our solution the velocity (20) is nonzero for any τ 0 > 0. Simulations with nonvanishing u ⊥ (τ 0 ) do exist (see, e.g., [30][31][32]), but so far they only studied v 2 in the final state. It would be interesting to see the effect of the initial velocity on the early-stage development of ǫp in these simulations.
viscous hydrodynamic simulations [12,27,34]. With the natural identification S ⊥ = L 2 , parametrically the agreement between Eqs. (41) and (42) is perfect. 9 Thus we have presented an explicit justification of the empirical formula Eq. (42) to linear order in the 'Knudsen' number K ∼ η 0 (L/τ ) 2/3 . At higher order in K, our result suggests a power series with alternating sign (remember that E ∼ (1 − K) 4 before the d 2 x ⊥ integration). Since this is not a geometric series, we expect a deviation from the formula Eq. (42) starting from O(K 2 ). Nevertheless, Eq. (42) does capture the main feature of the series and this is probably the reason of its successful applications in phenomenology.

V. CONCLUSION
In conclusion, we have constructed an approximate boost-invariant elliptic solution of the Navier-Stokes equation with a well-defined region of validity. In this region, the leading contribution to the momentum anisotropy (30) of the fluid is analytically calculated. Driven by the spatial anisotropy, a stronger flow develops in the direction of larger pressure gradient. While this scenario is well-known and well-tested numerically, in the presence of viscosity it has not been demonstrated with the level of analytical precision we have been able to present in this paper. We also pointed out the potentially large negative contribution to ǫ p due to the source anisotropy which, depending on models, can even flip the sign of ǫ p . Finally, by focusing on the flow anisotropy contribution in Eq. (39), we have discussed the connection with the empirical formula (42) previously suggested in heavy-ion phenomenology.