{'fraction_non_alphanumeric': 0.03384279475982533, 'fraction_numerical': 0.011512504962286622, 'mean_word_length': 4.630902486728137, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Non-equilibrated post freeze out distributions

We discuss freeze out on the hypersurface with time-like normal vector, trying to answer how realistic is to assume thermal post freeze out distributions for measured hadrons. Using simple kinetic models for gradual freeze out we are able to generate thermal post FO distribution, but only in highly simplified situation. In a more advanced model, taking into account rescattering and re-thermalization, the post FO distribution gets more complicated. The resulting particle distributions are in qualitative agreement with the experimentally measured pion spectra. Our study also shows that the obtained post FO distribution functions, although analytically very different from the Juttner distribution, do look pretty much like thermal distributions in some range of parameters.


Introduction
Since the very beginning of heavy ion collision physics fluid dynamical models were widely used for simulation of such reactions. Their advantage is that one can vary flexibly the Equation of State (EoS) of the matter and test its consequences on the reaction dynamics and the outcome.
This makes fluid dynamical models a very powerful tool to study possible phase transitions in heavy ion collisionssuch as the liquid-gas or the Quark-Gluon Plasma (QGP) phase transition. For highest energies, achieved nowadays at RHIC, hydrodynamic calculations give a good description of the observed radial and elliptic flows [1,2,3,4,5,6], in contrast to microscopic models, like HIJING [7] and UrQMD [8].
Freeze out (FO) process is an important and necessary part of the fluid dynamical modeling. Generally saying, fluid dynamical models are determined by the initial state, EoS and the way FO is realized, in addition to well known fluid dynamical equations. FO describes the last stages of the reaction, when the matter becomes so dilute and cold that particles do not interact anymore, their momentum distribution freezes out, and particles stream freely toward the detectors. FO models allow the evaluation of two-particle correlation data, different flow analyzes, transverse momentum and transverse mass spectra and other observables. There are many different theoretical approaches for the freeze out problem. Sometimes these approaches are oversimplified to the extent that conservation of energy and momentum is violated (see Ref. [20,22] for an overview).
A basic standard assumption in fluid dynamical modeling of heavy ion collisions is that freeze out happens across a hypersurface in space-time. Such a hypersurface is an idealization of a layer of finite thickness (of the order of at least mean free path or collision time) where the frozen-out particles are formed and the interactions in the matter become negligible. The dynamics of this layer is described in different kinetic models such as Monte Carlo models [9,10] or four-volume emission models [11,12,13,14,15]. Under such an assumption FO can be pictured as a discontinuity in relativistic flow, where the kinetic properties of the matter, such as energy density or flow velocity change suddenly.
It took a long time to develop general theory of discontinuities in the relativistic flow. The story goes like this. In 1948 Taub [16] discussed discontinuities across propagating hypersurfaces with space-like normal vectors. If one applies Taub's formalism to FO hypersurface with timelike normal vector, one gets a usual Taub adiabat, but the equation of the Rayleigh line yields imaginary values for the particle current across the front. Thus, these hypersurfaces were thought unphysical. However more recently, Taub's approach has been generalized to these hypersurfaces as well [17] and the imaginary particle currents arising from the equation of the Rayleigh line were eliminated.
Thus, it is possible to take into account conservation laws exactly across any surface of discontinuity in relativistic flow. The corresponding equations are where dσ µ is the unit 4-vector normal to the discontinuity hypersurface, is the energy-momentum tensor, is a baryon current (usually in heavy ion collision modeling we neglect electric charge current) These consist of local thermodynamical fluid quantities: the energy density ǫ, pressure p, baryonic density n and the collective Since 1974 people were using the Cooper-Frye formula [18] to calculate final particle spectra. Much later it was realized that the Cooper-Frye FO description has a conceptual problem of negative contributions, when it is applied to the FO hypersurface with space-like normal vector [19,20,21,22,23,24,25,26]. 1 This leads to the non-equilibrated cut-off post FO distributions: The simplest example of such cut-off distribution is the cut Jüttner distribution [28], proposed in Ref. [19]. In Refs. [21,22,23,24] it was shown that one can obtain cut Jüttner post FO distribution in a oversimplified kinetic model. If one improves this model, taking into account rescattering and re-thermalization, or if one uses more re-1 To our knowledge the first attempt to modify post FO distribution has been done by Yu. Sinyukov in 1989 [27]. For the critical review of this work see [25].
alistic four-volume emission model [22], then final distri- We generalize the simple model presented in Refs. [21,22,23,24] for the FO hypersurface with time-like normal vector (the model is also not 1D anymore). The philosophy of this work is also similar -first we will show that in the oversimplified model the obtained post FO spectra are really thermal ones. Then in more realistic model, taking into account rescattering in the still interacting gas, we will obtain post FO distribution in the more complicated non-equilibrated form, which is in qualitative agreement with experimental data, mentioned above. Actually the main goal of this work is not quantitative calculations, but a qualitative understanding of the phenomenon in simple models.
2 Simple freeze out model As it was mentioned above our basic model is in some sense similar to the one discussed in Refs. [21,22,23,24], but on the other hand it is much more realistic and also habitual for the people working in ultra-relativistic heavy ion collision field. Instead of very special geometry (1D, infinitely long tube with particle source on the left end and vacuum on the right one) used in Refs. [21,22,23,24], this model can be combined with any relativistic hydrodynamical model. Below we will use proper time coordinate keeping in mind the Bjorken hydrodynamical model, as the simplest but still reliable one, although the similar treatment could be done for any pre FO state.
Let assume that at some moment τ = τ 0 our matter is still completely equilibrated and described by the thermal distribution f Juttner (τ 0 , r, p, T, n, u µ ). Starting from this moment matter starts to gradually freeze out. At each particular moment τ i the FO hypersurface will be τ = τ i , and we will follow this hypersurface, i.e. we will describe our system from the reference frame of the front (RFF), where normal vector is dσ µ = (1, 0, 0, 0). We can describe the FO kinetics of such a system assuming that we have two components of our distribution [12,13,14,21,22,23,24]: f f ree (τ, r, p), describing frozen out particles, and f int (τ, r, p), describing particles, which still interact (below we will use shorter notations, f f ree (τ ) and f int (τ ), but the meaning is the same). At the initial moment, τ = τ 0 , the distribution f f ree vanishes exactly and f int (τ 0 ) = f Juttner (τ 0 , r, p, T, n, u µ ). In this work we will use an extended definition of the classical Jüttner distribution [28], which is adequate at ultra-relativistic energies and allow for particle creation. Thus the particle density, which depends on T , is not the same as the conserved baryon charge density, n. Then, f int gradually disappears and f f ree gradually builds up as τ increases. The most simple kinetic model describing the evolution of such a system is the following: where τ f is a parameter with dimension of time determining how fast particles freeze out. It is of the order of collision time. In general case τ f may depend on particle density and momentum of the particle, but in our simple model we neglect such effects and consider τ f = const.
Our assumption becomes particularly bad for the large times, when the particle density is very low -we will discuss this problem separately in section 2.3. If this approximation is not done, i.e. the relaxation time is assumed to be momentum and density dependent, that would be a large step towards the full solution of the BTE, which is then a much more involved problem, and in most cases cannot be performed analytically. Nevertheless even this strongly simplifying assumption leads to nontrivial results as we shall see later.
The interacting component, f int , will gradually freeze out in an exponential rarefaction with time. Then At τ −→ ∞ the distribution f f ree will tend to the original Jüttner distribution, while f int will vanish. The frozen out particles keep their velocities and consequently the volume occupied by the frozen out matter will gradually increase with time.
This is a highly unrealistic model, similar to the one resulting in the cut Jüttner distribution in Refs. [21,22,23,24] for the FO hypersurface with space-like normal vector.
We can improve this model taking into account rescattering in the interacting component in the same way as it was done in the References above.

Freeze out model with rescattering
The assumption that the interacting part of the distribution, f int , remains after some drain the same Jüttner distribution with decreased amplitude is, of course, highly unrealistic. As suggested in Refs. [21,22,23,24] rescattering within this component will lead to re-thermalization and re-equilibration of this component. Thus, the evolution of the component f int is determined both by the drain term and the re-equilibration. If we include the collision terms explicitly into the transport equations (6) this in general case leads to a combined set of integro-differential equations. We can, however, take advantage of the relaxation time approximation to simplify the description of the dynamics. In this framework the two components of the momentum distribution develop according to the modified coupled differential equations: 2 Thus, the interacting component of the momentum distribution, described by eq. (9), shows the tendency to approach an equilibrated distribution with a relaxation time τ r . Of course due to the energy, momentum and conserved particle drain, this distribution, f eq (τ ) is not the same as the initial Jüttner distribution, but its parameters, n eq (τ ), T eq (τ ) and u µ eq (τ ), change as required by the conservation laws.
If we do not have collision or relaxation terms in our transport equation then the conservation laws are trivially satisfied, because the change of the full distribution collision or relaxation terms are present, these contribute to the change of T µν and N µ , and this should be consid-2 Notice that in the ultra-relativistic case the total number of particles is much bigger than the baryon charge and depends on T . Thus the collision term in system (9) gives a nontrivial and nonvanishing contribution. ered in the modified distribution function f int (τ ). In this case from the conservation laws we get: and where we used kinetic definitions of N µ and T µν . The above equations can be easily solved: and 2.2 Immediate re-thermalization limit As a first approximation to this solution let us assume that τ r ≪ τ f , i.e. re-thermalization is much faster than particles freezing out, or much faster than parameters, n eq (τ ), T eq (τ ) and u µ eq (τ ) change. This leads to This assumption may be unrealistic for large times, when particle density becomes very low. We shall discuss the problem of the large times in the next section.
For f eq (τ ) we can assume the usual Jüttner form at any Thus, we can write eqs. (12,13) in the following forms: and Using the expressions derived in Ref. [22] we can easily find the solutions for n(τ ), ǫ(τ ), p(τ ) and u µ (τ where ∆ µν (τ ) = g µν − u µ (τ ) u ν (τ ) is a projector to the plane orthogonal to u µ (τ ). Thus, Then, From eqs. (15,18) we obtain the expression for the baryon density Similarly, from eqs. (16,18,19) we obtain This last expression means that from the moment our gradual FO has started the pressure in the interacting component is fixed by the model and conservation laws, and one does not need in addition the EoS as usual. The EoS during FO looks like:

Freeze out of massless pion gas
To proceed further we need to make some assumptions about the type of matter we are dealing with. First of all we assume that our system is homogeneous (this assumption will be used for the next sections as well). Let us take the simplest case of massless pion gas. Then, EoS is and it will stay like this during all the FO process. Using eq. (19) we find temperature as a function of τ : Thus, we can write the expression for the Maxwell-Boltzmann as it should. For such a distribution Now, based on the second equation in system (9), we can find the distribution function of the frozen out particles. where is the Exponential integral function. We would like to underline the following points about these graphs. One may notice that for the small times τ For the very large times our distribution becomes very strongly peaked at low momenta. This is not a surprise, since our temperature is going down exponentially and for large times we keep only particles with very low momenta in our interacting component, which keeps gradually freezing out and cooling down.
As it was mentioned before our model becomes unrealistic for large times, but there is staightforward way to solve this problem. After some time the particle density will be so low that the particles will have no chance for further interaction. So, we have to stop at some moment and transfer the rest of still "interacting" particles (probably a few percents) to frozen out component. In any case, since our FO is a very fast exponential process, the uncertainty generated at large times will be rather small.

Entropy condition
To obtain a physically realizable result, we have to check the condition for entropy increase: where, for both the equilibrium and the non-equilibrium distributions the entropy current is given as This condition is not necessary to obtain a solution of the freeze out problem, but it should always be checked to exclude non-physical solutions.
In our case the original entropy density is given by the well known expression for the Maxwell-Boltzmann gas: where we have used eq. (28). The entropy density for τ = ∞, when all the particles are frozen out, can be calculated through eq. (33) using the distribution function given by expression (31). We have done this integration numerically and obtained for the ratio of entropy densities: for any initial temperature, T (τ 0 ) 3 .
Thus, we see that the entropy density is slightly decreasing during our FO process. Since we have a gradual FO this means, that in order to avoid a decrease of total entropy, our FO model has to be accompanied by an expansion, which is natural, as mentioned before, and which should lead to a volume increase by at least 6% during FO. This does not sound unrealistic.

Freeze out of nucleons
Now we would like to study the situation for baryon rich matter -in this case the distribution functions will be different, since we have to take care about baryon density conservation in addition. If our system contains some nucleons in addition to pions (for simplicity we assume that there are no other mesons in the system), then the situation is the following. The EoS at the pre FO state was more complicated then eq. (24), but for the FO times it is fixed by the eq. (22). We assume that the temperature is high and the chemical potential is not too large -such that 3 It is not a surprise that R is independent on T (τ0), since for the massless gas we have only one dimensional parameter -T (τ0) and correspondingly s ∝ T 3 (τ0).
our system is pion dominated (this is reasonable assumption for ultra-relativistic heavy ion collisions). For such a system the temperature will be still basically defined by pions and we can connect ǫ and T through equation (23) and correspondingly eq. (25) is also valid. 4 The expression for the distribution function of the interacting component now becomes (see for example [31]): where m is the mass of the nucleon, K 2 (x) is the Bessel function of the second kind. Of course, for the proper treatment of the nucleons we would have to use Fermi statistics, but in this work we restrict ourselves to MB distributions first of all to get an analytical/semianalytical results, and secondly because MB distributions are nevertheless frequently used to fit experimental data, see for example [30]. Since for the FO in heavy ion collisions T ≪ m, we can use the expansion of the Bessel function: If we would use Bose and Fermi statistics we could formulate this approach in the following way: and From the last equation using eqs. (25,20), we can find that and, thus, µ is decreasing faster than T , and therefore if the inequality µ ≪ T , holds for τ = τ0, it will be even stronger for higher τ .
Thus we obtain Now, based on the second equation in system (9), we can calculate the distribution function of the frozen out particles. where is the incomplete gamma function. 5 The nucleon distribution functions are presented in Again f f ree looks similar to some thermal distribution. In this case f f ree is saturating very fast -already the f f ree 5 Please note that the Exponential integral function, eq.
(30), is actually a special case of incomplete gamma function: (△τ = 3) is a good approximation of the final f f ree (△τ = ∞). This is due to the fact that f int is freezing out much faster for baryons since, in addition to the exponentially decreasing temperature, it now has the baryon density, n(τ ), in front, which drops down even faster than T .
Then during FO process energy density is given by the expression:

Conclusions
Our analysis shows that, although the thermal post FO distributions are not restricted by general theory for the FO hypersurface with time-like normal vector, it may be too optimistic to fit particle spectra with thermal distributions. In our simple kinetic models we are able to calculate post FO distributions analytically, and we may indeed reproduce thermal post FO distribution, but only in highly simplified situation. If we take into account rescattering and consequently re-thermalization in the interacting component, the post FO distribution gets a more complicated, non-thermal form. The resulting particle distributions are in qualitative agreement with the experimentally measured pions spectra.
To clarify the principal place of our work in the field we shall make a comparison with the cut Jüttner approach of Bugaev [19], eq. 5, which is the only other post FO distribution published up to recently, which satisfies the FO requirements for space-like FO surfaces. Afterwards, Bugaev's approach was improved for space-like surfaces in Refs. [21,22] and in subsequent publications by this research group. The present work completes this improvement work and extends it to time-like FO hypersurfaces also. Thus, more precisely, this approach is an improvement compared to Ref. [19] just as the works [21,22] .
This work actually completes the approach introduced in [21,22] by extending them to time-like freeze-out surfaces.
The approach is not the best of all possible approaches, it includes a few simplifying assumptions for the more transparent presentation, but it is an important advance compared to the cut Jüttner approach.
Another important point one can learn from our study is that although our post FO distribution functions are analytically very different from the Jüttner distribution, they do look pretty much like thermal distributions in some range of their parameters. Thus, observed spectra may easily be misinterpreted.