Quasiadiabatic modes from viscous inhomogeneities

The viscous inhomogeneities of a relativistic plasma determine a further class of entropic modes whose amplitude must be sufficiently small since curvature perturbations are observed to be predominantly adiabatic and Gaussian over large scales. When the viscous coefficients only depend on the energy density of the fluid the corresponding curvature fluctuations are shown to be almost adiabatic. After addressing the problem in a gauge-invariant perturbative expansion, the same analysis is repeated at a non-perturbative level by investigating the nonlinear curvature inhomogeneities induced by the spatial variation of the viscous coefficients. It is demonstrated that the quasiadiabatic modes are suppressed in comparison with a bona fide adiabatic solution. Because of its anomalously large tensor to scalar ratio the quasiadiabatic mode cannot be a substitute for the conventional adiabatic paradigm so that, ultimately, the present findings seems to exclude the possibility of a successful accelerated dynamics solely based on relativistic viscous fluids. If the dominant adiabatic mode is not affected by the viscosity of the background a sufficiently small fraction of entropic fluctuations of viscous origin cannot be a priori ruled out.


Introduction
The first releases of the WMAP data [1] presented convincing evidence that the initial conditions of the Einstein-Boltzmann hierarchy are predominantly adiabatic. This conclusion follows from the relative position of the first anticorrelation peak of the cross-spectrum between the temperature and the polarization of the Cosmic Microwave Background (CMB in what follows). The subsequent WMAP releases and the Planck explorer results [2] confirmed (and refined) the early determinations of the first WMAP data [1] so that today we can say, with a fair degree of confidence, that single field inflationary models lead naturally to adiabatic initial data of the CMB temperature and polarization anisotropies. Nonetheless every deviation from this simple paradigm leads necessarily to various entropic solutions (see e.g. [3,4]). The entropic modes can be explicitly constrained using CMB physics but their presence is not mandatory for a consistent explanation of the observational data. Conversely the presence of an adiabatic mode is essential and cannot be overlooked at least in the framework of the standard concordance paradigm. According to the current data [1,2], a small fraction of anticorrelated entropic modes in the presence of a dominant adiabatic mode may even improve the fit of the temperature autocorrelations accounting for potential large-scale suppressions of the corresponding angular power spectra.
The conventional distinction between the adiabatic and the entropic solutions, going back to the pioneering analyses of the temperature and polarization anisotropies [5,6], assumes an ambient fluid that is thermodynamically reversible but this hypothesis is not necessary and it can be relaxed by making the plasma viscous. The gauge-invariant perturbations of the viscous coefficients lead to new fluctuations modes of the predecoupling plasma [7]. The physical features of these viscous solutions differ from the four conventional nonadiabatic modes 2 . Since the inhomogeneities of the viscous coefficients cause entropic fluctuations of the spatial curvature, the rôle of viscosity at large scales must either be constrained by the initial data of the Einstein-Boltzmann hierarchy [7] or totally absent. If correct this conclusion would threaten the possibility of an accelerated phase only driven by the viscous coefficients. In this paper we shall therefore analyze the fluctuations induced by the viscous coefficients both at the linear and at the nonlinear level. It will be shown that the large-scale fluctuations induced by inhomogeneous viscosities are not necessarily entropic, as argued in the previous paragraph, but they can be very close to adiabatic at large scales (hence the terminology quasiadiabatic) provided the viscous coefficients solely depend on the energy density of the relativistic plasma. The evolution equations of the gauge-invariant curvature perturbations in the case of a relativistic, irrotational and irreversible fluid differ substantially from the ones valid in the perfect fluid case and derived by Lukash [8] even prior to the formulation of the Bardeen formalism [9]. The equations for the normal modes of the curvature perturbations driven by the viscous inhomogeneities reduce anyway to the corresponding expression valid for a perfect and relativistic fluid [8] when the viscous coefficients vanish exactly. The results obtained in perturbation theory (within a gauge-invariant approach) are confirmed by a fully nonlinear analysis where the curvature inhomogeneities are treated within the expansion in spatial gradients. In this approach the curvature inhomogeneities are not required to be perturbative.
If the inflationary phase is solely driven by the viscous coefficients the quasiadiabatic scalar mode is more suppressed than in the conventional case where inflation is driven by a single scalar field. Consequently the tensor to scalar ratio exceeds the observational limits and turns out to be excessively large if compared with a bona fide adiabatic solution. Viscous stresses have been widely studied already in the mid seventies since the presence of bulk viscosity in the relativistic plasma can influence the character of the cosmological singularity [10]. Barrow [11] presented detailed studies on inflationary Universes driven by a bulk viscosity coefficient (see also [12,13]). While the possibility of an inflationary phase only driven by bulk viscosity received various critiques also in the past (see e.g. second paper of Ref. [11]) one of the byproducts of the present analysis is a systematic strategy for a more concrete phenomenological scrutiny of the large scale inhomogeneities induced by the viscous coefficients.
The plan of this paper is therefore the following. In section 2 we shall describe the viscous fluctuations of a relativistic plasma in an explicitly gauge-invariant language. Section 3 is devoted to the normal modes of the system; we shall also investigate if and when the perturbative fluctuations of the bulk viscosity can become quasiadiabatic. In section 4 the themes scrutinized in section 3 will be examined in a fully nonlinear perspective by expanding the geometry and the hydrodynamical variables in spatial gradients. In section 5 the tensor to scalar ratio will be computed when the quasi-de Sitter phase is solely driven by bulk viscosity. Finally section 6 contains some concluding remarks.

Gauge-invariant fluctuations viscous coefficients 2.1 General considerations
In what follows we shall present the full governing equations for an irreversible fluid where the viscous coefficients can fluctuate in space and time. If the total energy-momentum tensor of the plasma includes the viscous contributions, the adiabatic limit, in a strict sense, is recovered when the viscosities are neglected and the total entropy four-vector is conserved. The total energy-momentum tensor of the fluid shall then be written, for the present ends, as: where, as usual, P ν µ = δ ν µ − u µ u ν , ξ denotes the bulk viscosity coefficient, η the shear viscosity and four-velocity obeys g µν u µ u ν = 1. We shall consider the situation where the bulk viscosity coefficient ξ is the sum of a homogeneous part (denoted by ξ(τ )) supplemented by the inhomogeneous contribution (denoted by δξ(τ, x)): Note that the fluctuations of the bulk viscosity are defined in a slightly different manner 3 in comparison with Ref. [7]. The fluctuations of the bulk viscosity coefficient would read η(τ, x) = η(τ ) + δη (τ, x). In what follows we shall only consider the fluctuations of the bulk viscosity coefficient. A potential homogeneous component of the shear viscosity coefficient would not contribute the background equations (2.9); indeed σ ν µ of Eq. (2.1) does not have a homogeneous component, as it can be explicitly verified. The shear viscosity does not contribute to the background so that, without loss of generality, we can take η → 0. The shear viscosity is therefore important only over small scales. Indeed across the matterradiation transition the shear viscosity coefficient η determines the optical depth, the Silk damping scale and, ultimately, the shape of the visibility function [5,14]. Over large scales (possibly exceeding the Hubble radius) the shear viscosity only couples to the traceless part of the extrinsic curvature 4 .
The viscous energy-momentum tensor at large-scales can be evaluated in the Landau-Lifshitz frame [15]. The covariant conservation of the total energy-momentum tensor (i.e. ∇ µ T µν = 0) can be projected along u ν and along P α ν ; the two obtained equations together with the covariant conservation of the particle current are given hereunder: Using then Eqs. (2.3) and (2.5) together with the first principle of thermodynamics, the evolution of the entropy 5 can be easily derived: where the right hand side of Eq. (2.6) comes directly from ∇ α u β T αβ /T ; in Eq. (2.6) µ = µ/T is the chemical potential rescaled through the temperature, s is the entropy density and ν α 3 In the first paper of Ref. [7] the bulk viscosity coefficient includes a supplementary scale factor. More precisely we have that a ξ old = ξ current where ξ old and ξ current correspond, respectively, to the previous and to the current definitions of the bulk viscosity coefficients. 4 To be precise we refer here to the extrinsic curvature of the spatial slices called K ij ( x, τ ); see, in this respect, the discussion of section 4. 5 The explicit form of Eq. (2.6) has been obtained by trading the term u ν ∇ µ T µν for (∇ µ u ν )T µν since, in the Landau frame, ∇ µ (u ν T µν ) = 0. is given by: where χ denotes the heat transfer coefficient. The adiabatic limit is recovered when the viscous contributions are neglected and the total entropy four-vector is conserved.

Metric fluctuations induced by bulk viscosity
We shall now derive the gauge-invariant system of metric fluctuations with the purpose of computing the curvature perturbations induced by the fluctuations of the viscous coefficients.
In section 4 this problem will be addressed in fully nonlinear terms, i.e. without relying on the separation of the various quantities into a background value supplemented by the corresponding fluctuations. For the time being the metric, the total energy density and the total pressure will be split as where g µν will be taken to be conformally flat, i.e. g µν = a 2 (τ )η µν where a(τ ) is the scale factor, τ is the conformal time coordinate and η µν is the Minkowski metric with signature mostly minus. In general δg µν to linear order can always be written as the sum the scalar, tensor and vector fluctuations, i.e. δg µν = δ s g µν + δ t g µν + δ v g µν , where δ s , δ t and δ v denote, respectively, the scalar tensor and vector fluctuations of the perturbed metric. Using Einstein equations and recalling the notations of Eqs. (2.1)-(2.2) the evolution for the homogeneous expansion rate are 6 : where P t = p t − 3ξ H/a is the shifted background pressure. In Eq. (2.9) the prime denotes a derivation with respect to τ while the overdot will denote a derivation with respect to the cosmic time coordinate t. The relation between H and H =ȧ/a is given, as usual, by H = aH. It is appropriate to remark, at this point, that the perfect fluid contribution is characterized, in the simplest case, by the barotropic index w = p t /ρ t and by the related sound speed c 2 s = p ′ t /ρ ′ . The scalar fluctuations of the conformally flat metric are parametrized by four independent functions so that the entries of the perturbed metric can be written as: 10) 6 The Planck length and the Planck mass will be defined, respectively, as ℓ P = /1 √ 8πG = 1/M P where M P = M P / √ 8π. Recall also that M P = 1.22 × 10 19 GeV where, as already mentioned, δ s denotes the scalar fluctuation of the perturbed metric. The tensor modes are immediately gauge-invariant and are parametrized in terms of a divergenceless and traceless rank-two tensor in three dimensions: 11) where, as already mentioned, δ t denotes the tensor fluctuation of the perturbed metric.
Finally the vector modes are parametrized as 12) where, as already mentioned, δ v denotes the vector fluctuation of the various entries of the perturbed metric.

Evolution of the tensor and vector modes
The evolution of the tensor and vector modes is completely standard. In particular the tensor modes obey, as usual, Equation (2.13), as expected, does not include a direct contribution of the bulk viscosity coefficient. The contribution is only indirect (i.e. through the scale factor). From Eqs. (2.12) and recalling that u 0 δ v u i = V i , the (0i) and (ij) components of the perturbed Einstein equations together with the evolution of the velocity give 15) where, for simplicity, we just selected the gauge W i = 0. Equation (2.13) and its solution will play a rôle in the determination of the tensor to scalar ratio to be specifically discussed later on. Equations (2.14) and (2.15) determine the rate of decrease of the vector modes which do not play a major rôle in the present investigation, exactly as in the other more conventional cases. Equations (2.13), (2.14) and (2.15) share a common feature: the effect of the viscosity enters only through the evolution of the background. This is in sharp contrast with what happens in the case of the scalar modes.

Evolution of the scalar modes
The scalar fluctuations of the effective energy-momentum tensor 7 are given by: where u 0 δ s u i = v i , and ∂ i v i = θ. Concerning Eqs. (2.16) and (2.17) two comments are in order: first, as expected, the bulk viscosity and its fluctuations affect the spatial components of the perturbed energy-momentum tensor; second the bulk viscosity couples both to the peculiar velocity and to the metric fluctuations. Consequently the inhomogeneities of the bulk viscosity cannot be rescaled away with simple redefinitions of the metric perturbations as in the case of Eqs. (2.13), (2.14) and (2.15). Equations (2.16)- (2.17) are written in general terms since no particular gauge choice has been imposed so far. The same strategy will be also used in all the other relevant equations with the aim of deriving a consistent gauge-invariant evolution involving also the viscous coefficients and their fluctuations. Since this procedure is a bit lengthy but only the main steps will be swiftly outlined by focussing on the gauge-invariant meaning of the bulk viscosity coefficients. Without imposing a specific gauge choice, the Hamiltonian and the momentum constraints stemming from the (00) and (0i) perturbed Einstein equations in the presence of the viscous energy-momentum tensor are given by: The trace of the perturbed spatial components of the Einstein equations is: Similarly the traceless projection of the spatial components becomes: where L ij = (∂ i ∂ j − δ ij ∇ 2 /3) and Π ij is the total anisotropic stress. Note that Π ij is gauge-invariant and this implies that also the left hand side must be gauge-invariant 8 . This observation essentially determines the form of the so-called Bardeen potentials where Ξ denotes the gauge-invariant fluctuation of the bulk viscosity coefficient defined as: The gauge-invariant expression of the viscosity fluctuation will be particularly important in section 3 when discussing the quasiadiabatic mode. In Eq. (2.22) we also introduced the gauge-invariant counterparts of δp t and θ namely where the subscript g in Eq. (2.24) recalls that the corresponding fluctuation is is invariant under infinitesimal coordinate transformations (for short gauge-invariant). The gaugeinvariant counterpart of Eqs. (2.18) and (2.19) becomes respectively The evolution of the total energy-momentum tensor implies the following pair of equations where the two source terms F ξ and G ξ are defined, respectively, as: Both F ξ and G ξ contain the bulk viscosity and its fluctuations. As we shall see in the following section the presence of spatial fluctuations in the bulk viscosity coefficient induces further source terms in the evolution of the spatial curvature. As we shall see these terms play effectively the same rôle of an intrinsic nonadiabatic pressure fluctuation.

Quasiadiabatic modes
From the governing equations of the previous section the gauge-invariant evolution of the curvature perturbations can be easily obtained. We shall show, as anticipated, that the spatial perturbations of the viscous coefficients act as a source term of the evolution equations of the curvature perturbations. They are physically equivalent to nonadiabatic pressure fluctuations. However, if the viscous coefficients depend solely on the energy density of the fluid, the source terms induced by the viscous coefficients can be neglected for typical scales larger than the Hubble radius. Over large-scales the evolution equations of the normal modes of the system reproduce the ones of the adiabatic modes. Conversely, inside the Hubble radius their evolution equations are very different from their adiabatic counterpart. This is the reason why these modes have been termed here quasiadiabatic.

Adiabatic and nonadiabatic fluctuations of the pressure
In terms of the Bardeen potentials the gauge-invariant curvature fluctuations are where the second equality follows from the momentum constraint in its gauge-invariant form using the notation Θ = ∇ 2 V g . In the comoving orthogonal gauge (where both β and the three-velocity vanish) R coincides up to a sign (which is matter of conventions) with the fluctuations of the spatial curvature. While in the comoving orthogonal gauge, R is related to fluctuations of the spatial curvature, in a different coordinate system R will have the same numerical value but will not necessarily be related to curvature fluctuations. Equation (3.1) can be complemented with the definition of the curvature fluctuation on uniform density hypersurfaces: Using Eq. (2.25) the difference of Eqs. (3.1) and (3.2) is proportional to the Laplacian of the Bardeen potential and can therefore be neglected at large scales: This equation generalizes the relation between R and ζ to the case when the bulk viscosity coefficient is non-vanishing and it implies that, up to Laplacians of Ψ, R ≃ ζ.
Even if the conventional terminology might suggest otherwise, the nonadiabatic modes arise in a globally inviscid fluid, as the preceding considerations illustrate. In the present paper we want to drop this hypothesis since the total energy-momentum tensor of the plasma could include the contributions of the shear viscosity, of the bulk viscosity and of the heat transfer. The adiabatic limit, in a strict sense, is recovered when the viscous contributions are neglected and the total entropy four-vector is conserved. Equation (3.2) accounts for the curvature fluctuations on uniform density hypersurfaces. If the total fluid contains a number of different constituent components (for instance two, the a-fluid and the b-fluid) we can decompose the total ζ as ζ While the weighted sum of ζ a and ζ b is related to the total ζ, the difference between them gives what we normally define as the the entropy perturbations [3,4]. More specifically, the relative fluctuations in the specific entropy ς can be written as 11 : The the second equality in Eq. (3.4) holds in the case of two barotropic constituents. Equation (3.4) applies, for instance, in the discussion of CDM-radiation isocurvature mode. Let us now define the quantities relevant to the evolution of curvature perturbations. The gauge-invariant pressure fluctuation δp g can always be split into the adiabatic contribution (containing the total sound speed of the system) supplemented by the nonadiabatic contribution (containing the entropy fluctuations) where the two subscripts imply that the two relative variations at the right-hand side should be taken, respectively, at constant entropy and energy densities 12 . In the simple case of two fluids the total speed of sound and the nonadiabatic pressure density variation are: (3.6) where the speeds of sound in the two fluids of the mixture have been explicitly introduced.
Recalling the connection between ζ and the weighted sum of ζ a and ζ b it is also possible where the speeds of sound refer to the inviscid contribution to the total energy-momentum tensor. Thus c 2 s b = w b and c 2 s a = w a .

Decoupled evolution of the curvature perturbations
The decoupled evolution of the curvature perturbations is obtained in two steps. The first step is to derive the first-order equation obeyed by the gauge-invariant curvature perturbations. The second step involves some lengthy but straightforward algebra to pass from the first-order (but still coupled) system to a second-order decoupled equation. The first order equation obeyed by R ′ can be obtained at least in two different ways either starting from the evolution of the metric perturbations or from the total velocity field. The first derivation consists in taking the difference of Eq. (2.22) and of the Hamiltonian constraint (first equation of (2.25)). This combination will lead directly to a term containing a term proportional to δp nad . The same result obtained with this procedure can be derived from Eq. (2.28). In this case the observation is that, thanks to the momentum constraint (second equation of Eq. (2.25)), ∇ 2 (R + Ψ) = HΘ: using this relation to eliminate Θ from Eq.
(2.28), the wanted equation can be immediately derived. In both cases the final result can be expressed in the following manner:

7)
12 To perform the variation at constant (total) energy density means that δρ g a = −δρ g b . Similarly, to perform the variation at constant ς means that δς = 0, i.e. from Eq. (3.4 where the total source term Σ R is defined as In Eq. (3.7) and (3.8) we traded the difference of the Bardeen potentials for the total anisotropic stress as it follows directly from Eq. (2.26). Concerning Eqs. (3.7) and (3.8) few comments are in order. We first note, as already mentioned, that the fluctuations of the viscous coefficients play the same rôle of the nonadiabatic pressure fluctuations δp nad . In the limit ξ → 0 the viscous coefficients do not contribute to the background but the fluctuations always affect the curvature perturbations. This is the situation leading to the viscous modes across matter-radiation equality discussed in the first paper of Ref. [7]. Equation (3.7) leads to second-order equation for R which will be the basis for our subsequent considerations. By taking the time derivative of both sides of Eq. (3.7) various terms will arise: the terms containing the Laplacian of Ψ can will be eliminated through Eq. (3.7) while the term proportional to the time derivative of the Laplacian of Ψ (i.e. ∇ 2 Ψ ′ ) can be expressed via the constraint (second equation of Eq. (2.25)) and in terms of Eq. (3.1). After this straightforward algebraic procedure the explicit form of the second-order equation becomes: It is clear from the previous expressions that the presence of a gauge-invariant fluctuation of the bulk viscosity coefficient produces a computable source for the evolution of curvature perturbations which is fully equivalent to a nonadiabatic fluctuation of the pressure. As a cross-check it is useful to remark that, in the limit ξ → 0, Ξ → 0 and Π = 0, Eqs. (3.8) and (3.9) reproduce the well known results firstly obtained by Lukash [8] in the absence of nonadiabatic pressure fluctuations (i.e. δp nad = 0). In this limit all the terms at the right hand side of Eq. (3.9) disappear and Eq. (3.9) becomes, as expected 13 The case of a single scalar field is implicitly contained in Eq. (3.10). Formally 14 the case of single scalar field is obtained by requiring c s → 1, Σ R → 0 and ξ = 0. In this case (ρ t + p t ) → ϕ ′ 2 /a 2 (where ϕ denotes the scalar field) and z t → z ϕ = aϕ ′ /H. 13 In Eq. (3.10) we have that P t coincides with p t since, in this case, the bulk viscosity coefficient vanishes. 14 From the purely algebraic viewpoint the situation is slightly more complicated. Indeed the scalar field has an effective sound speed c 2 ϕ = 1 + 2a 2 V ,ϕ /(3Hϕ ′ ). In this case δp nad = δp ϕ − c 2 ϕ δρ ϕ = −4V ,ϕ ∇ 2 Ψ/(3ℓ 2 P ϕ ′ H) (where V ,ϕ = ∂V /∂ϕ and V is the potential of the scalar field). These two modifications combine in Eq. (3.9) and lead to the standard form of the evolution equation for the normal modes.

Quasiadiabatic normal modes
Whenever Ξ = Ξ(ρ t , H) the contribution of the bulk viscosity coefficient and of its fluctuations rearrange and cancel so that the source term is negligible at large scales. In spite of this observation the evolution equation for the curvature perturbations does not reproduce exactly the canonical result for the standard adiabatic modes, hence these solutions have been named quasiadiabatic. Since this is a particularly relevant point of the discussion we shall now present a more detailed analysis.
The term Σ R appearing in Eq. (3.7) contains two kinds of contributions: the terms that do not vanish in the limit ξ → 0 and those that do vanish in the same limit. By separating the terms of different origins Eq. (3.7) can therefore be rewritten as We will now show that all the terms appearing in the second line at the right hand side o Eq. (3.11) are O(∇ 2 R) and therefore negligible for typical length-scales larger than the Hubble radius.
For this purpose let us notice immediately that the last term at the right hand side of Eq. (3.11) can also be written as ξ∇ 2 (R + Ψ)/[a(ρ t + P t )]. We shall therefore focus the attention on the first two terms appearing in the second line at the right hand side of Eq. (3.11) and remark that in the linearized treatment Ξ(ρ t , H) is equivalent to the case Ξ(ρ t ) since, according to the background equations, H = ℓ P ρ t /3. Since Ξ defines the gauge-invariant fluctuation of the bulk viscosity given in Eq. (2.23) we can say, by definition, that where the second equality follows from the first one by recalling that ξ ′ = (∂ξ/∂ρ t )ρ ′ t and from the gauge invariant definition of δρ g of Eq. (2.24). Thanks to Eq. (3.13) (and using the Hamiltonian constraint to express δρ g ) we obtain: (3.14) Using now Eq. (3.1) we can immediately obtain, from Eq. (3.14) the following result Inserting now Eq. (3.15) into Eq. (3.11) and recalling the momentum constraint to eliminate Θ we finally arrive at the equation where all the terms containing the Laplacian of Ψ have been collected: the quantity appearing in the round bracket is actually P ′ t /ρ ′ t so that Eq. (3.16) becomes: where c 2 ef f = P ′ t /ρ ′ t is just an auxiliary variable 15 . From Eq. (3.17) the second-order form of the evolution equation of R can be written as: In the case δp nad → 0 (absence of nonadiabatic pressure fluctuations) and Π → 0 (absence of anisotropic stress) we also have Σ R → 0. When the preceding conditions are all met in the limot ξ → 0, the result of Eq. (3.18) reproduces the results of Lukash [8] since, in this case, c 2 ef f coincides with c 2 s . Whenever ξ = 0, however, the situation is totally different and Eq. (3.18) describes, as anticipated, the evolution of the quasiadiabatic modes.

Evolution of the quasiadiabatic normal mode
To analyze the evolution of the quasiadiabatic modes it is natural to set to zero both the anisotropic stress and the nonadiabatic pressure fluctuations. Going to Fourier space Eq. (3.18) can also be written as: Equation (3.19) is homogeneous but it is non-standard insofar as the pump fields get corrected both inside and (partially) outside the Hubble radius. We shall now argue that the coefficient of the third term of Eq. (3.19) is positive semi-definite if the inflationary phase is driven by the bulk viscosity coefficient. To demonstrate this point it is useful to reverse the question and demand that the coefficient of the third term in Eq. (3.19) is positive semi-definite; such a request implies: 15 It must be stressed that c 2 ef f does have the physical meaning of a sound speed only in the limit ξ → 0 since, in this limit, c 2 ef f coincides with c 2 s . Whenever ξ = 0 the effective sound speed is not only given by c 2 ef f , as we shall demonstrate in a moment.
Assuming now that the inflationary phase is triggered by a dynamical bulk viscosity coefficient we will haveḢ = − 3 2 (1 + w)H 2 + 3 2 ℓ 2 P Hξ. (3.21) From Eq. (3.21) we can argue that ξ is always expressible as where ǫ = −Ḣ/H 2 is the standard slow-roll parameter. During inflation ǫ ≪ 1 and therefore we shall have, in the first approximation, that ξ ≃ H M 2 P (1 + w). The inequality of Eq. (3.20) becomes immediately (3.23) Using now Eq. (3.13) and the fact that c 2 ef f is asymptotically constant during inflation it is possible to prove that Eq. (3.20) implies that 16 ǫ ≤ (1+w)/(1+w −2c 2 s ): all the terms at the right hand side of the previous inequality are O (1). Therefore the inequality simply requires that ǫ < 1 which is always true since, by definition, ǫ ≪ 1 during the slow-roll phase.
The term containing the first derivative of R can be eliminated and therefore Eq. (3.19) can be written, in Fourier space, as where the variable q is implicitly defined as . (3.25) If ξ is assigned with the constraint that asymptotically there is a slow-roll phase of quasi de-Sitter type we must always demand, The solutions pinned down by the conditions (3.26) are clearly not the most general ones compatible with the condition that the bulk viscosity coefficient depends solely on the energy density of the plasma. At the same time these solutions will be useful for a specific comparison of the quasiadiabatic mode with the genuine adiabatic solution to be discussed in section 16 As it ca ne explicitly verified from its definition and from the asymptotic expression of ξ we have that c 2 ef f → c 2 s − γ − 1 up to slow-roll corrections which are subleading.
5. Note in particular that there are solutions of Eqs. (2.9) where the slow-roll parameters of Eq. (3.26) are time-independent 17 . Even if the equations for the quasiadiabatic modes are completely general we shall find it convenient, for the sake of simplicity, to assume that the slow-roll parameters are constant at least approximately. Since to leading order in the slow-roll parameters ξ ∝ HM 2 P the explicit form of C(kτ ) becomes It is clear that during the inflationary expansion (or during the fully developed inflationary phase) the numerical value of C 2 (kτ ) is dominated by the third and the fourth terms of the previous expression. The equation to be solved becomes, therefore, This equation can be studied in two opposite regimes, namely The solution of Eq. (3.29) is standard and describes the regime when the relevant wavelengths are larger than the Hubble radius: .

(3.31)
For a proper normalization of the whole solution (and for a correct determination of the scalar power spectrum) we need also to solve Eq. (3.30) and, in this case, the solution can be solved in terms of parabolic cylinder functions. In the case of Eq. (3.30), however, the equation gets even simpler: 3.32) whose solution is given by The solution parametrized by the coefficient C 2 (k) is physically unacceptable. Taking into account that kτ > 1 and that ǫ < 1 this solution describes the situation where starting from perturbatively small initial data the curvature perturbations quickly become O(1) and rapidly jeopardize perturbation theory. We are therefore left with the solution parametrized by C 1 (k) describing a the exponential suppression of the curvature perturbations in the regime kτ > 1. Up to a phase we shall normalize |C 1 (k)| = 1/ √ 2k. This choice will be useful when comparing the quasiadiabatic power spectra with their adiabatic counterpart. 17 These solutions are well known (see e.g. [11]). In the power-law case a good example is the standard solution a(t) ∝ (t/t 1 ) α with α ≫ 1 and ξ ∝ HM 2 P . In this case, for instance, ǫ → 1/α.

Nonlinear curvature perturbations 4.1 ADM decomposition and normal coordinates
In section 3 we have shown that the gauge-invariant inhomogeneities of the viscous coefficients provide a supplementary source term in the evolution equations of curvature perturbations. If the bulk viscosity coefficient only depends on the energy density and on the Hubble rate the resulting curvature perturbations are approximately conserved over largescales even if the evolution of the corresponding normal modes differs substantially from the perfect fluid case. It is desirable to scrutinize the validity and the implications of this result without relying on the perturbative expansion. To achieve this goal we shall study the effects of the bulk and shear viscosity within the expansion in spatial gradients. This technique has been used for the first time in the analysis of general relativistic singularities [16] and subsequently exploited in a variety of contexts ranging from inflationary models [17,18] to large-scale structure [19].
In the Arnowitt-Deser-Misner formalism [20,21] (ADM in what follows) the line element is expressed in terms of the conventional (3 + 1)-dimensional decomposition: (4.1) where N = N(τ, x) denotes the lapse function, N i = N i (τ, x) is the shift vector and γ ij = γ ij (τ, x) is the three-dimensional metric. In the ADM variables of Eq. (4.1) the extrinsic curvature of the spatial slices will be denoted by K ij (τ, x) and it is defined as where ∇ i are the covariant derivatives defined with respect to γ ij . The other standard notations for the traces are K = γ ij K ij and TrK 2 = K j i K i j . The traceless part of the extrinsic curvature will be denoted by K j i = K j i − δ j i K/3. For the sake of simplicity the shift vectors will be assumed to vanish (i.e. N i = 0) and in this case the coordinate observers coincide with the normal observers 18 .

Nonlinear evolution, curvature modes and viscosity
In the viscous case the nonlinear generalization of the curvature perturbations on comoving orthogonal hypersurfaces is: where γ = det[γ ij ] and u i is the spatial component of u µ . The nonlinear generalization of the density contrast on uniform curvature hypersurfaces becomes instead where ρ and P = p + 3ξK are now nonlinear quantities; the expression of P and ρ holds to lowest order in the gradient expansion but are not necessarily homogeneous and this is why we distinguished them from their corresponding background values denoted, respectively, by P t and ρ t in sections 3 and 4. Equations (4.2) and (4.3) define a set of nonlinear variables which are also gauge-invariant. Both R i and ζ i do not depend on the choice of time hypersurfaces and are exactly invariant for infinitesimal coordinate transformations in the perturbative regime. Equations (4.2) and (4.3) correspond, in linear theory, to the variables R and ζ. Indeed in the conformally Newtonian frame where the gauge freedom is removed and the coordinate system is completely fixed and N 2 (τ, x) = a 2 (τ ) [1 + 2φ(τ, x)] and w ij (τ, x) = a 2 (τ ) [1 − 2ψ(τ, x)]δ ij . In the limit set by the two preceding expressions we have that and coincide therefore with the expressions of Eqs. (3.1) and (3.2). In the linearized approximation the Eqs.
where F k contains at least one spatial gradient 19 . Probably the first nonlinear generalization of inflationary curvature perturbations has been proposed in [22] after the pioneering analyses on the gauge-invariant treatment of linearised cosmological perturbations [9]. Similar variables have been subsequently scrutinized and rediscovered by different authors [23].
Although the nonlinear evolution of the curvature perturbations caused by the inhomogeneities in the viscous coefficients can be followed either in ζ i or in R i , it appears to be more useful in the latter than in the former since R is directly related to the normal mode of the system. In the presence of viscous stresses the nonlinear evolution of R i is given by: (4.4) where the ellipses stand for terms which contain, at least, three spatial gradients and are therefore of higher order in the gradient expansion. Equation (4.4) generalizes the results of Refs. [22,23] (see also [7] second paper). For the sake of conciseness the ellipses shall be neglected altogether in the subsequent discussions but it is understood that the forthcoming results hold to lowest order in the gradient expansion. Recalling the expression of the effective pressure the right hand side of Eq. (4.4) can be made more explicit; the result is: where the terms have been grouped in such a way that each of the two square brackets reproduces, respectively, the nonadiabatic and the viscous contributions in the perturbative limit. To further simplify the right hand side of Eq. (4.5) we can use first the evolution of ρ (i.e. ∂ τ ρ = K N (ρ + P)) and then rearrange the various terms. The result of this step is given by: (4.6) The first at the right hand side of Eq. (4.6) is the nonadiabatic pressure fluctuation; the second term at the right hand side of Eq. (4.6) vanishes since its contribution is of higher order in the gradients. More specifically this term can be rewritten as: but the term at the right hand side vanishes. Indeed the inhomogeneous Eintsein equations are 2ℓ 2 P ρ = K 2 − TrK 2 , 3Nℓ 2 P (ρ + P) = 2∂ τ K − 3NTrK 2 + NK 2 . (4.8) In the first equation TrK 2 = K 2 /9 + TrK 2 and TrK 2 is of higher order 20 . In summary, thanks to the results of Eq. (4.8) and in the absence of nonadiabatic pressure fluctuations, Eq. (4.6) becomes: (4.9) If the source term in Eq. (4.9) vanishes the curvature inhomogeneities will be conserved and the equations of motion will enjoy a further symmetry 21 since R i (τ, x) can be shifted by a a term constant in time (but not in space). When ξ = ξ(ρ) Eq. (4.9) implies that ∂ τ R i = 0: in this case the two terms at the right hand side simplify because ∂ i ξ = (∂ξ/∂ρ)∂ i ρ and ∂ τ ξ = (∂ξ/∂ρ)∂ τ ρ.

Bulk viscosity versus shear viscosity
To obtain Eqs. (4.7), (4.8) and (4.9) two results have been used, namely that TrK 2 is of higher order in the gradient expansion and that the shear viscosity does not contribute to the evolution to leading order in the gradient expansion. These two results also imply that bulk viscosity does contribute to the deceleration parameter while the contribution of the shear viscosity is of higher order. Given their relevance for the nonlinear discussion, these two points will now be discussed in some detail. For the purposes of this discussion we shall write explicitly the Einstein equations in their contracted form 22 where Π j i denotes the anisotropic stress (which is by definition a traceless rank-two tensor in three dimensions) and r ij are the components of the Ricci tensor of the spatial slices 23 . Both Π j i and r j i are of higher order in the gradients and we shall see that the traceless part of the extrinsic curvature is also of higher order in the gradients and it is the only component affected by the presence of shear viscosity. Indeed, after taking the the traceless part of Eq. (4.12) the following equation is obtained: (4.13) where r j i = (r j i − δ j i r/3) is traceless. Equation (4.13) shows that the shear viscosity (unlike bulk viscosity) completely decouples from the trace of the extrinsic curvature and only affects the traceless part. When η(τ, x) Eq. (4.13) can be easily solved and the result is 24 14) where τ * = τ * ( x) denotes some arbitrary integration time while, in the last line, τ 1 and τ 2 denote two generic times. Equation (4.14) demonstrates that the traceless part of the extrinsic curvature is determined by the anisotropic stress and by the traceless part of the 22 In other words we shall give, respectively, the (00), (0i) and (ij) components of the equations written in the form R ν µ = ℓ 2 3) Γ n im where the Christoffel connections are computed in terms of γ ij . 24 intrinsic curvature. Both quantities are of higher order in the gradient expansion. The last term at the right hand side of Eq. (4.14) shows that the shear viscosity suppresses the traceless part of the extrinsic curvature even further in comparison with the case η → 0.
Since the evolution of η decouples from the trace of the extrinsic curvature, it does not contribute to the inhomogeneous generalization of the deceleration parameter: the only accelerated solutions obtainable in the case of irreversible fluids are determined by the bulk viscosity coefficient. For the sake of comparison with the fully homogenous case we choose Gaussian normal coordinates and set N = 1; in this situation Eq. (4.10) can be written as: (4.15) where q( x, t) = −1 +K/TrK 2 is the inhomogeneous generalization of the deceleration parameter 25 . To discuss the sign of q(t, x) the following three remarks are in order: i) since γ ij is always positive semi-definite we have that u 0 u 0 = 1 + u 2 ≥ 1 (where, we remind, u 2 = γ ij u i u j ); ii) with the preceding observation it follows from the first two terms at the right hand side of Eq. (4.15) that q(t, x) is always positive semi-definite as long as (ρ + 3P) ≥ 0 up to correction O(u 2 ); iii) the terms multiplying the shear viscosity η all contain at least to gradients since each velocity field carries at least one gradient because of the momentum constraint of Eq. (4.11).
Thanks to the three previous observations the sign of the generalized deceleration parameter only depends on P = p − ξ∇ α u α ≃ p + Kξ (and hence on the bulk viscosity coefficient) while the shear viscosity does not play any rôle. According to Eq. (4.15) bulk viscosity only enters to second order in the gradient expansion where, however, the bulk viscosity also contributes through the term (ρ + P)u 2 implicitly contained in (ρ + P)u 0 u 0 . In conclusion the proof presented in the first part of this section is now complete and curvature perturbations are nonlinearly conserved if and when the bulk viscosity coefficient only depend on the energy density (or on the trace of the extrinsic curvature).

Tensor to scalar ratio
After the scrutiny of the nature of the quasiadiabatic modes both at the linear and at the nonlinear level, it seems useful to draw an explicit comparison between the quasiadiabatic solution and the genuine adiabatic paradigm. In this respect the simplest and most revealing quantity to estimate is the amplitude of the tensor to scalar ratio, i.e. the ratio between the scalar and tensor power spectra. As in section 3 we shall assume a slow-roll phase is only supported by the bulk viscosity (see, in particular, Eq. (3.26) and discussions therein). This case will then be compared with the conventional situation of a single field inflationary model. Even if the discussion could be conducted in fairly general terms thanks to the results of the previous sections, the attention will be focussed on the case where the slow-roll parameters are approximately constant in time.
According to Eq. (2.13) the tensor modes only couple to the curvature and their evolution equations are always the same both in the adiabatic and in the quasiadiabatic case even if the basic fields driving the slow-roll dynamics change completely from one case to the other. Equation (2.13) can be solved in the two relevant regimes namely for kτ ≫ 1 (when the relevant wavelengths are shorter than the Hubble radius) and for kτ ≪ 1 (when the wavelengths are larger than the Hubble radius). After counting properly the tensor polarizations, the power spectrum becomes, in terms of the rescaled variable µ k (τ ) The solution for µ k in the two regimes kτ ≫ 1 and kτ ≪ 1 can be written as: The values of A k and B k appearing in Eq. (5.2) are determined by demanding the continuity in τ (T ) ex of the solution and of its first derivative. We remind that, by definition, kτ (T ) ex ≃ 1. The full result, valid for τ ≥ τ ex can also be written as where a ex = a(τ (T ) ex ) and analog notation is used for H ex . In the adiabatic case the large-scale curvature perturbations are caused by a single scalar field. As already mentioned in section 3 the corresponding evolution equations can be recovered from our results by trivially setting ξ → 0, δξ → 0, δp nad → 0 and c s → 1 (see e.g. Eq. (3.9)). Ultimately the well known evolution equation of the curvature perturbations will be given, in Fourier space, by (5.4) In this case the power spectrum is given by where q (ϕ) denotes the normal mode in the scalar field case while q (without superscript) denotes the normal mode in the quasiadiabatic case (see, in particular, Eq. (3.24)). The same procedure described in the case of the tensors and leading to Eq. (5.3) can be applied in the case of q (ϕ) . The full solution analog to Eq. (5.3) but valid in the scalar case is given by: . (5.6) Note that, in the case of Eq. (5.6) the solutions valid for kτ ≪ 1 and kτ ≫ 1 have been matched for kτ (S) ex ≃ 1 where τ (S) ex coincides, in the genuine adiabatic case, with τ (T ) ex . From these expressions it is therefore possible to compute the tensor to scalar ratio denoted by r (ϕ) T which is given in the single scalar field case by: where z ex = z ϕ [τ (S) ex ]. As expected the tensor spectral index is n T = −2ǫ to lowest order in the slow-roll approximation 26 . When the fluctuations are induced by bulk viscosity the tensor contribution is exactly the same since what matters is the evolution of the extrinsic curvature which is only sensitive to theḢ/H 2 , at least during the slow-roll phase. The evolution of the scalar modes has been already discussed in detail in the last part of section 3 (see, more specifically, Eq. (3.24) and discussion therein). Therefore from Eqs. (3.29)-(3.30) and from Eqs. (3.31)-(3.33) we can obtain the tensor to scalar ratio in the quasiadiabatic case: . (5.8) Naively from Eq. (5.8) we could argue that r T gets much larger than 1 in the limit ǫ → 0: this means that, in the quasi-adiabatic solution, that the scalar modes are suppressed in comparison with the tensor modes. Even this qualitative argument is grossly correct the situation is a bit more subtle since, as indicated, the conditions for the horizon crossing are different for the scalar and tensor modes. In particular, recalling Eqs. (3.29)-(3.30) and Eqs. (3.31)-(3.33), during the slow-roll phase we have that According to Eq. (5.9) the evolution of the tensor modes implies that the horizon crossing occurs for kτ (T ) ex = 1 while for scalars it occurs for kτ (S) ex = 2 √ ǫ, i.e. |τ (S) ex | ≪ |τ (T ) ex | since, during slow-roll, ǫ ≪ 1. With these specifications we have that Eq. (5.8) can also be written as (5.10) If we compare Eqs. (5.7) and (5.10) we can observe that r T is at least O(100) times larger than r (ϕ) T . This estimate follows if we consider that, at most, |c ef f | 2 < O(1). Similarly we can take w = O(1) (but smaller than 1). Note that this estimate does not rely on a particular inflationary solution driven by bulk viscosity but just on the assumption that slow-roll dynamics is compatible with a background viscosity.
In summary, if the quasi-de Sitter phase is caused by the evolution of the bulk viscosity the quasiadiabatic scalar mode is subleading in comparison with the tensor mode. Conversely the bona fide adiabatic dominates against the tensor mode in the single field case. This result means that when the quasi-de Sitter phase is driven by bulk viscosity we are getting closer to the situation of the exact de Sitter space-time where the scalar fluctuations of the geometry should strictly vanish in comparison with the tensor mode.

Concluding remarks
Depending on the specific dynamical situation the large-scale inhomogeneities of the viscous coefficients can be classified as entropic (i.e. nonadiabatic) or quasiadiabatic. Whenever the bulk viscosity does not have a homogeneous background the resulting fluctuations are automatically gauge-invariant and their contribution to the evolution equations of the curvature perturbations reminds of the familiar source terms arising in connection with the (four) conventional entropic modes customarily constrained by means of the temperature and the polarization anisotropies of the CMB. Along this first perspective, the viscous modes are only tolerable as a subleading component of a dominant adiabatic solution whose presence is instead mandatory in the light of current large-scale observations. A second complementary possibility stipulates that the viscous coefficients have a spatial variation but in the presence of a homogeneous background. In such a situation the curvature perturbations inherit a source term depending on the fluctuations of the viscous coefficients, on the background viscosities and on the inhomogeneity of the expansion rate. The fluctuations of the bulk viscosity coefficient always act as a supplementary nonadiabatic pressure fluctuation but, this time, a quasi-de Sitter stage of expansion can be driven solely by the viscous coefficients. The present analysis demonstrates that the evolution of curvature perturbations is in general nonadiabatic. Moreover, if and when the bulk viscosity coefficient leads to quasi-de Sitter solutions and to large-scale curvature perturbations, the shear viscosity coefficient at large scales only couples to the evolution of the traceless part of the extrinsic curvature of the spatial slices and does not contribute to the accelerated expansion.
The potentially dangerous nonadiabatic source terms summarized in the previous paragraph disappear whenever the viscosity coefficients depend solely on the energy density of the relativistic fluid. In this case the curvature perturbations are effectively quasiadiabatic since they coincide with the standard adiabatic solution but only in the large-scale limit. In perturbation theory this conclusion follows from the analysis of the gauge-invariant fluctuations of the spatial curvature. The same result can also be obtained from a fully nonlinear analysis where the evolution of the curvature perturbations is studied within the gradient expansion appropriately extended to handle the viscous situation. Unfortunately the curvature power spectrum of the quasiadiabatic solution is parametrically smaller than the tensor power spectrum. Hence the corresponding tensor to scalar ratio turns out to be larger than in the standard adiabatic case where the scalar power spectrum dominates, over large scales, against its tensor counterpart. Taken at face value the obtained results show that the viscous coefficients alone cannot drive a phase of accelerated expansion and, at the same time, reproduce the standard adiabatic scalar mode. We can get very close to an acceptable phenomenological situation if the bulk viscosity coefficient only depends on the energy density of the plasma. Even in this case, however, there are serious drawbacks since the dominance of the tensors against the scalars is at odds with a pretty robust observational evidence. In this respect the obtained results suggest a novel strategy for a concrete phenomenological scrutiny of the large scale inhomogeneities induced by the viscous coefficients.
In a more optimistic perspective the large-scale fluctuations of the viscous coefficients remain a viable possibility only when the dominant adiabatic solution comes from a different physical origin. In this case the fluctuations of the bulk viscosity play the same rôle of a supplementary nonadiabatic solution in the space of the initial conditions of the Einstein-Boltzmann hierarchy. Such a component can be constrained prior to photon decoupling and across the matter-radiation transition with the same techniques customarily employed to bound the presence of the standard four nonadiabatic solutions.