Anisotropic escape mechanism and elliptic flow of bottomonia

We study the role of anisotropic escape in generating the elliptic flow of bottomonia produced in ultrarelativistic heavy-ion collisions. We implement temperature-dependent decay widths for the various bottomonium states, to calculate their survival probability when traversing through the anisotropic hot medium formed in non-central collisions. We employ the recently developed 3+1d quasiparticle anisotropic hydrodynamic simulation to model the space-time evolution of the quark-gluon plasma. We provide a quantitative prediction for transverse momentum dependence of bottomonium elliptic flow and nuclear modification factor for Pb$\,+\,$Pb collisions in $\sqrt{s_{\rm NN}}=2.76$ TeV at the Large Hadron Collider.


I. INTRODUCTION
Heavy quarks, such as charm (c) and bottom (b) quarks, and their quarkonium bound states (cc and bb) are very useful internal probes of the hot and dense medium created in collisions of heavy nuclei at high energies [1]. Heavy quarkonia created in high-energy collisions at the BNL Relativistic Heavy Ion Collider (RHIC) and the CERN Large Hadron Collider (LHC) have been found to be appreciably affected by the medium. This leads to distinctive features in their observed final yields, which encode information about the thermodynamic and transport properties of the medium. Therefore understanding the dynamics of heavy quarks and quarkonia in a deconfined medium is of great interest for the heavy-ion physics community [1][2][3][4][5][6][7][8].
In the classical picture [2,3], heavy quarkonia embedded in a static, equilibrated quark-gluon plasma (QGP) may survive at temperatures above the QGP crossover temperature due to their large binding energies. However, if the QGP energy density becomes sufficiently high, the resulting Debye screening of the quark-antiquark potential eventually leads to the dissociation of charmonia and bottomonia [2,3]. In this classical picture, the bound states with the largest binding energy, respectively the J/ψ and the Υ, have the highest dissociation temperatures.
In recent years, this simple picture has been challenged by a number of findings. Based on first-principles finitetemperature quantum chromodynamics (QCD) calculations, it was shown that the in-medium quark-antiquark potential contains an imaginary part which is associated with the in-medium quarkonium breakup rate [9,10]. This results in significant thermal widths for quarkonia, at variance with the older assumption of uniquely defined binding energies, and leads to the suppression of quarkonia at temperatures at which they would survive in the traditional scenario [11][12][13][14][15].
Another extension to the standard idea is the introduction of dynamics. In the dynamical picture, quarkonium states are both dissociated and (re-)associated with the rate for each process depending on the open heavy flavor density and temperature of the system [17,18]. In particular, in an evolving QGP with decreasing temperature, a quark and an antiquark pair may combine into a stable bound state which was until then unstable [19,20]. For the fireball created in heavy-ion collisions at RHIC or LHC, recombination seems to be marginal for bottomonia [21,22] while playing a more important role in the observed yields of charmonia [23]. In this paper we will investigate bottomonia, whose dynamical evolution is not significantly affected by regeneration, due to the fact that bb pairs are less abundantly produced by initial hard scatterings than cc pairs. An additional important aspect of the dynamics is that the dissociation process is not instantaneous, but depends on how long the quark-antiquark pairs experience a high medium energy density. This is irrelevant in a static infinite QGP, but becomes important in an expanding finite-sized fireball. In that case, and for pairs in motion relative to the QGP, this translates into a dependence on the in-medium path length of the pairs. Thus, the bound states created in a high energy heavy-ion collision that survive the ensuing immersion in the fireball are those that quickly reach a region of low enough temperature. This path-length dependence of the bound-state survival probability is naturally described within an "escape mechanism" scenario [24][25][26][27]. Since a generic heavyion collision gives rises to a spatially anisotropic overlap in the plane transverse to the beam axis, the path-length dependence of states propagating through the generated QGP results in an anisotropic emission pattern of the quarkonia, measured in momentum space, as was first predicted for J/ψ [28]. For bottomonia, the anisotropic escape probability should constitute the major ingredient to the observed momentum anisotropy, quantified as usual in terms of Fourier harmonics v n .
In this paper, we provide a quantitative prediction for bottomonium elliptic flow v 2 generated by the anisotropic escape mechanism; for a qualitative discussion, see Ref. [29]. In Sec. II we give the details of the model used to simulate the hydrodynamic evolution of the fireball created in Pb + Pb collisions at √ s NN = 2.76 TeV at the LHC. For the bottomonium states, we implement temperature-dependent decay widths and calculate their resulting survival probability in the hot and dense anisotropic medium in Sec. III. Including the feed down contribution to the bottomonium ground state from higher excited states, we find that the elliptic flow of the Υ(1S) is at the level of a few percent and investigate its dependence on the parameters governing the medium expansion in Sec. IV. We also study the dependence of the bottomonium elliptic flow on transverse momentum and centrality.

II. HYDRODYNAMICAL MODEL
For the spatiotemporal evolution of the quarkgluon plasma we make use of the recently developed 3+1d quasiparticle anisotropic hydrodynamics (QaHydro) framework [30]. This framework has been successfully used to describe identified-particle spectra, charged-particle multiplicity vs pseudorapdity, identifiedparticle mean transverse momentum, identified-particle and charged-particle elliptic flow, and Hanbury-Twiss-Brown radii in both LHC 2.76 TeV [31,32] and RHIC 200 GeV heavy-ion collisions [33]. The framework includes both shear and bulk viscosities in addition to an infinite number of transport coefficient.
In the version of the anisotropic hydrodynamics code used, the shear viscosity to entropy density is assumed to be constant and all other transport coefficients are self-consistently determined within a quasiparticle model with the temperature-dependent quasiparticle mass extracted from lattice QCD results for the entropy density [34]. For details of the framework we refer the reader to Ref. [31][32][33] and Ref. [35] which presents a comprehensive review of the approach. More importantly this approach allows one to extend the dissipative hydrodynamical evolution to early times due to an infinite-order resummation in the inverse Reynolds number [36].
We begin the dissipative hydrodynamical evolution at τ = 0.2 fm/c. We use optical Glauber model to construct the initial energy density profile in the transverse plane and take the initial central temperature to be T 0 = 600 MeV and the shear viscosity to entropy density ratio to be η/s = 0.2. These values were determined in prior studies by making fits to ALICE identified-particle spectra obtained in 2.76 TeV Pb-Pb collisions [31,32]. The effective temperature profiles in different centrality classes were obtained from the anisotropic hydrodynamics code and were exported to disk. These exported profiles were then used to construct 4-D interpolating functions which provided the spatiotemporal evolution of the effective temperature in each centrality class. For more details concerning the hydrodynamical framework, fits, and comparisons to 2.76 TeV experimental data, we refer the reader to Ref. [32].

III. DECAY WIDTH AND SURVIVAL PROBABILITY OF BOTTOMONIA
Let us now turn to the description of the production of bottomonia and their behavior in the QGP. The bottomonium states are produced in initial hard scattering processes during the very earliest stages of the heavyion collision. The spatial distribution of the production points in the transverse plane is assumed to follow that of the number of binary collisions, N coll (x, y). We assume a power-law transverse momentum (p T ) distribution of the Υ's obtained from PYTHIA simulations for p + p collisions, scaled by the mass number of the colliding nuclei. Note that this kind of scaling implicitly assumes that the bottomonia do not "flow" with the medium and any v 2 that we obtain in our model will be purely due to the anisotropic escape mechanism. The initial Υ distribution for p + p collision is assumed to be given by [37] with Y being the longitudinal rapidity in momentum space. Here is the most forward rapidity of the bottomonia. Eventually, the momentum rapidity density follows a Gaussian distribution: For our calculations, we integrate over Y and consider the resulting p T distribution. The formation of each bound bottomonium state requires a finite formation time τ form . The value τ 0 form of the latter in the bottomonium rest frame is assumed to be proportional to the inverse of the vacuum binding energy for each state. For the Υ(1S), Υ(2S), Υ(3S), χ b (1P ) and χ b (2P ) states we use τ 0 form = 0.2, 0.4, 0.6, 0.4, 0.6 fm/c, respectively [14]. In the laboratory frame, relative to which a bottomonium state with mass M has transverse momentum p T , the formation time becomes After they are formed, since the bound bb states are color-neutral, their elastic scatterings on the color charges in the QGP are fewer and, because of their large rest mass, they propagate quasi freely, following nearly straight-line trajectories. We use the simulation results from QaHydro framework [30] to obtain the temperature of the medium along the bottomonium trajectory. With this temperature, we compute the thermal decay widths of the bottomonium states, adopting the recent state-ofthe-art estimations of in-medium dissociation of different bound bb states [11], which are here considered for the case of a locally momentum-isotropic medium for simplicity. For a given bound state, the 3-D Schrödinger equation is solved numerically with a temperature-dependent complex heavy-quark potential [16]. The in-medium breakup rates for each state are then computed from the imaginary part of the binding energy as a function of the temperature T /T c . The real part of the binding energy, on the other hand tells us when the state becomes completely unbound. We then set the temperature scale to T c = 160 MeV. Below T c , we assume that plasma screening effects are rapidly diminished due to transition to the hadronic phase and the widths of the states becoming approximately equal to their vacuum widths which in comparison to the in-medium widths are essentially zero.
The potential used for the solution of the 3-D Schrödinger equation is based on a generalized Karsch-Mehr-Satz potential [3] obtained from the internal energy. It includes an imaginary part which emerges due to Landau damping of the exchanged gluons in the hardthermal-loop framework [9][10][11] but does not include effects of singlet to octet transitions [38] and is explicitly based on the high-temperature limit of quantum chromodynamics. By making this assumption we do not include any explicitly non-perturbative contributions to the imaginary part of the in-medium quarkonium potential. In general, the imaginary part of the potential encodes the break-up rate of heavy-quarkonium bound states and can be properly understood in the context of open quantum systems in which the heavy-quark system is quantum mechanically coupled to a thermal heat bath allowing for states to transition from bound states to the environmental sector of the in-medium density matrix [38][39][40][41][42][43][44][45][46][47].
For a bottomonium with transverse momentum p T along the azimuthal φ p direction, which is created with transverse co-ordinates (x, y) at the formation time τ form , the position at any future time τ is given by where v T = p T /E T is the bottomonium transverse velocity and τ = τ − τ form . The thermal decay width Γ(T (x , y , τ )) is then modeled as Γ(T (x , y , τ )) = 2 for > 0, 10 GeV for ≤ 0, where ≡ Im[E bind (T (x , y , τ ))] and ≡ Re[E bind (T (x , y , τ ))] denote the imaginary and real parts of the in medium binding energy, respectively, at a transverse position (x , y ) at time τ . The value of 10 GeV in the second case is chosen to rapidly dissociate states which are fully unbound and, in practice, the results are not significantly sensitive on this value as long as it is large enough to quickly melt the state under consideration. The in-medium decay width for a given state, so obtained, thus determines its survival probability as it propagates in the medium. The final transmittance for a bb bound state labelled by j is given by [11][12][13][14] T j (x, y, p T , where Θ is the usual step function. The final time τ f in the above equation is self consistently determined in the simulation as the proper time when the local effective temperature of the medium becomes less than the freeze out temperature T f = 130 MeV. From the above equation, we can obtain the transmitted spectra as a function of transverse momentum and azimuthal direction of all the produced bottomonium states, dN j p T dp T dφ p = dx dy n coll (x, y)

IV. RESULTS AND DISCUSSIONS
In this section, we numerically evaluate Eq. (6) to calculate the elliptic flow of bottomonia due to escape probability through a medium having an anisotropic shape. We consider Pb + Pb collisions at √ s NN = 2.76 TeV in 40−50% centrality class which corresponds to an average impact parameter of 10.41 fm [48]. In order to compare to the experimental results, we integrate over the entire temperature profile in the transverse plane, to obtain the weighted average "raw" spectra for each bottomonium state, as shown in Eq. (6). To account for the post-QGP feed down of the excited states, we use a p T -averaged feed down fraction obtained from a recent compilation of p + p data at LHC. The inclusive spectra for Υ(1S) is then calculated from a linear superposition of the raw spectra for each state: (GeV/c) .1% and f 3P →1S = 1.5% as adopted from [14]. Since contribution from 3S and 3P states are small, we include their percentage contributions in 2S and 2P states, respectively. The inclusive spectra so constructed is then used to calculate the elliptic flow v 2 of Υ(1S). Note that, while considering feed down, the transverse momentum of the mother and daughter bottomonium states are assumed to be the same. This assumption can be justified by considering the average p T value of the mother excited states and the daughter 1S states. Due to the large mass of the bottomonium states, we find that the mean p T value of the mother excited states and the daughter 1S states are almost identical. Our main objective is to find to v 2 (p T ) of Υ(1S) states generated from anisotropic dissociation in the plasma. However no such experimental measurement is available for the same, with which our model results can be contrasted. Hence to check the viability of our model, in addition to v 2 (p T ) we also calculate the p T dependence of nuclear modification factor, R AA (p T ) that has been widely studied at LHC.
Using Eq. (6) and following standard formulae, we thus calculate R j AA (p T ) and v j 2 (p T ) of j-th bound state. For inclusive Υ(1S) states we make use of Eq. (7). While calculating v 2 , we have not accounted for the reaction plane angle as we have considered smooth optical Glauber model initial conditions.
The p T dependence of v 2 for different directly produced bottomonium states for Pb + Pb collisions at √ s NN = 2.76 TeV in 40 − 50% centrality class is shown in Fig. 1. For each state, v 2 initially increases with p T and gradually tends to saturate. A reverse hydrodynamic mass ordering between the different bottomonium states is observed with more massive state having larger v 2 . At a given p T , owing to weaker binding and hence broader decay width, v 2 is larger for the excited states, as was found in Ref. [21]. However this ordering according to the binding energies gets broken for the χ b (2P ) state, which is seen to acquire smaller v 2 than Υ(2S) state. Even though χ b (2P ) remains practically unbound during the entire evolution of the plasma, due to large intrinsic formation time, the suppression effects remain operative for a shorter period, resulting in lesser flow. The nonmonotonic nature of v 2 (p T ) for the excited states can be attributed to the melting of those states beyond which the decay width is abruptly set to 10 GeV. Note that such non-monotonocity is not seen for Υ(1S) states, for which the melting temperature is around 900 MeV. 1 In Fig. 1, we also show the inclusive v 2 of the 1S state which takes into account feed-down. As the largest contribution to inclusive v 2 of 1S states comes from the decay of 1P states, it is close to the direct v 2 of 1S states. At very high p T 10 GeV, v 2 tends to decrease with increasing p T , an effect which is more prominent for excited states. This reduction can be attributed to the competition between the dynamics of the plasma and that of the bottomonium state under consideration. Due to both its velocity as well as its dilated formation time in the plasma frame, a bottomonium with large p T escapes faster from the plasma, spending less time inside and therefore the suppression effect is less important, leading to a reduced v 2 .
Before moving forward, it might be interesting to look at the p T dependence of R AA within the same model framework. In Fig. 2, we show R AA (p T ) for different bottomonium states with same set of input parameters. For any given state, as expected R AA is seen to increase with p T due to a faster escape from the dense medium. The suppression pattern for different states exhibit a monotonic trend. Owing to weakest binding the 2P states are seen to have largest suppression. For inclusive 1S states the results are in line with the previous calculations [15].
It might be interesting to check the sensitivity of the results to different model inputs. Instead of showing different bottomonium states separately, we only report the inclusive 1S states for this purpose. To check the effect of intrinsic resonance formation time on resulting flow, we calculate v 2 (p T ) assuming similar intrinsic formation time, τ 0 f orm = 0.2 fm/c, for all states. The result is contrasted with our previous calculation (different intrinsic formation times) in Fig. 3. If all the bottomonium states are assumed to be formed at same time, they experience screening effects for a longer time, resulting in larger v 2 , an effect more prominent at higher p T .
Next we check the effect of decay widths of the un-  bound states. In Fig. 4, we show the inclusive v 2 of 1S states for two arbitrary values of Γ ∞ 1 GeV and 10 GeV. Since 70% of the measured Υ(1S) states are direcly produced which are always bound inside fireball, the resulting v 2 values are practically same.
Finally, in Fig. 5, we show the transverse momentum dependence of v 2 of Υ(1S) including feed down contributions from higher excited states for different centrality for Pb + Pb collisions at √ s NN = 2.76 TeV. We see that there is a non-monotonic behaviour of v 2 with respect to centrality and maximum v 2 is obtained for 50 − 60% centrality. While eccentricity increases with centrality, the breakup rates decreases with temperature and goes to zero if the bottomonium has escaped. As we increase the impact parameter, the central temperature decreases and the resulting bottomonium momentum anisotropies decrease because of the decrese in thermal decay width Γ(T ). However spatial anisotropy increases at a much faster rate than the decreasing of thermal decay width leading to the monotonic behaviour of bottomonium v 2 with respect to centrality as seen in Fig. 5.
At this point, it is worth mentioning that for bottomonia that are in motion relative to the expanding quarkgluon plasma, the in-medium dissociation depends on the effective local temperature that the quarkonia experiences due to the relativistic Doppler effect [49,50]. As a consequence, bottomonia with a higher velocity relative to the QGP are expected to be more affected. For effective hydrodynamic expansion, that we consider in the present case, the QGP is initially at rest so the high-p T bottomonia would indeed be more suppressed. However, these high-p T bottomonia have a larger formation time and hence the bottomonium experiences smaller temperatures due to rapid cooling of the medium. Since the dissociation width is smaller at lower temperatures, the suppression is lesser leading to smaller v 2 . Moreover, the averaging over redshifted and blueshifted regions can potentially wash out the influence of Doppler corrected temperature and therefore we expect this effect on v 2 to be small. We postpone this analysis for future work.
Before closing we note that the ALICE Collboration has very recently reported the first measurement of v 2 of inclusive Υ(1S) states in Pb + Pb collisions at √ s NN = 5.02 TeV, at forward rapidity [51]. Due to paucity of statistics v 2 (p T ) has been calculated for a large centrality interval 5 − 60%. The measured v 2 values is consistent with zero and with small positive values predicted by the transport models within large uncertainty. Of course we can not make a one-to-one correspondence with our present calculations, due to difference in beam energy and rapidity interval.

V. SUMMARY AND CONCLUSION
In this paper, we have provided a quantitative prediction for the elliptic flow of bottomonia produced in mid-central collisions in √ s NN = 2.76 TeV Pb + Pb collisions at LHC via an anisotropic escape mechanism. We employed the Glauber model to generate initial distribution of energy density in the plane transverse to the beam axis. Using temperature-dependent decay widths for bottomonium states, we calculated their survival probability when traversing through the hot and dense anisotropic medium formed in non-central collisions. For the expansion of the fireball, we have used results from the recently developed 3+1d quasiparticle anisotropic hydrodynamic simulation. We also accounted for the feed down contribution to the bottomonium ground state from higher excited states. We found that the transverse momentum dependence of the elliptic flow of bottomonia is of the level of few percent, consistent with the finding of Ref. [21]. We also found a monotonic behaviour of v 2 with respect to centrality. For completeness, we have also calculated the p T dependence of R AA which is consistent with the existing theoretical and experimental estimates. Looking forward, it will be interesting to consider the effect of medium-induced transitions between bound states which is predicted from the open quantum system approach [38,46,47,[52][53][54]. In this "state reshuf-fling" scenario, transitions between various bound states become possible which counteracts the usual suppression picture by allowing for the re-formation of otherwise suppressed states even above the hadronization temperature. Since the excited states of bottomonium acquire more elliptic flow due to anisotropic escape mechanism, one may expect to generate larger flow owing to the feed-down from excited states. We leave these questions for future work.