{'fraction_non_alphanumeric': 0.04404029597931711, 'fraction_numerical': 0.014442364268520995, 'mean_word_length': 4.162933415158024, 'pattern_counts': {'":': 0, '<': 6, '<?xml version=': 0, '>': 11, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Analytical approach to relaxation dynamics of condensed Bose gases

The temporal evolution of a perturbation of the equilibrium distribution of a condensed Bose gas is investigated using the kinetic equation which describes collision between condensate and noncondensate atoms. The dynamics is studied in the low momentum limit where an analytical treatment is feasible. Explicit results are given for the behavior at large times in different temperature regimes.


INTRODUCTION
A distinctive feature of the quantum kinetic theory of the condensed Bose gas arises from the correlations between the superfluid component and the normal fluid part corresponding to the excitations. This causes the occurrence of number-changing processes determining the relaxation to equilibrium when excitations collide frequently. As a consequence, in the hydrodynamic regime, a collision integral C 12 describing 1 ↔ 2 splitting of an excitation into two others in the presence of the condensate is needed. The specific form of such a term depends on the dispersion relation ω(k) for the energy of quasiparticles and on the matrix element M of the effective Hamiltonian describing the interaction between them.
In a series of papers before the experimental realization of BEC, Kirkpatrick and Dorfman [1,2] have derived a kinetic equation in a uniform Bose gas which includes these processes, and they have computed explicit values for the transport coefficients of the twofluid hydrodynamics. More recently, Zaremba, Nikuni and Griffin [3] have extended the treatment to a trapped Bose gas by including Hartree-Fock corrections to the energy of the excitations, and they have derived coupled kinetic equations for the distribution functions of the normal and superfluid components. Many implications of their approach concerning transport coefficients and relaxation times have been thoroughly reported in Ref. [4], mainly when k B T ≫ gn c , where g = 4πam −1 is the interaction coupling constant in terms of the s-wave scattering length a, and n c is the condensate density. In this regime the dispersion law approaches the free particle law ω(k) ≈ k 2 /(2m). In the opposite low temperature limit k B T ≪ gn, the relevant part of the excitation spectrum is ω(k) ≈ ck, where c = gnm −1 is the sound speed near zero temperature and n is the particle density. The kinetic equation in this regime has been derived in Ref. [5], but its features have been less considered. Due to the absence of studies on the eigenvalue problem posed by the linearized kinetic equation, little is known about the time behavior of the solutions of kinetic equations in both regimes. This paper is devoted to deriving some results on the free evolution of initial perturbations of thermal equilibrium as described by the three-excitation collision term near the critical and near zero temperature. By making appropriate approximations of the kinetic equation we shall see that the problem becomes analytically tractable in some limits, and hence the asymptotic behavior of the solutions for large time may be evaluated. In the first regime, k B T ≫ gn c , we consider the small momentum limit where the equilibrium distribution function obeys the condition n 0 (k) ≫ 1. By approximating the full linearized kinetic equation in this limit we obtain an integro-differential equation which, after integration by parts, agrees with linearized evolution equation of wave turbulence [6,7] with the appropriate indices. This result is not a priori obvious when one starts from the full kinetic equation, because of the strong singularity of the integrand. At low temperature, k B T ≪ gn, we consider two opposite regimes n 0 (k) ≫ 1 and n 0 (k) ≪ 1. By performing a low momentum approximation, ck ≪ k B T , we will show that the wave turbulence framework still emerges but in a form not consistent with energy conservation, so an additional improvement is needed in order to restore this conservation law. In the opposite regime near zero temperature, where ck ≫ k B T , we shall consider an approximation based on the dominance of Beliaev damping processes.
A notable feature of the linearized collision terms in the wave limit is the homogeneous dependence on the momentum. Such a property furnishes a systematic way to compare the relative orders of different collision integrals at low momentum. This is achieved by simply considering the degrees of homogeneity, since their values depend on some indices related to the dispersion law, the scattering amplitude and the number of bodies in the collision. Thus one can see the dominance of C 12 over the binary collision term C 22 .
The plan of the paper is the following. In Sec.II we discuss the general form of the linearized kinetic equation as written in terms of irreducible components of the perturbation. Sec.III is first devoted to deriving the low momentum approximation of the evolution equation below the critical temperature. Then, once we show that the wave turbulence picture is recovered, we briefly review the general procedure for solving this kind of equation, we apply these techniques to the present context and we determine the asymptotic behavior of the perturbation for large time. Anisotropic perturbations are also discussed. In Sec. IV we derive the evolution equation in the thermal regime, with particular attention to the requirement of energy conservation, and we deal with the issue of the dominance of C 12 over C 22 . In Sec. V we consider the case of very low temperature and derive some asymptotic results specific for this regime. Section VI contains some concluding remarks.

II. LINEARIZED COLLISION INTEGRAL IN THE FREDHOLM FORM
The form of the three-excitation collision integral depends essentially on the dispersion law ω(k) and the matrix element M of the threeexcitation interaction. These quantities determine R by means of In this work we will consider two extreme regimes of the Bogoliubov dispersion law near the critical temperature and near zero temperature. In the first limit k B T ≫ gn c the dispersion law and the squared amplitude become [3] where n c is the superfluid density. Since in this regime the average value of the energy of a quasiparticle is of order k B T the term gn c is subleading in Eq. (4). In the opposite limit k B T ≪ gn we have [5] ω(k) ≈ ck + αk 3 , The amplitudes M can be derived heuristically by golden rule arguments from the effective Hamiltonians in these regimes [8]. Note that no distinction is made between n and n c in the scattering amplitude valid at T → 0, because we neglect the small low temperature depletion. Since we are interested in the low momentum region, we shall retain the first term of the dispersion law. Here we will consider the homogeneous regime where the distribution function is independent of the position. The linearization of the kinetic equation ∂ t n = C 12 [n] proceeds by the insertion of where n 0 (ω(k)) is the Bose distribution function without chemical potential with the energy given by the first term of Eqs. (4) or (6); the dimensionless function χ(k, t) parametrizes the departure from equilibrium for the noncondensate distribution function. The linearized kinetic equation adopts the form which must be completed with the linearized equation for the fluctuation of the condensate density δn c , The kernel L is symmetric in k and k ′ and, due to the Bose functions, falls off sufficiently as k, k ′ → ∞. We may then introduce the following scalar product and pose the eigenvalue problem The solution of the linearized equation would be formally written as the following series where the coefficients c j would be determined by the initial condition. In this work we do not intend to follow this procedure since it is uncertain that there exists a spectrum of discrete eigenvalues for this problem. This doubt is based on the observation that the lifetime of a long-lived well-developed sound mode of energy ω = ck varies continuously with the momentum, τ (k) ∝ k −δ , where δ = 1 at T = 0 [9] or δ = 5 at T = 0 [10]. So, this property does not seem to be consistent with the discrete character of the spectrum. Instead, we will try to find an expression for the action of the operator L[χ] at low momentum, and then solve the evolution equation found in this approximation. Formally, the operator L is decomposed in two types of contributions where the function U is given by and the explicit form of the coefficient N (k) will be not required. The conservation laws of energy and momentum imply the presence of zero modes of L proportional to the energy and the momentum, χ(k, t) ∝ ω(k), h(t) · k. In particular, Thus, the linearized collision operator can be written as which suggests to introduce a dimensionless variable A(k, t) defined as In terms of A, the linearized kinetic equation has the form The rotational invariance of the kernel U(k, k ′ ) may be exploited by expressing the perturbation as a superposition of angular momentum eigenstates Using the addition theorem for spherical harmonics one arrives at the equation where the operator K l adopts the Fredholm form Here U l (k, k ′ ) is the l-coefficient in the series of Legendre polynomials of U: Next we derive two approximations to K l [A](k, t) at low momentum valid near the critical temperature and low temperature. They must formally satisfy energy conservation which is achieved if Clearly, in view of Eq. (22), the symmetry of the kernel U 0 (k, k ′ ) is sufficient in order to accomplish this requirement, so such a property must not be lost in the approximated kernel.

A. The linearized equation at low momentum
In this section we consider the regime k B T ≫ gn c near the critical temperature. Here the dispersion is well approximated by ω = k 2 /(2m), and the squared amplitude by |M| 2 = 8n c a 2 m −2 . The low momentum approximation to U 0 (k, k ′ ) is obtained by the rescaling k → ǫk, k ′ → ǫk ′ and the expansion around ǫ = 0. This gives The kernel U l (k, k ′ ) has a non-integrable singularity proportional to |k − k ′ | −1 as k ′ → k (see Fig. 1). In order to manage carefully this singularity, it is convenient to integrate by parts in the exact expression of K l [A] given by Eq. (22). It helps to write and to perform the change of the order of integration. This yields where and θ(x) is the unit step function.
The behavior of the first term in the right hand side of Eq. (27) as k → 0 is derived by determining the leading behavior of J < l (k, q) as both arguments k, q → 0. With respect to the second term in Eq. (27) it is far from obvious that the limit q → 0 of J > l (k, q) must be incorporated since q runs from k to ∞. The key idea emerges by noting first that the considered scattering amplitude involves the assumption that all three momenta are of the same order. It is then expected to be adequate to take the limit q → 0 when k → 0 in J > l (k, q). Another more formal justification arises from the fact that J > l (k, q) diverges as ln(q − k) as q → k. Therefore, by assuming that ∂ q A(q, t) falls off sufficiently to insure integrability, the main contribution to the second term in Eq. (27) as k → 0 comes from the region where q is also close to zero. The leading behavior of the integral in the last term of Eq. (27) as k → 0 is easily determined since the singularities of (2l + 1)U 0 and U l cancel each other.
The evaluation of the integrals in that limit can be carried out analytically, and they assume the scaling form Explicit expressions of H <,> l (x) will be given below. In this approximation the action of the operator K l is written as where Note that I 0 = 0. Since the right side of Eq. (21) behaves as ∂ t A/(βω), the kinetic equation will be approximated by The integral operator that we have obtained in Eq. (32) is slightly different from those that appear in the linearization of the kinetic equation of wave turbulence. However, by assuming the validity of integration by parts and noting that n 0 ∼ (βω(k)) −1 at low momentum, the kinetic equation (21) adopts the form where h = 1. The kernel V l (x) may include a term proportional to δ(x − 1), but in our case this does not occur. Essentially, such an equation is the final result of the linearization procedure of the nonlinear equation of wave turbulence [7], which may be obtained from the original (1) by keeping the leading part of the integrand when n ≫ 1 [11]. Nevertheless, in the absence of the large momentum part of the distribution functions assuring the fall off of the integrand at large momentum, the naive justification of this procedure may be inappropriate and a more delicate approach is required. That is the reason why we have addressed the linearization of the full quantum kinetic equation preserving suitable integrability properties.
In the remainder of this section we shall focus on the solutions of Eq. (34) with K l approximated by Eq. (32).

B. Evolution of perturbations close to T c
In order to obtain the characteristic time scale of the kinetic equation, it is convenient to introduce the dimensionless momentum variable k = k β/m; the scale (m/β) 1/2 corresponds to the average value of the momentum carried by a quasiparticle in the weakly interacting gas. It turns out that the right side of Eq. (34) reads where τ −1 sets the inverse time scale corresponding to the scattering rate involving three excitations in the presence of the (small) condensate. It is given by Here β may be replaced by the inverse critical temperature of the ideal gas β c = mζ(3/2) 2/3 /(2πn 2/3 ) since the above expression is linear in the small depletion n c /n. We set t = t/τ for the reduced variable time. The explicit form of H 0 (x) is derived in the Appendix (the computations for the cases l = 1, 2 proceed in an analogous way). The results for l = 0, 1, 2 are Hereafter in this section we will drop the lines over the reduced momentum variables.
The most efficient way to manage the convolution integral of the kinetic equation is to use the Mellin transform with respect to momentum and the Laplace transform with respect to time. If F (s, λ) denotes the image of A(k, t), the evolution equation becomes or alternatively where Here W H l (s) and Ψ(s) are the Mellin transforms of H l (x) and the initial condition A(k, 0), while W V l (s) denotes the Mellin image of the kernel V l (x) in Eq. (35). We will assume that A(k, 0) has compact support to assure that Ψ(s) is an entire function. Let us first concentrate in the l = 0 case. From Eq. (38) we see that the Mellin function is given by where γ E is the Euler constant and ψ(z) = Γ ′ (z)/Γ(z) is the digamma function. Note the reflection property The fundamental strip where W V 0 (s) is analytic is Re s ∈ (−2, 4), and contains two simple zeros at s = 0, 2. It is important for what follows to introduce the winding number of the Mellin function W V l (s) as the increment divided by 2π of the argument of W V l (s) as s moves from σ − i∞ to σ + i∞ along the line Re s = σ. In the l = 0 case, one can check numerically that the winding number κ(σ) of the Mellin function vanishes when σ ∈ (0, 2). To proceed further one must solve Eq. (46) for F (s, λ). (In this paragraph we will assume an arbitrary non-negative value for h.) This can be accomplished by following the method of Ref. [12], briefly reviewed in the sequel. These authors have shown that when an interval of zero winding exists within the fundamental strip, (σ − , σ + ) ⊂ (a, b), then the equation (46) has a unique solution. It can be constructed using a particular solution B(s) of the homogeneous equation whose requisite properties have been described in Ref. [12]. The most important features of B(s) are the meromorphic property in the strip Re s ∈ (a, b + h), and its analyticity and absence of zeros in the strip Re s ∈ (σ − , σ + + h). The location of the poles and zeros of B(s) outside that interval is dictated by the set of zeros of the Mellin function within the fundamental strip. They determine the asymptotic properties of the solution. In particular, it turns out that the pole of B(s) with maximal real part determines the asymptotics of A(k, t) as k → 0. Now, by expressing the Mellin image of the perturbation as the function f satisfies the inhomogeneous difference equation with constant coefficients which has an inverse Mellin image given by The end step is to write the convolution corresponding to Eq. (50) and to invert the Laplace transform. This yields the solution in terms of the inverse Mellin transform of the base function Notice that when h > 0 the asymptotic expansion of the solution (53) as t → ∞ requires the determination of the leading behavior of Q(k) as k → ∞. According to Eq. (52), this comes from the zero of B(s) with minimal real part. Next we turn to the l = 0 case. Due to the presence of the zero of W V 0 (s) at s = 0, the function B 0 (s) has a sequence of simple poles at s = 0, −1, −2, . . ., with the behavior then where B(1) is an arbitrary constant. This leads to a behavior independent of k at low momentum Note that B(1) must cancel another similar factor arising from Q(q). On the other hand, the zero of B(s) with minimal real part is at s = 3, and it comes from the zero of W V 0 (s) at s = 2. Using B(s + 2) = W V 0 (s + 1)W V 0 (s)B(s) when s → 1 we obtain This zero determines the high momentum limit of Q(k) in the form so, the long time behavior of the solution is given by which possess the self-similarity property Finally, it is easy to derive in closed form the asymptotics in the more restricted regime k → 0, t → ∞ either from Eq. (57) combined with (59), or from Eq. (60) combined with (56). The result is where The source term δΓ 12 only receives contribution from the l = 0 component, or using the above low momentum approximation which contains the partial Mellin transform of A 00 (k, t) for s = 1. From the convolution in Eq. (53) one finds The substitution of the large momentum behavior of Q(k) given in Eq. (59) yields the asymptotics Combining this with Eq. (65) then gives the long-time behavior of the source term Since the energy of the perturbation is proportional to F (s = 3, t), we see that in accordance with energy conservation. What about the non-isotropic perturbations? The Mellin functions corresponding to the l = 1, 2 harmonics are Now W V 1 (s) has a double zero at s = 1, and the winding number κ does not vanish on the entire strip of analyticity; the corresponding values are 1, (−1) if Re s ∈ (1, 3), (Re s ∈ (−1, 1)). For l = 2, the winding number does not assume the value zero either, but now κ = 1 on the entire interval (-2,4), and the Mellin function does not vanish in the strip of analyticity. In these situations of absence of an interval of zero winding, Balk and Zakharov [12] have shown that the initial value problem for A(k, 0) is not well posed, signaling the non-uniqueness of the solution when κ > 0 on the entire interval (as in the l = 2 case), or a strong instability of exponential growth if κ ≷ 0 (as in the l = 1 case). In terms of the Mellin function, we have W V 1,2 (0) > 0, which at least is a sufficient condition [12] for the instability of the Rayleigh-Jeans spectrum n ∝ k −2 under anisotropic perturbations. Thus, it seems praiseworthy to assign this property to the Bose-Einstein equilibrium solution of the kinetic equation with a collision term given by Eqs. (4) and (5). Finally, it may be instructive to compare these findings with those of Ref. [13]. The application of the methods of wave turbulence theory to the study of the binary collision term C 22 , approximated when n 0 ≫ 1, showed the Mellin function W 22 (s) to be symmetric with respect to s = 3/2, reflecting the symmetry of the corresponding kernel U 0 . At the same time, it was shown that W 22 (s = 0) > 0 while W 22 (s = 1) = 0. As a consequence, such approximated form of C 22 destroyed energy conservation but was consistent with particle number conservation. By simply shifting the minimum of W 22 (s) one can get W 22 (s = 0) = 0 and W 22 (s = 1) = 0. So that another approximated form of C 22 , which conserves the energy but destroys the conservation of the particle number, appears as possible. However, given the reflection property of W 22 (s) and the location of its zeros, it does not seem possible to get a low momentum approximation maintaining both conservation laws.

IV. EVOLUTION OF PERTURBATIONS IN THE THERMAL REGIME βck ≪ 1
In this section we consider the regime ω(k) ≪ k B T ≪ gn, where the energy of quasiparticles is ω = ck, with c = gnm −1 . Their low-energy interactions are well described by the squared amplitude |M| 2 given by Eq. (7). By ignoring the next positive contribution (8m 2 c) −1 k 3 in the dispersion law all U l become equal to (2l + 1)U 0 , due to the collinearity of wave vectors in the collision integral. Therefore in this approximation, the kinetic equation is the same for all values of l. The low momentum limit of this kernel is A comparison between the exact and the approximated kernels is shown in Fig. 2. Our first task is to find an evolution equation consistent with energy conservation. When Eq. (72) is substituted into Eq. (22), the latter turns into where we have introduced an ultraviolet cutoff Λ ≫ k in the integral which gives the coefficient of A(k, t). We left for the moment this term undetermined. The well-defined convolution term in the above equation determines the (provisional) class of admissible solutions of the approximated evolution equation. Such a class would be formed by functions A(k, t) which grow as k → 0 more slowly than k −4 , and fall off faster than k −3 as k → ∞.
In accordance with that, the integral for the Mellin transform of the gain term is convergent for 3 < Res < 4. Note however that the analytic continuation of this integral to s = 0 is the value of the gain term formally evaluated for A(k ′ , t) = 1. Based on energy conservation, we may insist on considering the functions A = constant within the admissible class of solutions, so the loss term must cancel out this value 1/12. Therefore, the replacement Λ 3 /(3k 3 ) − 1/12 → −1/12 produces such a cancellation. The approximated evolution equation obtained in this way reads which again has the form of a linearized kinetic equation of wave turbulence, where the degree of homogeneity is now h = −4. In fact, we will find that the above equation still introduces a conflict with energy conservation. This defect is repaired by adding a source term independent of k to the right-hand side of Eq. (75). It has the form and produces a partial cancellation of the integral term of Eq. (75) with the result This implies that the class of solutions of the improved equation must be changed. Such a class is formed indeed by functions which fall off faster than k −4 as k → ∞. In order to understand the role of the source term and how it arises, it is convenient first to discuss the solution of the wrong evolution based on Eq. (75). By introducing the function we define the Mellin function corresponding to Eq. (75) as which has zero winding in the entire strip of analyticity, and W V (s + 7/2) = W V (7/2 − s).

The corresponding difference equation for the Laplace-Mellin image of A becomes
where v = 3/(32πmnβ 2 c 2 ). The method used to solve Eq. (80) is the same as before. We seek an appropriate solution of B(s − 4) = −W V (s)B(s) for Re s ∈ (3, 4), and we write the solution in the factorized form Re s+i∞ sin π 4 (s + 1) .

(83)
The other solutions of the homogeneous difference equation in terms of products and fractions of Γ-functions, such as are not appropriate for the present problem since they have a pole at s = 0 within the strip Re s ∈ (−1, 4). In accordance with the factorization of F (s, λ), one finds the solution as the convolution where Z(k) and Q(k) are given by the integrations in Eqs. (54) and (52) with σ ∈ (−1, 4). Thus, to derive the asymptotics of the perturbation as t → ∞, we must evaluate the low momentum behavior of Q(k). This is dictated by the pole of B(s) −1 at s = −4 which leads to If F (s, t) denotes the Mellin image of the solution A(k, t), the above result permits the calculation of the leading behavior of the quotient F (s, t)/B(s) as t → ∞. It reads We can see now the violation of energy conservation. In the approximation where n 0 ≫ 1 this requires for any time that or F (s = 3, t) = F (s = 3, t = 0) = Ψ(3), since Ψ(3) corresponds to the energy of the initial perturbation. But, in accordance with the above asymptotics, the quantity F (s = 3, t) assumes a value proportional to Ψ(−4)B(3)t −7/4 for large time, so dE/dt = 0. Another indication that Eq. (75) is conflicting comes from the fact that, in order to check energy conservation, the integration of both sides after multiplication by k 2 leads to which does not converge for the class of functions considered up till now.
To obtain an improved evolution equation, let us consider a factorized solution F (s, t) for a different base function respecting energy conservation. If we assume that Q(k) ∼ Dk α as k → 0, then the limit behavior of F (s, t) as t → ∞ is Since energy conservation requires that F (s = 3, t = ∞) = Ψ(3), we see that α = −3, and B(s) must possess a simple zero at s = 3 in order to cancel the pole of the gamma function. Furthermore, to accomplish the low momentum behavior of Q(k), the prescription for Re s in Re s+i∞ and the analog of (54) must be Re s > 3. Clearly, if the base function assumes the value which was discarded before, the requirement F (s = 3, t) = Ψ(3) is accomplished for large time. Some features of the solution generated by this base function are the following. The inverse images of Ψ(s)/B(s) and B(s) are simply and the behavior of Q(k) as k → 0 reads Hence the corresponding solution is written as which for large t assumes the value Similarly, the behavior or the source term as t → ∞ becomes so the insertion of Eq. (96) produces We can check directly energy conservation for any time, Applying the inverse Mellin transform we obtain The origin of the extra term proportional to k −4 is the pole of B(s) at s = 0. With this result at hand, the integral term of Eq. (75) evaluated for A(k, t) in Eq. (97) becomes The last term can be understood as a self-consistent source which restores energy conservation. In consequence, if one uses the partial cancellation produced by this term, the corrected evolution equation adopts the form where which gives rise to the same Mellin function as before , Re s > 4, but with a different strip of analyticity. The relevant interval of zero winding is now Re s ∈ (7, ∞) and the base function B(s) has neither zero nor poles in the corresponding strip Re s ∈ (3, ∞). We expect that the evolution equation obtained in this way captures the essential features of the low momentum regime at low temperature. Dominance of C 12 over C 22 . To conclude this section, it is interesting to observe that the collision term that we have considered in this and the previous section is the prevailing summand in the sum C 12 +C 22 in the low momentum regime, n 0 (k) ≫ 1. This is based on the rules of power counting for the degrees of homogeneity h of the linearized collision integral which have been given in Ref. [12]. These authors have shown that when the dispersion law has the form ω(k) ∝ k δ , and the scattering amplitudes of C 12 and C 22 scale according to the index h is given by where d is the spatial dimension, and ν is the exponent in a solution n 0 (k) ∝ k −ν of C[n 0 ] = 0. Near the critical temperature, the scattering amplitude involving four excitations comes from the term g(ψ † ψ) 2 in the effective Hamiltonian, where ψ † (ψ) creates (destroys) excitations, while terms as g √ n c ψ † ψψ give rise to the scattering amplitude involving three excitations. In both cases κ = η = 0 since M 4,3 ∝ g, g √ n c . The insertion into Eq. (111) of δ = 2, κ = 0, d = 3, ν = 2 produces h[C 22 ] = 0, while h[C 12 ] = 1. Thus, we can neglect the contribution of C 22 in the linearized description at low momentum.
At very low temperature, when the excitations are phonons, the effective description is made in terms of a Hamiltonian H[n(x), θ(x)] depending on the particle density n and the Goldstone mode or phonon field θ conjugate to n. Now, since each θ-field carries the In this section we consider the kinetic equation and its solutions when the term U 0 (k, k ′ ) of the collision integral (22) is approximated in the quantum regime k B T ≪ ck ≪ gn. This proceeds by replacing n 0 (1 + n 0 ) → e −βck on the right side of Eq. (21), while U 0 (k, k ′ ) is symmetrically approximated by its leading behavior as βck → ∞, (see Fig. 3) The relative size of the second term is O(e −βc(k ′ −k) ). The corresponding subleading contribution to the collision term may be evaluated by using the Watson's lemma [14] 9 16πmn Although the approximation that arises by ignoring this contribution is necessarily nonenergy conserving, it may be interesting to explore its consequences, since we would expect this drawback to have a negligible effect as T → 0. Therefore we consider the evolution equation where, for convenience, we have written the gain term as a convolution. Unsurprisingly we see that the loss term reduces to −Γ B (k)A(k, t) where Γ B (k) = 3k 5 /(320πmn) is precisely the decay width of Beliaev damping at T = 0 [10]. Although in this regime, n 0 ∼ e −βck ≪ 1, we are not in the vicinity of a power law spectrum as those appearing in the wave turbulence framework, the above linear operator is homogeneous of degree h = −5. Therefore, we can analyze the evolution equation with methods similar to the ones we have described before. We define the Mellin function corresponding to Eq. (114) , Re s < 4.
It remains to find the second factor Q(k) of the convolution which gives the solution of the initial value problem In particular, if the initial perturbation is A(k, 0) = C δ(k − k 0 ), the inversion of Ψ(s)/B(s) yields Q(k) in the form Since all singularities of B(s) and 1/B(s) are located in the half-plane Re s ≥ 0, the Meijer function G 0, 3 3,3 (z) vanishes if 0 ≤ z < 1. Some explicit limiting forms are where M ≈ 0.48. Consequently, the solution may be rewritten as It is interesting to find the leading behavior of the last integral as t → ∞. To apply the Watson's lemma, we use Eq. (121) and thus we obtain We see that this integral corresponds to a subleading term of O((Γ B (k 0 )t) −2 ) in the expression for A(k, t). Noting that when k > k 0 the second term in Eq. (123) proportional to Q 1 (k) is subleading with respect to the first one, we obtain the final result for the long time asymptotics of the solution or A(k, t) ∼ 6C 10 − 2 √ 5 k 0 cosh( √ 71π/10) exp −v B k 5 0 t , t → ∞, k ≫ k 0 . VI.

CONCLUDING REMARKS
In this paper we have studied the temporal evolution of a perturbation of the equilibrium distribution of noncondensate atoms in Bose gases as described by a kinetic equation including only collisions between condensate and noncondensate atoms. Due to the difficulty to address this problem without any simplification, we have considered some approximations to the kinetic equation in different temperature regimes. The evolution equations that arise by approximating the full kinetic equation when the energy of an excitation is small compared to the thermal energy turn out to have definite homogeneity properties, and may be regarded as kinetic equations for waves in the limit of large occupation numbers. This occurs near the critical temperature when the dispersion law is quadratic, and near zero temperature for Bogoliubov excitations of very low energy. Notably, by using the explicit expression that determines the degree of homogeneity, one may easily show the dominance of C 12 over the binary collision term C 22 at low momentum. This provides a justification to our approach (i. e. dealing only with C 12 ). In the opposite quantum regime where the occupation number of noncondensate atoms decreases exponentially with the momentum, we recover unexpectedly a kinetic equation of the same type of the wave-turbulence equations in the previous situations.
To obtain the solution of the initial value problem in all cases, we have applied the method developed by Balk and Zakharov [7,12] to analyze the behavior of weak-turbulent media near Kolmogorov spectra. This is a scarcely used technique due to the difficulties in the analytical calculation of the Mellin function and the base function, mainly for binary processes 2 ↔ 2 in turbulent media. Here, however, this task is made easier as we deal with splitting processes 1 ↔ 2 and the background distribution corresponds to thermal equilibrium.
Just below the critical temperature when the condensate is small, the isotropic perturbations and the time derivative of the condensate density vanish for long time according to power laws, while anisotropic perturbations are unstable. In the thermal regime at low temperature, the requirement of energy conservation leads to a modification of the naive form of the evolution equation, and indeed an additional source term is needed. All solutions of the correct equation decay exponentially for large time and the time derivative of the condensate density tends to zero according to a power law. In the quantum regime at very low temperature the large time behavior of perturbations for the non-condensed dynamics turns out to be exponentially decaying due to Beliaev damping processes.
Finally, one may contemplate the possibility of using a similar analysis to that we have presented in this paper to evaluate estimations of the transport coefficients.