Backreaction from inhomogeneous matter fields during large-scale structure formation

We study how inhomogeneities of the cosmological fluid fields backreact on the homogeneous part of energy density and how they modify the Friedmann equations. In general, backreaction requires to go beyond the pressureless ideal fluid approximation, and this can lead to a reduced growth of cosmological large scale structure. Since observational evidence favours evolution close to the standard growing mode in the linear regime, we focus on two-component fluids in which the non-ideal fluid is gravitationally coupled to cold dark matter and in which a standard growing mode persists. This is realized, e.g. for a baryonic fluid coupled to cold dark matter. We calculate the backreaction for this case and for a wide range of other two-fluid models. Here the effect is either suppressed because the non-ideal matter properties are numerically too small, or because they lead to a too stringent UV cut-off of the integral over the power spectrum that determines backreaction. We discuss then matter field backreaction from a broader perspective and generalize the formalism such that also far-from-equilibrium scenarios relevant to late cosmological times and non-linear scales can be addressed in the future.


I. INTRODUCTION
The evolution equations used in cosmology are usually obtained under the assumption of a spatially homogeneous and isotropic distribution of matter. Although these symmetries remain preserved at all times in a statistical sense, the deviations from homogeneity and isotropy become locally sizeable during cosmological large-scale structure formation at late times. Here we ask whether this can have any influence on the cosmological evolution equations.
Because Einsteins gravitational field equations are non-linear, this question amounts to asking whether the Einstein equations for spatially averaged fields are affected by so-called backreaction effects, i.e., averages of second and higher orders in spatial inhomogeneities. Arguments that such backreaction effects are always irrelevant for cosmology have been made [1] and have been contested [2,3], the discussion focussing mainly on inhomogeneities of the metric. The smallness of the gravitational coupling constant G N makes it plausible that terms non-linear in inhomogeneities of the metric remain negligible for cosmological evolution [4,5].
However, the Einstein equations contain also the energy-momentum tensor of the matter distribution which develops large inhomogeneities. Starting from an approximation of matter fields as a non-ideal fluid, it was demonstrated in Ref. [6] how such inhomogeneities give rise to backreaction effects in late-time cosmology where non-linear terms in metric inhomogeneities and their temporal variations are negligible compared to those of inhomogeneities in the matter fields.
The purpose of the present manuscript is two-fold. First, we provide within the formalism of Ref. [6] the first explicit calculations of backreaction effects in dynamical model scenarios. Second, we generalize the formalism of Ref. [6] to systems in which a Navier-Stokes fluid description cannot be taken for granted.

II. BACKREACTION FORMALISM
As derived in Ref. [6], backreaction effects in a Friedmann-Robertson-Walker (FRW) universe can be encoded in a time or scale-factor dependent source term D(a) for the evolution of the spatially averaged energy densityε, 1 aε + 3H (ε +p eff ) = D(a) . (1) Here, the Hubble constant H =ȧ/a 2 is expressed in terms of the scale factor a(τ ) of the FRW metric ds 2 = a 2 −dτ 2 + δ ij dx i dx j and the dot denotes a derivative with respect to conformal time τ . The effective pressurē p eff =p +π bulk contains thermal pressure and a possible bulk viscous modification. Eq. (1) holds irrespective of whether a fluid description of the matter distribution with inhomogeneities is applicable (see section IV for details), but the explicit form of D depends on it. The spatial average of the trace of the Einstein field equations can be shown to be free of backreaction effects, Combined with the energy conservation (1), this closes the evolution equations for background fields. The two equations combine to arXiv:2106.10184v1 [gr-qc] 18 Jun 2021 Integrating over τ , one finds the modified Friedmann equation where we have assumed that the backreaction D(τ ) vanishes for sufficiently early times τ < τ i . We emphasize that (4) holds irrespective of whether a fluid approximation applies.
Further consequences of (4) are explored in Appendix A. In particular, we show that if backreaction arises in a purely radiative sector, the D-dependent contribution toε(τ ) in eq. (4) cancels exactly the second term on the right hand side, so that H(τ ) is in this sence D-independent. For other equations of state, however, the D-dependence does not cancel in eq. (4), i.e., the Friedmann equation is modified.

III. LINEARIZED INHOMOGENEITITES AND BACKREACTION IN COSMOLOGICAL FLUIDS
In this section, we first consider cosmological matter fields in dynamical scenarios for which fluid dynamics applies and where inhomogeneities can be followed by linearized fluid equations. The backreaction formalism of section II can then be elaborated on in detail.

A. The backreaction term for a fluid
For the case that matter fields can be described as non-ideal fluids with finite pressure and/or non-vanishing shear-and bulk viscositiesζ andη, the backreaction term D(a) can be written explicitly as [6] Here, θ ( q 1 )θ ( q 2 ) = δ (3) ( q 1 + q 2 ) P θθ ( q 1 ) defines the power spectrum P θθ in terms of a spatial average . . . of the expansion scalar θ(x) = d 3 qθ(q)e iqx . The spectra for vorticity w and pressure p are defined analogously.
To understand how to arrive at (5), one may consider the case of a fluid with pressure but with vanishing viscosities. The fluid dynamic equation for the energy density takes then the forṁ Here, the flow field u µ = (γ, γ v) with γ = 1/a √ 1 − v 2 allows for small inhomogeneities | v| 1 on top of the Hubble flow. The combination v · ∇ε + ε ∇ · v is a total derivative with vanishing spatial average and therefore, The spatial average δp ∇ · v over second order perturbations can then be expressed in terms of an integral over the power spectrum P θp ( q). This is the first line of (5). One obtains the last two lines in (5) by including the shear and bulk viscous contributions in the Navier-Stokes approximation in eq. (6) and following the same derivation.
B. Non-ideal fluids with standard growing modes The main aim of the present section is to calculate D explicitly for model scenarios of large-scale structure formation in which cosmological matter is treated in the fluid approximation. To this end, we expose first some qualitative considerations which inform our choice of model scenarios: The growth of large scale structure after photon decoupling is thought to follow a growing mode approximately ∝ a. This growing mode is realized, e.g., when cold dark matter (CDM) is described as an ideal fluid. However, the term D in (5) vanishes for a pressure-less and nonviscous ideal fluid. On the other hand, one-component fluids with pressure and/or viscosity are known to show a modified growth of cosmological large scale structure [7][8][9][10][11]. For instance, for a one-component fluid with pressure, the Jeans criterion implies that over-densities oscillate rather than grow if their mass is below the Jeans mass. This shows that it is not easy to generate any backreaction without modifying the late-time power-law growth of cosmological large scale structure.
In building cosmological models with backreaction, one bold approach could be to start without taking into account any constraint from the approximately linear growth of large scale structure; one may hope that at least some of the models constructed this way evade existing constraints or indicate modifications that will be supported by future observations. Alternatively, a very conservative approach would be to restrict the study to cosmological models which show a linear growing mode at all scales and which are therefore designed to pass observational constraints on large scale structure. In the following, we adopt this conservative approach.
Two-component fluids in which a fluid with pressure and/or viscosity is coupled gravitationally to a pressureless fluid show a standard linear growing mode in linear perturbation theory. One well-known case is a fluid with pressure coupled to a pressure-less fluid. This is realized, e.g., in the ΛCDM model where baryonic matter retains after recombination a small ionization and thus a small pressure. It is then the gravitational coupling to pressureless cold dark matter that ensures the growth ∝ a of the density contrast of baryonic matter on all scales, while the baryonic pressure does some work during the same epoch. We calculate the very small backreaction effect for this ΛCDM scenario in section III D before investigating other two-component models that allow for somewhat larger backreaction effects.

C. Gravitationally coupled two-component fluids
For a one-component fluid with pressure and finite viscosity, the equations of motion for the expansion scalar θ k and the energy density contrast δ k ≡ δε k /ε read at linear order and for subhorizon fluctuationṡ Here, the subscript k denotes Fourier modes of wavelength 2π/k and eqs. (10) and (11) are the Poisson equations for the two Newtonian potentials Φ k and Ψ k . We have parametrized the equation of state in terms of where H = H/a and where c s denotes the velocity of sound. We assume in the following that w, c s , c ad , ηH/(εa) are much smaller than 1, so that we can set them equal to 0, apart from terms in which they are multiplied by k 2 . In this approximation, the shear viscous difference between the two Newtonian potentials is neglected, Φ k = Ψ k . We now consider two sectors: one with pressure and/or shear viscosity that gives rise to a backreaction D (we denote it therefore by a subscript D), and one which is standard pressure-less and non-viscous cold dark matter (subscript C). The spatially averaged energy densities of both sectors evolve independent of each other, but it is the sum of their gravitational attractions that enters the trace of Einstein's field equations and that determines the expansion history. This is seen in the background equations a We denote matter inhomogeneities in the two sectors by δ k , θ k and d k , ϑ k , respectively. To linear order, these inhomogeneities are coupled gravitationally via the Poisson where Ω D ≡ ε D /(ε D + ε C ) is the energy fraction of the sector with pressure and/or viscosity, and Ω D + Ω C = 1.
Within this section, we consider a matter-dominated Einstein Universe with H ∼ a −1/2 and approximately constant energy fractions Ω D , Ω C . To linear order in inhomogeneities, this two-component system does not develop vorticity and with δp k = c 2 s δε D k = c 2 sεD δ k , its backreaction reads (21) It is helpful to rewrite the evolution of the linear perturbations (16)- (20) in the standard two-component fields One finds where the prime denotes a derivative with respect to ln a, and where we used (1 + H /H) = 1/2 for an Einstein Universe. In eq. (24), the velocity of sound enters in the combination c 2 s k 2 Similarly, the viscous corrections are written in terms of In the absence of pressure and viscosity, c 2 s = c 2 η = 0, the evolution equations have the growing mode φ 1 = φ 2 = φ 3 = φ 4 ∝ a. In more general two fluid models discussed below, this generalizes to As a consequence of the relative minus sign between the growth of δ k and the growth of θ k , it is then clear that the pressure-induced backreaction term in the first line of (5) is always positive if evaluated on the growing mode. This can be made more explicit by rewriting

D. Baryons coupled to CDM
The solutions of the linear evolution equations (23)-(26) depend on the material properties c 2 s and c 2 η and their scale-dependence. In this subsection, we consider first the case of normal baryonic matter gravitationally coupled to cold dark matter (CDM) shortly after the era of recombination. The square of the baryon sound velocity is [13] where µ mol is the mean molecular weight (including electrons and all ions of H and He), T b is the temperature of the baryon fluid and k B is Boltzmann's constant. For a rough estimate, we use µ mol = 1 g/mol and we assume a matter-dominated universe with T b = T • a −1 and T • = 2.7 K. This yields c 2 s a = 4 kBT• 3µ mol = 3 · 10 −13 which corresponds to c s ≈ 2·10 −5 at recombination (a ≈ 0.001). Consistent with the approximations leading to (8) - (11), this sound velocity is small enough to be neglected in the equations of motion except where c 2 s is enhanced by a factork 2 . The shear viscous contribution ∝ c 2 η of the baryonic fluid is much smaller, so that we can set c 2 η = 0 in the following.
Since T b ∼ a −1 in a matter-dominated Universe, eq. (31) allows us to set in the evolution equations (23)-(26) where α is constant in time but depends on the wavenumber k. The resulting system has been considered in chapter 8.3 of Weinberg's Cosmology book [14]. For Here, we have renamed the energy fraction Ω D in (23)-(26) by Ω B , since it is the baryonic fluid that exhibits pressure in the example of this subsection. For Ω B 1, the growing-mode solution of (34) is µ = 1, i.e., the small pressure contribution leaves the a-dependent growth of the density contrast (almost) unaffected, while the amplitude of the pressure-full fluid component is somewhat reduced by [14] To calculate from this solution the backreaction (30), we need to know the power spectrum P φ1φ2 (q). Since φ 2 = φ 1 = φ 1 on the growing mode, we have P φ1φ2 (q) = P φ1φ1 (q) = P δδ (q). The solution (35) relates the energy density fluctuations of the baryonic fluid to those of cold dark matter, δ k = ξ d k , and therefore P δδ (q) = ξ 2 P dd (q). The power spectrum P dd (k) of cold dark matter is given by . A very simple parametrization that is sufficient for our purposes is with q eq = 0.02 h Mpc and N = 10 2 Mpc h 3 . This allows us to write (30) as where the denominator of the integrand arises from the factor ξ 2 in P δδ (q) = ξ 2 (q) P dd (q). This factor ξ 2 is scaleindependent but it depends on the wavenumber q, see eq. (32). This is important: without the factor ξ 2 (q), the integral (37) would be logarithmically UV-divergent. The pressure-induced correction ξ 2 (q) is the physics effect that regulates this UV-divergence. To evaluate (37), we introduce the shorthand and we find with For the Hubble constant H 0 = 1 3000 h Mpc , the prefactor in eq. (39) is 6π N q eq H 2 0 = 4.2 · 10 −6 and A = 7.9 · 10 −10 . This yields F(A) = 10.2 and Given the smallness of the pre-factor λ, this backreaction is negligible and cannot affect the evolution (13) of cosmological background fields in a measurable way.

E. An explicit example of a modified Friedmann equation
How would the background fields and the Friedmann equation be modified in the presence of a more sizeable backreaction of the form (41)? As this would not be the physically realized baryonic matter, we switch notation which has the analytic solution Here,ε i,D is the energy density at the initial time τ i at which a i = a(τ i ). Since a ≤ 1, the factor e λ (a−ai) is O(1) throughout the evolution and the correction to the characteristic power-law decay of cold dark matter becomes negligible for very small λ. The factor e λ (a−ai) is always larger than unity which is consistent with the general expectation that the work done by backreaction increases the energy density of the corresponding matter component.
In the modified Friedmann equation (4), the energy density (43) enters in the combination For a a i , this allows us to write (4) in the form  To address this question, let us consider first a twocomponent system of the kind described in section III D, but where the sound velocity c 2 s of the pressure-full matter component takes a value different from that of the baryonic fluid. Since the backreaction D in (30) is proportional to c 2 s , one may try to increase the backreaction D by increasing c 2 s . However, there is a tight upper bound to any such effort: the amplitude ξ of the growing mode in (35) is reduced by increasing c 2 s . This reduction is larger for higher wavenumber q and it thus regulates the integral (37) that determines D. The resulting c 2 sdependence of D is encoded in the dependence of the factor A F(A) in (39). Plotting this factor in Fig. 1, we see that backreaction reaches its maximum for A ≈ 1 which corresponds to c 2 s a = 4 · 10 −4 and that it decreases for larger c 2 s . The value c 2 s a = 4 · 10 −4 still lies within the range of validity of our approximation scheme (c 2 s 1) in which only terms proportional to c 2 s k 2 were kept. For the baryonic fluid, we had determined A F(A) = 8.1 · 10 −9 . According to Fig. 1, max [A F(A)] = 0.2 is a factor 2.5 · 10 7 larger. Therefore, a dark matter component with minimal pressure (c 2 s 1) that couples gravitationally to cold dark matter and that participates in the growth of large scale structure with standard growing mode could exhibit a backreaction term as large as While this is much larger than the value (41) obtained for the baryonic fluid, it still falls short of being detectable.

Backreaction for pressure and viscosity with arbitrary scale dependence
So far, we have considered only dark matter components which exhibit negligible viscosity and whose sound velocity has a particular a-dependence c 2 s a = const. We next ask to what extent other material properties could yield backreactions that are even larger than eq. (46). To this end, we explore now cosmological two-component fluids in which one of the two components shows material properties of arbitrary scale κ s , κ η and with arbitrary power-law scale-dependence We introduce the shorthands φ i ≡ φin ainφ i for i = 1, 2, 3, 4. Eqs. (23)-(26) can be cast into two second-order differential equations for the growth of the density contrasts φ 1 andφ 3 , We consider again the case Ω D 1 in which the growing mode is (almost) unperturbed. In the limit Ω D → 0, eq. (50) has the solutionφ 3 = a, so that eq. (49) becomes For the viscosity-free case (κ η = 0) with constant c 2 s a (i.e. ν s = 0), the solutionφ 1 = ξ a of this equation reduces to eq. (35). A numerical solution of the general case (51) is possible. We find it more instructive, however, to consider two special cases that illustrate the generic properties of the solution.
Case I: κ s = 0, κ η = finite. For this case, it is convenient to use the evolution equationφ 2 =φ 1 = a dφ 1 /da and to rewrite eq. (51), The solution is The first term in the above expression can be neglected forφ 2i = a i 1. For ν η > 0, a useful approximate solution can be obtained if one observes that the evolution is characterized by two regimes. During the early evolution, for small a, the last term in the rhs of eq. (52) is negligible and the solution isφ 2 = a. At late times, for a ∼ 1, this last term becomes dominant and the solution isφ 2 3a 1−νη /(2κ ηk 2 ). A good fit of both regimes, giving also the correct order of magnitude for the short intermediate region, is given by the relatioñ Case II: κ η = 0, κ s = finite. For this case the analytical solution is complicated and not very useful because of the presence of strong oscillations for largek. However, if one averages over the oscillations, the average field φ 1 can be described by a simple expression. In analogy with the viscous case, for small a, the solution isφ 1 = a, while for a ∼ 1, the last term in the rhs of the averaged eq. (51) dominates and φ 1 3a 1−νv /(2κ sk 2 ). For ν s > 0, a relation of sufficient accuracy for our purposes is

Numerical Results
The denominators of (54) and (55) are scale-dependent generalizations of the factor ξ in (35) that regulates the UV-divergence of the integral over the power spectrum (37). In close analogy to the discussion in subsection III F 1, this limits the possible growth of backreaction with pressure or viscosity. To be specific, let us consider for the case I the integral that determines the viscous contribution to the backreaction in (30), Here, the prefactor 6πN H 2 0 q eq = 4.2 · 10 −6 is the same as in (39), and the solution ∝ C F(C) of the integral is of the same functional form as (39), but with the a-dependent argument We know from Fig. 1 that max [C F(C)] ≈ 0.2 and therefore, the viscous backreaction in (30) has a tight upper bound for the entire class of two-component fluid models considered here. This is consistent with the κ η -dependence of I η plotted in Fig. 2. We finally write the corresponding contribution to the backreaction (30) for the case II when viscosity vanishes but sound velocity has an arbitrary a-dependence, c 2 s ∼ a νs−1 .
While this integral differs somewhat from (56), the UV-cut on the spectrum P dd (q) is set by a similar qdependent denominator. This explains why the numerical evaluation of (59) in Fig. 3 yields also in this case a tight upper bound of the backreaction D ≤ O(10 −6 )Ω Dε H a.

IV. BACKREACTION FOR COSMOLOGICAL MATTER BEYOND THE LINEAR FLUID APPROXIMATION
Our discussion so far was based on describing matter in terms of ideal or Navier-Stokes fluid dynamics, and on following only scalar perturbations that were assumed to propagate in a linear way. These assumptions become questionable at late times and on small scales when cosmological structure formation becomes non-linear. Here, we provide a formulation of backreaction that remains valid beyond this linearized fluid-dynamic regime and that reduces to the formalism of section III in the limit in which Navier-Stokes fluid dynamics applies.
The need to go beyond a fluid dynamic formulation is particularly clear for the modelling of cold dark matter at late times and sufficiently small scales. At early times and very large scales, CDM is well-described by the single stream approximation which is equivalent to a description in terms of an ideal and pressure-less fluid. At late times, however, deviations from local equilibrium grow large, and the velocity dispersion induced by shell crossing modifies the dynamics [15][16][17]. In addition to scalar perturbations, also vector perturbations (e. g. vorticity) are generated by non-linear terms. Also the effect of dark matter self-interactions would show up in this regime. At even smaller scales, our current understanding of cosmological evolution is incomplete.
To understand cosmological evolution, we are particularly interested in the energy density ε(x). In a particle picture, it receives contributions from rest masses, from interactions and from the kinetic motion of particles (internal energy). For a one-component fluid with a single (at least approximately) conserved particle quantum number, the differential of this energy density is with temperature T , entropy density s, chemical potential µ and particle density n. Eq. (60) holds even if the fluid evolution is highly non-linear or turbulent and can be seen as a definition of T and µ. Such more general dynamical scenarios also allow for the definition of an entropy current s µ , as well as a particle number current N µ , so that one can write with fluid velocity u µ , heat current q µ and diffusion current ν µ . The heat current q µ and diffusion current ν µ are orthogonal to the fluid velocity, u µ q µ = u µ ν µ = 0.
Depending on the precise definition of the fluid velocity, there is an additional relation for q µ and ν µ (e.g., q µ = 0 in the Landau frame or ν µ = 0 in the Eckart frame), but we keep this frame definition open. Combining eqs. (60) and (61) and using ε+p = T s+µn, we obtain with Entropy production ∇ ρ s ρ and the change ∇ ρ N ρ of particle number vanish in global equilibrium, and so does the heat current q µ and the diffusion current ν µ . Therefore, the termσ is proportional to gradients and it vanishes in equilibrium or for an ideal fluid. Amongst all possibilities for gradient terms to contribute to (62), the combination ∇ µ u µ is distinguished in the cosmological context, because it is non-vanishing already for a homogeneous and isotropic fluid with Hubble expansion. We define therefore where π bulk is defined such that σ vanishes for a homogeneous and isotropic but expanding or contracting fluid. Alternatively, and equivalently, it is defined through the tensor decomposition of the energy-momentum tensor with respect to the fluid velocity u µ , T µν = εu µ u ν +(p+π bulk )∆ µν +π µν +q µ u ν +u µ q ν . (65) Eq. (62) becomes then By construction, the right hand side of eq. (66) vanishes for a homogeneous and isotropic fluid with FLRW expansion. In that case, eq. (66) gives the standard energy conservation law. However, for an inhomogeneous fluid for which homogeneity and isotropy are only symmetries in a statistical sense, the right hand side of eq. (66) can be, and will be, non-vanishing, for example due to terms quadratic in perturbations.
Assuming that metric perturbations are negligible, we can rewrite eq. (66) in terms of deviations v from the local Hubble flow. Using u µ = (γ, γ v), γ = 1/a √ 1 − v 2 , we finḋ As the combination v · ∇ε + ε ∇ · v is a total derivative with vanishing spatial average, the spatial average of eq. (67) reads where We stress again that D vanishes by construction in an exactly homogeneous and isotropic universe and that it starts quadratically in fluctuations. The form of (69) is more general than the derivation given in Ref. [6] and used in section III, since it does not assume that the fluctuations are small or that D is dominated by terms that are second order in fluctuations. For illustration, we consider some special cases for this general formula: (i) For an ideal fluid with pressure one has σ = π bulk = 0, and This is the case of our illustrative introductory example (7) that describes a modification of energy density due to work done in an imhomogeneous fluid by contraction against local pressure gradients. In Fourier-space, δp(x) = d 3 q δp(q)e iqx , it can be written as an integral over the power spectrum D = − 1 a d 3 qP θp ( q), which is the first term of (5).
(ii) In a first order gradient expansion (Navier-Stokes approximation) in the Landau frame one has q ρ = 0 and the local entropy production (71) For ∇ µ N µ = 0, this leads to the simple expression For exactly homogeneous and isotropic fields without inhomogeneities but with FLRW expansion, the only non-vanishing term in ∇ µ s µ is proportional to the spatially averaged bulk viscous pressureπ According to the definition (64), this term cancels in σ and therefore D = 0 in the absence of inhomogeneities. However, (73) appears also on the left hand side of (68). The resulting bulk viscous cosmologies without inhomogeneities have received much attention in model studies [7][8][9][18][19][20][21][22][23][24].
Keeping in (71) terms up to second order in inhomogeneities, one recovers the dissipative terms derived in [6] and given in the second and third line of eq. (5), These terms correspond to the contribution to energy evolution from the entropy production (dissipation) in fluid velocity perturbations.
(iii) The backreaction term D in eq. (69) can also be evaluated in far-from-equilibrium situations. One particularly explicit example is given by a system without conserved particle number, N µ = 0, in the Landau frame q µ = 0, when (63) reduces tõ In this case, the term σ √ 1 − v 2 in (69) is the spatial average over local entropy production times temperature and an inverse gamma-factor. In any kinetic theory description, even far-fromequilibrium, local entropy production ∇ ρ s ρ ≥ 0 follows from Boltzmann's H-theorem and is nonvanishing as a result of collisions away from detailed balance. By working in kinetic theory, or in non-equilibrium quantum field theory, one can therefore evaluate D beyond the simple first-order fluid approximation entering eq. (5).
(iv) From the decomposition (65) and the covariant conservation law ∇ µ T µν = 0, one obtains the alternative relation This is formulated in terms of the energymomentum tensor only and it does not make reference to entropy. For any form of out-of-equilibrium dynamics that follows the energy momentum tensor and that associates local flow fields (such as e.g. a kinetic theory formulation), all terms entering (76) are then given.
These remarks illustrate the general point that working within the framework of a fluid-dynamic evolution with linearized treatment of inhomogeneities does not exhaust scenarios in which backreaction can occur. We would find it particularly interesting to understand whether increased backreaction effects can be realized in far-fromequilibrium scenarios that are accompanied by significant entropy production or significant violations of particle number conservation.

V. CONCLUSIONS
This work analysed the combined evolution equations for the homogeneous part of the cosmological energy density and the growth of inhomogeneities, accounting for backreaction from the latter to the former. The evolution equations for the spatial averages (see eqs. (1) and (67)) and for the linearized inhomogeneities (see eqs. (8) to (11)) combine to a set of integro-differential equations that are difficult to solve in general. However, for perturbatively small inhomogeneities, one can solve the linearized evolution equation for perturbations in the absence of backreaction, and then determine the backreaction from the solutions. Based on this approximation, we gave in the present paper novel, explicit and compact expressions for how backreaction modifies Friedmann's equation in late-time cosmology.
A basic conclusion that can be drawn is that backreaction on the cosmological evolution is not necessarily linked to large metric perturbations, such as (Post-) Newtonian potentials. The original proposal of Ref. [6] had already emphasized this point. We discussed here an explicit realization, which has the appealing property that it makes use of two fundamental ingredients of the cosmological model: baryonic and dark matter. In section III we saw that the growth of structure proceeds in the gravitationally coupled systems of baryons and dark matter, despite the pressure that develops in the baryonic sector because of the small amount of ionization that it retains even after recombination. The simultaneous presence of pressure and significant structure growth makes the first contribution to the backreaction term of eq. (5) nonvanishing. The second important conclusion is that the sign of the backreaction term is such that the effect leads to the enhancement of the averaged energy density that drives the expansion. For backreaction through pressure this growth arises through the work done against perturbations. A similar conclusion can be reached for the backreaction through viscosity.
In the above example, the known baryon sound velocity is too small for the effect to be detectable. Replacing baryons by a hypothetical second dark-matter sector with nonzero pressure or viscosity enlarges the parameter range and leads to the enhancement of the effect by several orders of magnitude. However, the standard cosmological expansion is still not modified significantly. The analysis of section III reveals the reason by providing a refined general picture of what limits the maximal size of backreaction effects from inhomogeneities of matter fields. The backreaction terms in eq. (5) are integrals over measured or measurable power spectra. They are multiplied by the velocity of sound or by a viscous transport coefficient respectively. One may thus have expected that the maximal size of backreaction is set by the current phenomenological constraints on warm dark matter or viscous matter. But this is not the case. The growth of large scale structure is certainly known to put significant constraints on non-ideal fluid alternatives to CDM, but our choice of non-ideal fluids gravitationally coupled to CDM was designed to evade possible constraints from LSS. In doing so, we have encountered an even tighter and more generic constraint: any pressure or viscous correction gives rise to a physical UV-cut-off of the (logarithmically divergent) integral over the power spectrum that determines backreaction, see eq. (37), (56) and (29). This reduces backreaction below the conceivably detectable level for arbitrary choices of sound velocity or viscosity.
The analysis in section III has concentrated on cosmological times and regimes of scales where the evolution of inhomogeneities is close to linear. This allowed to make progress with analytical methods, but also implies -almost by construction -that backreaction terms remain small (and perturbation theory remains consistent). On the other side it is also clear that cosmological perturbation theory breaks down at late times and on small scales. The evolution of matter fields is then modified for example by the generation of velocity dispersion generated by stream crossing, or by non-linear terms in the evolution equations for scalar perturbations. Also vector and even tensor perturbations, like vorticity or shear stress, are generated by non-linear terms and can play an interesting role at small scales. Dark-matter self-interactions could modify the dynamics in the non-linear regime as well. Backreaction terms can be written as integrals over equal-time power spectra of inhomogeneities, and it is well possible that these integrals are dominated by small scales where non-linear effects are large. It would therefore be particularly interesting to extend our analysis into that regime. Of course this is technically challenging and needs a good understanding of non-linear physics.
In order to prepare for a future extension into the nonlinear and out-of-equilibrium regime, we have formulated the backreaction terms for the evolution of energy density in very general terms in section IV. This new formulation does not rely on the Navier-Stokes fluid approximation and can also be applied for example in the context of far-from-equilibrium quantum field or kinetic theory. It would be interesting to establish in the future an understanding of backreaction effects from the non-linear regime of matter field inhomogeneities to the overall cosmological expansion.
In this appendix, we further explore general consequences of the modified Friedmann equation (4). To this end, we first define the critical total energy densityε c by a(τ )D(τ ) .
(A1) When the integral over D is positive, the total average energy density todayε 0,c , for vanishing spatial curvature, is larger than the value usually taken for the critical energy density,ε If the dissipative term D is non-zero during some time interval, a part of the total background energy densityε c has been produced by dissipative processes. It might be useful to separate that part and to decompose the total background energy density additively as whereε is the conventional contribution of matter and radiation that does not result from backreaction andε D is the part of the internal energy that results from the dissipation of perturbations. The partε D has the evolutionε D + 3ȧ a (ε D +p D +π bulk,D ) = aD, for suitably defined pressurep D and bulk viscous pressurē π bulk,D . We assume now for simplicity the relation whereŵ D is some constant. After integrating, one finds ε D (τ ) = One first notes that one recovers the conventional Friedmann equation for D(τ ) = 0 or ifŵ D = 1/3. The latter condition would correspond toε D consisting of pure radiation. Ifŵ D is smaller than 1/3, the Friedmann equation is modified and in particular, today's value of the energy densityε not produced by backreaction becomes which for D > 0 is smaller than the conventional value because the integral is positive.