Kinetic description of particle emission from expanding source

The freeze out of the expanding systems, created in relativistic heavy ion collisions, is discussed. We combine kinetic freeze out equations with Bjorken type system expansion into a unified model. The important feature of the proposed scenario is that physical freeze out is completely finished in a finite time, which can be varied from 0 (freeze out hypersurface) to infinity. The dependence of the post freeze out distribution function on the freeze out time will be studied. As an example, model is completely solved and analyzed for the gas of pions. We shall see that the basic freeze out features, pointed out in the earlier works, are not smeared out by the expansion of the system. The entropy evolution in such a scenario is also studied.


I. INTRODUCTION
In the ultrarelativistic heavy ion collisions at RHIC the total number of the produced particles exceeds several thousands, therefore one can expect that the produced system behaves as a "matter" and generates collective effects. Indeed strong collective flow patterns have been measured at RHIC, which suggests that the hydrodynamical models are well justified during the intermediate stages of the reaction: from the time when local equilibrium is reached until the freeze out (FO), when the hydrodynamical description breaks down. During this FO stage, the matter becomes so dilute and cold that particles stop interacting and stream towards the detectors freely, their momentum distribution freezes out. The FO stage is essentially the last part of a collision process and the main source for observables.
In simulations FO is usually described in two extreme ways: either FO on a hypersurface with zero thickness, or FO described by volume emission model or hadron cascade, which require an infinite time and space for a complete FO. At first glance it seems that one can avoid troubles with FO modeling using hydro+cascade two module model [1], since in hadron cascades gradual FO is realized automatically. However, in such a scenario there is an uncertain point, actually uncertain hypersurface, where one switches from hydrodynamical to kinetic modeling. First of all it is not clear how to determine such a hypersurface. This hypersurface in general may have both time-like and space-like parts. Mathematically this problem is very similar to hydro to post-FO transition on the infinitely narrow FO hypersurface, therefore all the problems discussed for FO on the hypersurface with space-like normal vectors will take place here. Another complication is that while for the post FO domain we have mixture of non-interacting ideal gases, now for the hadron cascade we should generate distributions for the interacting hadronic gas of all possible species, as a starting point for the further cascade evolution. The volume emission models are also based on the kinetic equations [2,3,4,5] defining the evolution of the distribution functions, and therefore these also require to generate initial distribution functions for the interacting hadronic species on some hypersurface.
In the recent works [6] it has been shown that the basic assumptions of the Boltzman Transport Equation (BTE) are not satisfied during FO and therefore its description has to be based on Modified BTE, suggested in [6]. When the characteristic length scale, describing the change of the distribution function, becomes smaller than mean free path (this always happens at late stages of the FO), then the basic assumptions of the BTE get violated, and the expression for the collision integral has to be modified to follow the trajectories of colliding particles, what makes calculations much more complicated. In fact, in cascade model one follows the trajectories of the colliding particles, therefore what is effectively solved is not BTE, but Modified BTE. Once the necessity of the Modified BTE and at the same time the difficulty of its direct solution were realized, the simplified kinetic FO models become even more important for the understanding of the principal features of this phenomenon. Such simplified models were developed and studied in Refs. [5,7,8,9,10], but all these were missing expansion. Thus, the important question to be studied is whether the important freeze out features, pointed out in the earlier works, are not smeared out by the expansion of the system?
In this paper we present a simple kinetic FO model, which describes the freeze out of particles from a Bjorken expanding fireball [11]. Thus, this is a more physical extension of the oversimplified FO models without expansion [5,7,8,9,10]. Taking the basic ingredients of the FO simulations from Refs. [9,10] we can make physical freeze out to be completely finished in a finite time, which can be varied from 0 (freeze out hypersurface) to ∞. In the other words our freeze out happens in a layer, i.e. in a domain restricted by two parallel hypersurfaces τ = τ 1 and τ = τ 1 + L (τ is the proper time).
In Ref. [8] authors have also adopted kinetic gradual FO model to Bjorken geometry, but combined it with Bjorken expansion on the consequent, not on the parallel basis: system expands according to Bjorken hydro scenario, but when it reaches beginning of the FO process system stops expansion and gradually freezes out in a fixed volume. However, it was shown in [8] that such a scenario is not physical: the simultaneous modeling of expansion and freeze out is required in order to avoid decreasing total entropy. Now we propose such a generalized model.

II. FINITE LAYER FREEZE OUT DESCRIPTION
Let us briefly review gradual FO model, which we are going to generalize including expansion. Many building blocks of the model are Lorentz invariant and can be applied to both time-like and space-like FO layers, so at the beginning we will write these in general way. Starting from the Boltzmann Transport Equation, introducing two components of the distribution function, f : the interacting, f i , and the frozen out, f f ones, (f = f i + f f ), and assuming that FO is a directed process (i.e. neglecting the gradients of the distribution functions in the directions perpendicular to the FO direction comparing to that in the FO direction) we can obtain the following system of the equations [6,9]: The FO direction is defined by the unit vector dσ µ . FO happens in a layer of given thickness L with two parallel boundary hypersurfaces perpendicular to dσ µ , and s = dσ µ x µ is a length component in the FO direction. We work in the reference frame of the front, where dσ µ is either (1, 0, 0, 0) for the time-like FO, or (0, 1, 0, 0) for the space-like FO. The τ F O is some characteristic length scale, like mean free path or mean collision time for timelike FO. The rethermalization of the interacting component is taken into account via the relaxation time approximation, where f i approaches the equilibrated Jüttner distribution, f eq (s), with a relaxation length/time, τ th . The system (1) can be solved semi-analytically in the fast rethermalization limit [9].
The basis of the model, i.e. the invariant escape rate of the particles within the FO layer of the thickness L, for both time-like and space-like normal vectors is given as (see Refs. [9,10,12] for more details) where p µ is the particle four-momentum, u µ is the flow velocity. In fact the model based on the escape rate (2) is a generalization of a simple kinetic models studied in Refs. [5,7,8], which can be restored in the L → ∞ limit.
Here we will concentrate on the time-like case only, where the above Θ function is unity. The important feature of the escape rate in the form (2) is that physical freeze out is completely finished when s = L, i.e. it requires finite space/time extent. Furthermore, now we can vary this layer thickness, L, from 0 (freeze out hypersurface) to ∞ and study how the post FO distribution depends on the layer thickness. Interesting and unexpected result was found in [9,10,12], for both space-like and time-like FO layers, namely that if L is large enough, at least several τ F O , then post FO distribution gets some universal form, independent on the layer thickness.
Simple semianalytically solvable FO models studied in [5,7,8,9,10] are missing an important ingredient -the expansion of the freezing out system. The open question is whether the features of the FO discussed above will survive if the system expansion is included. In this work we are going to build a model, which includes both gradual FO and Bjorken-like expansion of the system, and answer this question.

III. BJORKEN EXPANSION WITH GRADUAL FREEZE OUT
Let us first remind the reader the basics of the Bjorken model [11]. Bjorken model is one-dimensional in the same sense as discussed before eq. (1) -only the proper time, τ = √ t 2 − x 2 , gradients are considered. Here the reference frame of the front, where dσ µ = (1, 0, 0, 0), is the same as the local rest frame (comoving frame), where u µ = (1, 0, 0, 0). The evolution of the energy density and baryon density is given by the following equations: where P is the pressure. The initial conditions are given at some τ = τ 0 : e(τ 0 ) = e 0 , n(τ 0 ) = n 0 . This system can be easily solved: where P = c 2 o e is the equation of state (EoS) with constant velocity of sound, c o .
It is important to remember that if we want to have a finite volume fireball, we need to put some boarders on the system. Here we assume that our system, described by the Bjorken model, is situated in the spacial domain is a space-time rapidity). Within these boarders system is uniform along τ = const hyperbolas due to model assumptions, while outside we have vacuum with zero energy and baryon densities as well as pressure. Thus, we have a jump, a discontinuity on the boarder, which stays there during all the evolution. Certainly, to prevent matter expansion through such a boarder (due to strong pressure gradient) some work is done on the boarder surface [13]. One can think about it as putting some pressure to the surface with the vacuum, exactly the one which would remove discontinuity, then work is done by the expanding system against this pressure. As the system expends the volume of the fireball increases as where A xy is the transverse area of the system. It has been shown in [13] that for such a finite Bjorken system the total energy is given by while the total entropy and the conserved charge appear to be Work done by the expanding system, W , is given by (8) One can easily check then the energy conservation: Applying our FO model to such a system, we obtain: where FO begins at τ = τ 1 and τ ′ = τ − τ 1 . Taking the fast rethermalization limit, similarly to what is done in [8], we can obtain simplified equations for f i , which is a thermal distribution f i (τ ) = f eq (τ ), and f f , and consequently for e i , n i and e f , n f : Now the idea is to create a system of equations which would describe a fireball which simultaneously expands and freezes out. Let us put our two components (e = e i + e f ) into the first equation of (3) and do some simple algebra: where last two terms add up to zero; the free component, of course, has no pressure. So far our eq. (14) is completely identical to the first equation of (3). The main assumption of our model is that our system evolves in such a way that eq. (14) is satisfied as a system of two separate equations for interacting and free components [15,16]: Similarly we can obtain equations for baryon density [15,16]: In all these eqs. (15-18) on r.h.s. there are two terms: the first one is due to expansion and the second one is due to freeze out.
As it was discussed above our volume is restricted in η space by two sharp boarders, |η| ≤ η R = const. All our interacting matter is maintained within these limits. However, the frozen out particles, which do not interact with system anymore, will not necessarily know about such a restriction and can, in principle, move faster than the boarder and thus will leave the volume under consideration. The momenta of the frozen out particles are distributed according to the distribution function, f f , which (as we shall see later) falls down exponentially with momentum, while only the particles strongly directed along the beam line can have high rapidity. So, from the particles emitted close to midrapidity only a very tiny fraction will have momentum rapidity higher than η R , and thus will leave our system volume. Particles freezing out from the system elements with higher space-time rapidity will have better, but still small, chance to leave the volume of consideration, and only the frozen out particles from the volume elements close to the |η| = η R boarders will disappear from our modeling with probability up to 50%, because f f is spherically symmetric in the comoving frame. For the reasonable η R boarders the amount of such particles can be neglected. It is, however, clear that our approach is best justified for the reaction region near midrapidity. In this work we will do only illustrative calculations, aiming for a qualitative understanding of the general features of the FO and checking the relation to the former models of this process [8,9,10], while for the more quantitative calculations and comparison with the experimental data we will restrict our modeling to midrapidity region [14].
Thus, finally, we have the following simple model of fireball created in relativistic heavy ion collision. Initial state, τ = τ 0 e 0 , n 0 Phase I, Pure Bjorken hydrodynamics, Phase II, Bjorken expansion and gradual FO, τ 1 ≤ τ ≤ τ 1 + L Solving Eqs. (15,17) we obtain: The difference with respect to the pure Bjorken solution (19) is in the last multiplier, and we see that, as expected, the interacting component completely disappears when τ reaches τ = L + τ 1 . With these last equations we have completely determined evolution of the interacting component [15,16]. Knowing e i (τ ) and EoS we can find temperature, T i (τ ). Due to symmetry of the system u µ i (τ ) = u µ (τ 0 ) = (1, 0, 0, 0). Finally, f i (τ ) is a thermal distribution with given T i (τ ), n i (τ ), u µ i (τ ). However, for us the more interesting is the free component, which is the source of the observables. Eqs. (16,18) give us the evolution of the e f and n f , and one can easily check that these two equations are equivalent to the following equation on the distribution function: The measured post FO spectrum is given by the distribution function at the outer edge of the FO layer, i.e. by f f (L + τ 1 ). Most of the observables will depend only on this momentum distribution. However the two particle correlations depend also on the space-time origins of the detected particles. Thus, in order to calculate two particle correlations we have to keep track on both the momentum and the space-time coordinates of the freeze out point of particles (this can be done, for example, in the source function formalism), i.e. we can use the full information regarding the space-time evolution of free component, eq. (22), and correspondingly can get some restrictions on it from the data.

IV. RESULTS FROM THE MODEL
Aiming for a qualitative illustration of the FO process we show below the results for the ideal massive pion gas with Jüttner equilibrated distribution [17]: in the comoving frame, where u µ = (1, 0, 0, 0). Here the degeneracy of pion is g = 3, while the baryon chemical potential in the case of pions is zero.
Contrary to the illustrative example in [15] here we do not neglect the pion mass. During FO the temperature of the interacting component decreases to zero, so at late stages of the FO process this new calculation is better justified. We will see below that T i falls below m π quite soon, and so the Jüttner distribution is a good approximation of the proper Bose pion distribution [16].
For our system we have the following EoS: where K n is Bessel function of the second kind, and a = m/T i . Eqs. (15) and (24) result into the following equation for the evolution of the temperature of the interacting component: .
(25) Furthermore, we have used the following values of the parameters: η R = 4.38, A xy = πR 2 Au , where R Au = 7.685 fm is the Au radius, τ 0 = 0.05 f m, T i (τ 0 ) = 820 M eV , During the pure Bjorken case the evolution of the temperature is govern by eq. (25) without the second (freeze out) term on the r.h.s.
In Fig. 1 we present the evolution of the temperature of the interacting matter, T i (τ ), for different values of FO time L, and Figs. 2 present the evolution of the energy densities for the interacting and free components. These satisfy the energy conservation E f (τ ) + E i (τ ) + W (τ ) = E(τ 0 ), where E f (τ ) and E i (τ ) are given by eq. (6) with corresponding energy density.
As it was already shown in Ref. [8,10], the final post FO particle distributions, shown on Fig. 3, are nonequilibrated distributions, which deviate from thermal ones particularly in the low momentum region. By introducing and varying the thickness of the FO layer, L, we are strongly affecting the evolution of the interacting component, see Fig. 1, but the final post FO distribution shows strong universality: for the FO layers with a thickness of several τ F O , the post FO distribution already looks very close to that for an infinitely long FO calculations, see Fig. 3 left plot. Differences can be observed only for the very small momenta, as shown in Fig.  3 right plot. So, the inclusion of the expansion into our consideration does not smear out this very important feature of the gradual FO.
Please note, that if one would look on our post FO distributions only in the medium momenta regions (with 0.5 GeV < |p| < 2.5 GeV ), then one could fit these spectra reasonable well with equilibrated distribution with temperature T = 160 M eV , see dashed line on Fig. 3. However, for low and high momenta such a fit would strongly disagree, this has been shown already in earlier FO model in Ref. [7b].
For our illustrative calculation we assumed η R = 4.38 what is one unit smaller than the original rapidity of the colliding nuclei for 100 GeV * A. The experimental data show a smaller plato in rapidity, but we used a classic onedimensional Bjorken model without transverse expansion which correspondingly requires a larger longitudinal size.
It is important to check the non-decreasing entropy condition [15,16,18] to see whether such a process is physically possible. Figs. 4 present the evolution of the total entropy, S(τ ), calculated based on the full distribu-tion function, f (p) = f i (p) + f f (p): The total entropy is then given by eq. (7). From this figure we can make an important conclusion, that gradual freeze out with rethermalization produces entropy. This is not new -in Ref. [3] the authors reached the same conclusion based on longitudially boost invariant hydrodynamical model with continuous emission (what means infinitely long FO in our language) with different escape probability [2]. In Ref. [2,3] P esc is given by the probability to escape without further collisions, calculated based on Glauber formula. In the case of complete rethermalization of the interacting component, which we also assume, their entropy increase can go up to 70%! This leads then to the corresponding increase in the pion production, which is a rough measure of entropy. Certanly such a hydro+continuous emission model leads the nonthermal final distributions, and authors found a reasonable agreement with NA35 transverse mass spectra. This model naturally explains the large number of the produced pions, which was a trouble for models with sharp FO hypersurface.
In our model, as well as in [2,3], there are two main sources of entropy production, namely the separation process (into interacting or free component) and the rethermalization process. Separation process leads to en- tropy increase because particles have now more accessible states being either interacting or free, not just interacting. The entropy increase in the rethermalization process is a manifistation of the H theorem. At the late stages of the long FO this entropy increase can be approximated analytically in our model.
At the late stages of the FO, when τ is close to the end of FO L + τ 1 , for long FO, i.e. L >> τ F O , the interacting component, f i , can be neglected next to free component, f f , (see eqs. (20,21) and Fig. 2 ), and so, using eq. (22) we obtain Then, putting eq. (28) into eq. (27) we get the following equation for the entropy evolution: If there is a baryon charge (or other conserved charges) in the system in general case the last term in the eq. (29) can then be related to n f : and for the baryon density, eq. (18), in this limit we have where n F O is the free baryon density at the end of FO. So, for the total entropy, eq. (7), we obtain: where N f (τ ) is the total number of frozen out baryon charges, which is in our limit, eq. (31), equal to the total number of frozen out baryon charges at the end of FO, N F O , which is, due to conservation law, the same as initial number of interacting baryon charges.
where S F O is the final entropy. Thus, we see that the total entropy increases during simultaneous expansion and gradual FO, as it is clearly seen on Fig. 4 for the L → ∞ curve. If there is a baryon charge in the system, then at the very end of the FO process entropy increases logarithmically, according to eq. (33).
On the other hand the fast FO in the narrow layer, which can be studied in our model, contrary to [2,3,4], does not lead to entropy increase. Already for the FO layer as thin as 2τ F O the entropy production is negligible, around 1%, see Fig. 4.

V. CONCLUSIONS
In this paper we presented a model, which describes simultaneously freeze out and Bjorken expansion, and thus, it is more a physical extension of the oversimplified FO models without expansion [5,7,8,9,10], which allowed us to study FO in a layer of any thickness, L, from 0 to ∞. Another important feature of the proposed model is that it connects the pre FO hydrodynamical quantities, like energy density, e, and baryon density, n, with  the post FO distribution function in a relatively simple way. As a zero order approach to the real physical system we studied the FO of the pion gas [16].
The results show that the inclusion of the expansion into FO model, although strongly affects the evolution of the interacting component, does not smear out the universality of the final post FO distribution, observed already in Refs. [8,9,10].
Another important conclusion of this work is stressing once again the importance to always check the nondecreasing entropy condition [15,16,18], since long gradual freeze out may produce substantial amount of entropy, as shown in Fig. 4.
This may have important consequence for QGP search. For many qualitative estimations it was assumed that all the entropy is produced at the early stages of the reaction and that the expansion, hadronization and FO go adiabatically, and thus number of pions can serve as a rough measure of entropy. However, if a non-negligible part of the entropy, say 10%, is produced during FO, then some estimations, for example of strangeness vs entropy (pion) production have to be reviewed.

VI. AKNOWLEDGEMENTS
Authors thank E. Molnar for fruitful discussions. This work was partially supported by Grant No. FIS2005-03142 from MEC (Spain) and FEDER.