{'fraction_non_alphanumeric': 0.04043309084757232, 'fraction_numerical': 0.013421304911746462, 'mean_word_length': 4.212815990593769, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Elliptic flow from final state interactions in the distorted-wave emission-function model

The effect of final state interactions on the elliptic flow coefficient measured in relativistic heavy ion collisions is investigated within the DWEF formalism established by Miller and Cramer [Phys. Rev. Lett. 94, 102302 (2005), J. Phys. G34, 703 (2007), Phys. Rev. C78, 054905 (2008)]. It is found that the optical potential previously found to give the best fit of particle multiplicity and Hanbury Brown and Twiss radii in RHIC events has a moderate effect on the resulting elliptic flow coefficient v_2. This indicates that final state interactions should be taken into account to confidently predict v_2 to better than ~20% accuracy.


I. INTRODUCTION
Hydrodynamic models have been quite successful at describing single particle hadron spectra measured in experiments at the Relativistic Heavy Ion Collider (RHIC) [1,2,3,4,5].
Recently, progress has been made in going beyond simple ideal hydrodynamic models and using viscous hydrodynamics to determine the value of various transport coefficients in the RHIC "fireball," e.g. the shear viscosity [6,7,8,9]. Of particular interest is the azimuthal anisotropy of produced particles. The existence of strong elliptic flow at RHIC is one of the best indications of thermalization at RHIC [10,11,12,13], and its precise value is sensitively connected to the viscosity of the medium [6].
Not every aspect of these hydrodynamic models, however, is completely physically justified. Due to uncertainty in especially early and late time dynamics in the evolution of a heavy ion collision, the initial and final conditions for hydrodynamic evolution are treated somewhat simplistically. In particular, the so-called Cooper-Frye freezeout algorithm [14] is typically implemented to describe how the medium transitions from hydrodynamic behavior to the essentially free particles that enter the detectors. Using this procedure, the medium, once it reaches some defined freezeout temperature, instantaneously "freezes out" from a hydrodynamic fluid to completely free particles that stream to the detectors. Knowing that this is not entirely physical, it makes sense to investigate how much this treatment might affect the conclusions reached from using these hydrodynamic models.
To that end, Miller, Cramer et al. investigated the effect of introducing some final state interactions to this freezeout behavior in the form of a one-body optical potential that the emitted particles interact with as they "fight" to escape the medium [15,16,17]. The main motivation behind this previous work was the notorious inability of hydrodynamic models to fit two-particle correlation data (so-called Hanbury Brown and Twiss radii) while simultaneously fitting single-particle data [18]. (For a review see Ref. [19].) Thus only multiplicity and HBT radii data were calculated, while the present work extends that to now investigate the effect of these final state interactions on measured elliptic flow. For simplicity only pions-the dominant hadron produced in a RHIC event-are considered.
The distorted-wave emission-function (DWEF) formalism is briefly reviewed in Sec. II, with the calculation of v 2 in this formalism described in Sec. III. Sections IV and V contain the results and conclusions. Details of the calculation are included in Appendix A for those interested, along with a semi-analytic derivation of a simple test case used to test the numerics, in Appendix B.

II. DWEF FORMALISM
The DWEF formalism was established in Ref. [15] and described extensively in Ref. [16], and the relevant parts are briefly summarized as follows.
The main quantity that we are interested in is the detected particle momentum spectra from which v 2 is defined as the second moment in the azimuthal momentum angle with p the momentum of the detected particle and Y the particle rapidity. φ p is the angle of the particle momentum with respect to the collision plane.
S(p, x) is known as the emission function. In conventional hydrodynamical models, it resides in space on a freezeout hypersurface defined by a surface of constant temperature (or other thermodynamic quantity) in the hydrodynamic simulation, with a momentum distribution at each point on the surface given by the appropriate equilibrium (or off-equilibrium in the viscous case) distribution at that given temperature. This freezeout hypersurface represents the surface of last scattering, from which free particles are emitted and travel directly to the detector.
Here, instead of running a full hydrodynamic simulation, we follow Ref. [20] and use an analytical parametrization of the freezeout surface, similar to one typically found in numerical hydro simulations, but with tunable parameters. In addition it is allowed to be a more general volume, with finite width in all dimensions, rather than the infinitely thin surface obtained in a conventional Cooper-Frye prescription.
Secondly, instead of freely streaming, particles that are emitted from this surface are then made to interact with an optical potential representing interactions with the medium from which the particles are escaping.
p (x) are the aforementioned distorted waves (as opposed to plane waves appropriate in the absence of interactions). They obey the equation of motion for pions interacting with optical potential Note that although the medium is time-dependent in principle, for simplicity the optical potential is taken here to be time-independent and can be interpreted as a time-averaged quantity.
Whereas in the original DWEF formalism only a rotationally symmetric transverse density ρ(b) was needed (corresponding to central collisions), here we are interested in azimuthal anisotropy and so we need to consider a more general form. Specifically we take the modified Woods-Saxon profile from Ref. [20] with R ws = 1 2 (R 2 x + R 2 y ). Thus lines of constant density in the transverse plane form ellipses with semimajor to semiminor axis ratio Ry Rx . Lastly we must specify the fluid velocity u, for which we again defer to Ref. [20]. It is parametrized using a transverse fluid rapidity η t (b) The transverse direction is taken to be perpendicular to lines of constant density. It can be shown that the angle of such a fluid velocity, φ b , obeys [20] The transverse fluid rapidity η t (b) is first taken to have the same elliptic symmetry as the density, increasing linearly with the "radial" coordinateb ≡ (b cos(φ)) 2 . Then added to this is a term proportional to cos(2φ) representing the amount of elliptic flow built up before freezeout The momentum in these coordinates takes the form We choose to focus on data at midrapidity, Y = 0, and so In all, then, the parameters involved in this model are: ∆η, ∆τ , τ 0 , µ π , T , w 0 , w 2 , R x , R y , a ws , η f , and a 2 . We are interested in the effect of an optical potential like the one found to give the best fit in Ref. [16] and so we will keep all of these parameters fixed to those best-fit values, and only adjust Ry Rx and a 2 to give reasonable results for non-central collisions. It should be noted that the formalism developed is not strictly correct when the optical potential is complex. (See the discussion in Ref. [17].) We therefore also investigate the best fit values of Ref. [17] for a vanishing imaginary part of the optical potential.

III. CALCULATING v 2
This section outlines how the calculations are carried out. A set of coupled differential equations must be solved numerically to obtain the wavefunctions ψ (−) p , and then a fivedimensional integral must be performed (two of which can be done analytically with suitable approximations.) and Eq. (5) becomes or Decomposing ψ results in So the term in brackets vanishes identically for each m, and we must solve a set of coupled differential equations. In practice, every f m above a certain m max is set to zero, and a finite set of coupled equations is solved numerically.
The boundary conditions are the same as for the cylindrically symmetric case-far outside the medium one should have a canonically normalized plane wave plus an outgoing wave, i.e.
with J m and H (1) m Bessel functions and Hankel functions of the first kind, respectively. Details of this calculation can be found in Appendix A. The program used to calculate the wavefunctions was tested in part by comparing to a semi-analytic solution described in Appendix B.

B. Integration
Once the wavefunctions are found, a five-dimensional integral must be performed: The τ integral can be done analytically The η integral can also be done analytically with the following approximations (as in Ref. [16]) where the Bose-Einstein distribution is approximated by a sum over Boltzmann distributions truncated at some j max , and so Finally, then, for the numerator we have and similarly for the denominator. The final three integrals are done numerically.
More details of this part of the calculation can also be found in Appendix A.

IV. RESULTS
As previously mentioned, we would like to determine the effect of adding final state interactions to hydrodynamic fits. To gain insight into this, we consider an emission function  [17]. with parameter values taken from Refs. [16,17], which give the best description of the single particle data in general, and also with the imaginary part of the optical potential held at zero (see Table I. Also note that in both fits the chemical potential was fixed at the pion mass).
We must make alterations to this central collision model to approximate a more peripheral First, as in Ref. [16], we scale down R ws , a ws , and τ 0 by the number of participants to the 1/3 power, with N part taken from the Glauber model (with the same parameters used in Ref. [6]) for an impact parameter of 0 and 7 fm (N part = 377.5 and 171.544). Specifically R ws → 0.7688R ws . Then we adjust the ratio Ry Rx such that the spatial eccentricity has a value of 0.035. This is a reasonable value corresponding to the spatial eccentricity at freezeout of hydrodynamic fits of peripheral collisions with impact parameter ∼ 7 fm. Note that the brackets in Eq. 26 indicate a spatial average with weight given by Eq. 7, while the spatial eccentricity in hydrodynamic simulations are typically given with respect to , e.g., energy density. We nevertheless keep the eccentricity from Eq. 26 fixed at this value with an understanding that it is only a rough but still realistic guide to the shape.
Lastly we must specify how much elliptic fluid flow is built up in earlier stages of the collision, represented by the value of a 2 (recall Eq. 10). First we set a 2 = 0 and see what v 2 is generated by interactions with the optical potential in the absence of significant elliptic fluid flow ( Fig. 1(a)). The calculated elliptic flow coefficient v 2 is plotted as a function of momentum, along with the relevant experimental data. (Note that p in our calculation is the momentum of an asymptotically free pion detected far outside the medium, not the momentum of a particle as it is emitted inside the medium, and can therefore be compared directly to experiment.) Although we are only able to calculate up to a limited momentum, it is clear that final state interactions alone do not generate an appreciable value for v 2 for either the general best-fit parameters (Fit 1) or those with a vanishing imaginary part of the optical potential (Fit 2).
Next we increase a 2 such that the experimental value for v 2 is roughly obtained (Fig. 1(b)).
A value of a 2 = 0.11 was required for the parameters from Fit 1, while a 2 = 0.10 was sufficient to bring the emission function from Fit 2 into the physical regime. One can see that the optical potential has a small but non-negligible effect-it decreases v 2 on the order of 10-25% of its zero-interaction value with a slightly smaller effect as momentum increases. The integral over the azimuthal angle of the pion momentum, φ p is done as a sum using a simple trapezoid rule. This is because for each different value of φ p , a new set of differential equations must be solved. This also allows for the numerator and denominator of Eq. 20 to be solved simultaneously, with just a factor of cos(2φ p ) multiplied to the numerator when adding terms to the sum.
For each term in the sum, then, first the wavefunctions ψ (−) p are obtained. They obey a set of coupled differential equations of the form for all integers m. This set is truncated, since large m moments (f m for m > p ⊥ R ws ) contribute little to the wavefunction. Therefore, all f m for m greater than some m max are set to zero, leaving a finite (2m max + 1) number of coupled ordinary differential equations.
These are solved by calling a GSL solver. Using an embedded Runge-Kutta-Fehlberg method seemed to give the best performance. For these solutions, Eq. 7 is integrated numerically to find the moments U n . This is done with the GSL adaptive integration routine for oscillatory functions.
To match to the proper boundary conditions, one must find (2m max + 1) linearly independent solutions to this set of equations and take the correct linear combination of these solutions that matches the desired boundary conditions. The straightforward choice for these linearly independent solutions is to sequentially solve for the case where only one of the partial waves is non-zero near the origin. For example, for the n'th solution let: and then solve the set of differential equations up to some arbitrarily large b max far outside the potential. We can then match each partial wave in this n th solution to the form: The final wavefunction is then given by the linear combination of these solutions that matches the form of Eq. 19 at b max : This part of the program was tested with the trivial case of zero optical potential, in addition to comparing to a separately written program that calculates only the cylindrically symmetric case, as well as to the results of the semi-analytical test case described in Appendix B.
Once these wavefunctions are obtained and stored in memory, the integral over b and φ in Eq. 25 can be performed in addition to the sum over Boltzmann factors. The integrations are done with two GSL adaptive integration routines, one embedded in the other. The sum is done inside the argument of the integrals.

APPENDIX B: SEMI-ANALYTIC TEST CASE
To test the numerics, the case of a pion moving through an elliptically-shaped stepfunction potential was solved (semi-)analytically making use of elliptic coordinates. This can be compared to the case of a ws → 0 (see Sec. II).
We want to solve Eq. 5 with U(b) an elliptically shaped step function-a finite potential inside an ellipse in the transverse plane, with zero potential outside.
It is useful to change to elliptic (cylindrical) coordinates, denoted u and v. Think of u as a 'radial' coordinate that runs from 0 to ∞ and v as an 'angular' coordinate that runs from 0 to 2π Note the major and minor axes of the resulting confocal ellipses are reversed from the shape of the density used in the main calculation (which is larger in the y direction). This is to maintain consistency with the conventional definition of elliptic coordinates. At the end one can simply take φ p → (φ p + π) to match the usual convention in RHIC papers.
Consider the case The sharp boundary at u = u 0 is an ellipse with major and minor axes In this coordinate system the Laplacian is and so Eq. (5) becomes 1 a 2 sinh 2 (u) + sin 2 (v) or equivalently On the inside of the potential and on the outside separately, q(u) does not depend on u and these cases can be solved with separation of variables and the solutions patched together at Start by expanding ψ p (b) in terms of so-called elliptic sines and cosines of the 'angular' variable v. They are solutions of 'Mathieu's equation' [26]: The general solutions are called 'Mathieu functions,' usually denoted C(α, q, v) for solutions even in the coordinate v and S(α, q, v) for odd. Demanding periodicity of the variable v allows only certain discreet eigenvalues α (denoted here α n for the even functions and β n for the odd functions). This (complete) set of periodic solutions is commonly called elliptic sines and elliptic cosines: C(α n , q, v) ≡ ce n (v, q) S(β n , q, v) ≡ se n (v, q).
The general solution of Eq. (B6) can be written in terms of these elliptic sines and cosines: [f cn (u)ce n (v, q) + f sn (u)se n (v, q)] .
and the plane wave coefficients p n and s n are 1 p n = 1 π 2π 0 dv e ip·x ce n (v, q out ) (B20) 1 s n = 1 π 2π 0 dv e ip·x se n (v, q out ).
The coefficients Ce out n and Co out n , along with the analogous 'inside' coefficients are determined by matching boundary conditions.
To match at the u = u 0 boundary, project the 'inside' angular functions (e.g. ce n (v, q in )) in terms of the 'outside' ones (e.g. ce n (v, q out )). Ce in j Je j (u, q in ) B c jn ce n (v, q out ) + Co in j Jo j (u, q in ) B s jn se n (v, q out ) .
The coefficients (Ce in n , Co in n , Ce out n , Co out n ) can then be determined by demanding that ψ and its gradient be continuous at u = u 0 , which gives the following relations: j Ce in j Je j (u 0 , q in )B c jn = 1 p n Je n (u 0 , q out )ce n (φ p , q out ) + Ce out n He (1) (u 0 , q out )ce n (φ p , q out ) (B27) j Co in j Jo j (u 0 , q in )B s jn = 1 s n Jo n (u 0 , q out )se n (φ p , q out ) + Co out n Ho (1) (u 0 , q out )se n (φ p , q out ) (B28) j Ce in j Je ′ j (u 0 , q in )B c jn = 1 p n Je ′ n (u 0 , q out )ce n (φ p , q out ) + Ce out n He ′(1) (u 0 , q out )ce n (φ p , q out ) (B29) j Co in j Jo ′ j (u 0 , q in )B s jn = 1 s n Jo ′ n (u 0 , q out )se n (φ p , q out ) + Co out n Ho ′(1) (u 0 , q out )se n (φ p , q out ). (B30)