Scalar and vector perturbations in a universe with nonlinear perfect fluid

We study a three-component universe filled with dust-like matter in the form of discrete inhomogeneities (e.g., galaxies) and perfect fluids characterized by linear and nonlinear equations of state. Within the cosmic screening approach, we develop the theory of scalar and vector perturbations. None of the energy density contrasts associated with the distinct components is treated as small. Consequently, the derived equations are valid at both sub- and super-horizon scales and enable simulations for a variety of cosmological models.


Introduction
Unraveling the physics behind the accelerated expansion of the late universe is one of the greatest challenges in modern cosmology. Going through a number of models that have been proposed to explain this phenomenon, one may first refer to Einstein's legendary cosmological constant as the underlying energy component. However, though it provides the sought-for behavior, it also brings along the question of what could be the origin of such a constant with the value required for the observed acceleration, yet unanswered. Perfect fluids with the linear equation of state (EoS) p = ωε, ω = const are also known to drive the acceleration provided that ω < −1/3. Nevertheless, ω should be very close to −1 [1] to be in agreement with observations. Alternatively, there exist the so-called quintessence [2,3] models in which the corresponding source is a scalar field. These include cases leading to constant ω [4], shown to impose severe restrictions on the form of the scalar field potential [5,6]. Perfect fluids with the nonlinear EoS p = f (ε) represent yet another class of models in this direction, and the Chaplygin gas model [7][8][9][10] may be studied as an example of those.
As one can see, there exists an extensive zoo of models that can explain the accelerated expansion of the universe. Obviously, the models to survive are those with predictions closest to the observational data. Among such criteria is the observed large-scale structure of the universe described on the basis of the widely studied perturbation theory [11][12][13][14][15][16]. In this direction, a cosmic screening approach has been proposed in the paper [17] for the scalar and vector perturbations, with the distinctive feature that the gravitational potential satisfies a Helmholtz-type equation and not a Poisson-type one. Consequently, at large cosmological distances from individual sources, the potential undergoes an exponential cutoff. Matter sources in [17] have the form of discrete point-like masses, and it is also important to emphasize that no assumptions are made regarding the smallness of the associated energy density contrast, ensuring the validity of the model both at sub-and super-horizon scales. This approach has been further developed in the papers [18][19][20][21][22][23][24][25][26][27]. Particularly in [18,19], generalizations to the case of perfect fluids with the linear p = ωε and nonlinear p = f (ε) equations of state have been performed, i.e. aside from point-like masses, additional components of two distinct types of perfect fluids have been included in the matter sources. Obviously, the linear component has the background EoS p = ωε. Nevertheless, for the nonlinear component, it is possible to write p = f (ε) only in the case of small fluctuations where the expansion p = f (ε) + (∂f /∂ε) ε δε + (1/2)(∂ 2 f /∂ε 2 ) ε δε 2 + . . . works well. This point is disregarded in some papers (see e.g. [28][29][30]) whereas in [31][32][33][34][35], for instance, it is clearly stated that the relation p = f (ε) does not hold in general.
In [18], the authors considered the case where energy density fluctuations of the nonlinear perfect fluid are small quantities. Meanwhile, density contrasts of the pressureless matter and perfect fluid with linear EoS were set arbitrary. In the present article, we investigate all three matter components with arbitrary density contrasts. Therefore, the considered model applies both to small/astrophysical scales, where matter fluctuations are large, and to large/cosmological distances, where the density contrast is small. We develop the theory of scalar and vector perturbations for this model within the cosmic screening approach and obtain a system of equations which enables the cosmological simulation for arbitrary forms of the function f (ε).
The paper is structured as follows. In Section 2, the basic equations are presented and the theory of scalar and vector perturbations is constructed for the considered model. The main results are summarized in Section 3. Appendix is devoted to showing that the auxiliary equations employed in the main proof are satisfied within the adopted accuracy for arbitrary density contrasts.

Scalar and vector perturbations in the cosmic screening approach
We investigate a universe which contains perfect fluids with nonlinear EoS p J = f J (ε J ), J = 1, 2, . . . Particular cases include the pressureless perfect fluid with p = 0 and perfect fluids with linear EoS p = ωε, ω = const. The space-averaged distribution (denoted by the overbar) of these components determines the dynamics of the homogeneous and isotropic universe described by the Friedmann equation in which a (η) denotes the scale factor, the Hubble parameter H ≡ a ′ /a ≡ (da/dη)/a, η is the conformal time and the constant κ ≡ 8πG N /c 4 (G N and c are the Newtonian gravitational constant and the speed of light, respectively). The total averaged energy density ε in the above equation has been split into its constituent parts with respect to their types of EoS: index "M" corresponds to pressureless matter (continuous as well as discrete) and indexes "I" and "J" correspond to perfect fluids with linear and nonlinear EoS, respectively. Fluctuations in the energy densities generate metric perturbations which, in the following, will be studied in terms of their associated scalar and vector components. The perturbed metric in the first-order approximation and in the conformal-Newtonian gauge reads [11,13,15,16] The only approximation in our approach is that the metric corrections Φ and B α as well as the peculiar velocitiesṽ α =ṽ α ≡ dx α /dη are considered small: Φ, B α ,ṽ α ≪ 1. On the other hand, the smallness of the energy density and pressure fluctuations is not demanded, i.e. the density and pressure contrasts may exceed unity: δε/ε, δp/p > 1. This serves as an indicator that our model works both on astrophysical and cosmological scales. Potentials Φ and B satisfy the linearized Einstein equations [18] △Φ − 3H(Φ ′ + HΦ) = 1 and its averaged energy density ε M = ρ M c 2 /a 3 . For such a component, the energy density fluctuation reads [36,37] with δρ M ≡ ρ M − ρ M . As for the linear perfect fluid, the energy density can be considered in the form [18][19][20] where A I ≡ A I + δA I and A I = const. Since each matter component separately satisfies the energy conservation equation (see Eq. (A.20) in [18]) in which ε represents any of the individual components and Φ and B are the total potentials produced by the combination of components, one can easily show that the function A I fulfills For the averaged quantities, Eq. (2.8) yields and evidently, ε I = A I /a 3(1+ω I ) satisfies this equation.
Let us now turn to the nonlinear perfect fluid with EoS p J = f J (ε J ), where f J represents some nonlinear function. Since the fluctuations in the energy density and pressure are not restricted to small values, it is no longer possible to substitute p J = f (ε J ) for the background pressure; we need to proceed in a rather different way. Similar to Eqs. (2.6) and (2.7), we consider the energy density in the form where F J is an unknown function for which we will derive an equation subsequently. Provided that |Φ| ≪ 1, expanding the quantities ε J and p J accordingly, we may write It is naturally demanded that the energy density (2.11) satisfies the conservation equation (2.8). In this connection, we substitute (2.12) and (2.13) into (2.8), which yields (2.14) In obtaining the above expression we have neglected the terms quadratic with respect to Φ and used the relation ∇B = 0 as well as f ′ (F ) = ∂f /∂ε| ε=F F ′ . Now, we decompose the functions F J and f J (F J ) into average values and fluctuations as where both quantities F J and f J (F J ) depend only on time. Obviously, ε J = F J and p J = f J (F J ) and hence the background equation (2.10) for the "J"-component reads Substituting the decomposed functions (2.15) into (2.14) and taking into account (2.16), we get Before proceeding further, we pause to make a few important comments. First, on large/cosmological scales the quantities Φ, B,ṽ are of the same order of smallness ǫ, i.e. Φ ∼ B ∼ṽ ∼ ǫ ≪ 1. Meanwhile at small/astrophysical distances, we have Φ ∼ ǫ and B ∼ṽΦ. Second, the perfect fluid is considered to behave the "normal" way. By this we mean that the squared speed of sound c 2 s = δp/δε ∼ ∂f J /∂ε J 1 and, additionally, the ratio of pressure fluctuations over pressure is of the order of its energy density counterpart: Keeping these in mind, the terms multiplied by the scalar perturbation Φ can be safely neglected in comparison to the rest of the expression in (2.17). For δf J B ≪ f JṽJ , we find that the left-hand side (LHS) is further simplified to yield where we have taken into account that the terms proportional to the products Φṽ α I,J are to be neglected in the first order. Moving further, we decompose the terms withṽ α n,I,J on the RHS into their longitudinal and transverse parts as (see Eqs. (2.24) and (2.25) in [18]) (2.23) Here the function Ξ has the form [17] whereas ξ I and ζ J are to be determined numerically. Using this system of equations, we split Eq. (2.20) into scalar and vector parts so as to obtain the equations we get the equation for the potential Φ: It is possible to reformulate (2.26) and (2.28) so that we have With only a single "J"-component 2 in the universe, Eqs. (2.29) and (2.30) are reduced to respectively, accompanied by the evident redefinition of λ.
We also need the momentum conservation equation for the "J"-component. Using (A.23) from [18] and eliminating the terms proportional to Φ 2 , ΦB, Φṽ, Bṽ, which are of the higher order of smallness, we obtain This equation may readily be used to define the time derivative of peculiar velocity and to this end, one may replace B ′ with the help of the equation B ′ = −2HB [17,18].
To conclude this section, we briefly consider the Chaplygin gas model [7][8][9] as a particular example of the "J"-type nonlinear perfect fluid. The EoS of the modified generalized Chaplygin gas has the form [10] for which the background quantities read These satisfy the conservation equation (2.16)

Conclusion
In our paper we have studied a universe consisting of three types of components. It contains dust-like matter (denoted by the index "M ") in the form of discrete inhomogeneities (e.g., galaxies, galaxy groups and clusters), perfect fluids with linear EoS p I = ω I ε I (ω I = const) and perfect fluids with arbitrary nonlinear EoS p J = f J (ε J ). The background spacetime geometry is defined by the Friedmann-Lemaître-Robertson-Walker metric. All three components have been considered to have arbitrary energy density contrasts. Within the cosmic screening approach, we have developed the theory of scalar and vector perturbations and obtained a system of equations that enables the numerical simulation of models with an arbitrary form of the function f (ε). Since we have not assumed the smallness of energy density fluctuations, these equations are valid on both small/astrophysical and large/cosmological scales. Additionally, we have checked some of the important auxiliary equations to demonstrate that they are indeed satisfied (up to the adopted accuracy) for arbitrary density contrasts.

A Check of the equation (2.25)
The cosmic screening approach [17,18] applies to all scales down to astrophysical distances where energy density fluctuations of the perfect fluids are no longer small quantities. In this appendix we intend to demonstrate that Eq. (2.25) holds for arbitrary values of density contrasts. For this purpose, we will be preserving the terms of the type δρΦ, δAΦ, δF Φ which can exceed ρΦ, AΦ, F Φ at small scales.
We have already noted that "M "-and "I"-components can be considered as particular cases of the "J"-component. Hence, it would indeed be sufficient to limit ourselves to the "J"-components in (2.25). Nevertheless, our reasoning will be clearer if we first prove the equation for the "M "-component, then for the "I"-component, and only after that for the "J"-component. Considering the "M"-component separately from the other two is especially justified by the fact that there exist exact expressions for the potentials Φ and B in the corresponding model.

A.1 Pressureless matter
We consider the "M "-component in the form of discrete point-like masses with the comoving mass density (2.5). For this case, the potential Φ in Fourier space reads [17] Above and hereafter the hat denotes the Fourier-transforms.
In the current arrangement, Eq. (2.25) takes the form where [17] Differentiating the expression (A.1) with respect to conformal time and neglecting the terms quadratic inṽ n , we obtain Now, taking into account the relation [17] The term −3i H k 2 k [δρ∇Φ] ∼ Hδρ·ǫ in the square brackets can be neglected in comparison with the term −H δρ ∼ Hδρ (their ratio gives precisely the order of smallness ǫ). Similarly, the term 3i H ak 2 k [δρ(aB) ′ ] is much less than the term −3i H 2 a n ρ n (kaṽ n ) k 2 . Their ratio is of the order of δρB/ (ρṽ) ∼ Bṽ/ ρṽ 2 /δρ ∼ Bṽ/Φ ∼ ǫ. Here, we use the relation B ′ = −2HB, the estimate (2.18) and take into account that at small/astrophysical scales B ∼ṽΦ while at large/cosmological scales B ∼ Φ (see also footnote 1).
Substituting the expressions for Φ and Φ ′ into the left-hand side (LHS) of (A.3), we get and this exactly coincides with the RHS of (A.3) (see Eq. (A.4)).

A.2 Perfect fluid with linear equation of state
In this section we consider perfect fluids with linear EoS p I = ω I ε I , ω I = const. From Eq. (2.25) we obtain where ζ I ≡ ξ I a 3(1+ω I ) . (A.10) Our intention now is to prove that Eq. (A.9) holds for arbitrary values of energy density fluctuations. From Eq. (2.29) and taking into account the definition (2.31) for λ, we get The time derivative of this equation is and substitution of (A.11) and (A.12) into the LHS of (A.9) results in where the terms ∝ κa 2 HδA I /a 3(1+ω I ) have been cancelled.
In order to proceed, we need to determine the quantities δA ′ I and ζ ′ I . For this purpose, we employ (2.7) as the energy density ε I = ε I + δε I of the "I"-component and re-express Eq. (2.8) in terms of this component only, which reads (A.14) Here we have taken into account the definition (2.22) as well as the transverse gauge condition ∇B = 0. It is worth noting that neither of the cases |ω I | ≫ 1 or |1 + ω I | ≪ 1 are being considered in the current configuration. Now, the second term in (A.14) can be neglected in comparison with the first one, since δA Further on, in order to find ζ ′ I , we use the momentum conservation equation (A.23) in [18] applied to the "I"-component (which is the analogue of (2.34) for the "J"-component): Using the definition (2.7), eliminating the second order terms, replacing B ′ with the expression −2HB and then acting by ∇, we obtain Since △ζ I = ∇ A IṽI /a 3(1+ω I ) , all terms in the first line can be neglected in comparison with the first and second terms in the second line (the ratios reduce to ∼ (δA I /A I )(B/ṽ I ) ∼ Bṽ I /Φ). Hence, and we arrive at Here, the third line represents a sum of summands of the order ∼ δA I Φ (with factor ∼ (1 + ω I )), and each summand of this order can be ignored in view of the fact that the terms of the superior order ∼ δA I (with factors ∼ 1 (see Eq. (A.13)) and ∼ ω I ) have already been cancelled in the previous steps. Our reasoning is that since we keep the terms in our equations up to only a certain order of smallness, previous cancellations of terms do not correspond to identical zero. To respect the adopted accuracy, it is necessary to discard also the terms with higher orders of smallness than these. Therefore, from (A.21) we obtain and this is an identity. Hence, Eq. (A.9) proves adequate within the considered accuracy. Now, to compare some of terms above with the previously cancelled ones, we should remember that, first, all terms in (A.31) were divided by κa 2 and second, the expression was acted on with △. This means that either we should apply △ to the cancelled terms or, vice versa, apply △ −1 to the remaining ones. For example, the cancelled term ∝ κa 2 H 2 where we have used (2.23). Then, employing the conditions δF J ∼ δf J , |∂f J /∂ε J | ∼ 1 and taking into account that time derivatives with respect to η are proportional in order to the Hubble parameter ∼ H, it can be deduced that all three terms in the second line behave as H 2 J ∇(BδF J ). Thus, the ratio of these terms and the cancelled term is of the order of (δF J /F J )(B/ṽ J ) ∼ Bṽ J /Φ ∼ ǫ (see the estimate (2.18)). Similarly, the very last term ∝ HδF J Φ may also be neglected provided that terms of higher order, i.e. ∝ κa 2 H J (δF J + δf J ) have already been cancelled in the above steps. Eventually, (A.31) is reduced to and taking into consideration that ∂f J ∂ε J ε J =F J is of the order of unity, we are left with the term ∝ ΦH (δF J + δf J (F J )), which is again to be neglected in comparison to previously cancelled terms in the expression. Therefore, the LHS of this equation is equal to zero up to adopted accuracy, and this serves as the proof of the equation (A.24).
To conclude this appendix, it is worth noting that in a similar way we can also prove the relation B ′ + 2HB = 0.