Bulk viscous Zel’dovich fluid model and its asymptotic behavior

In this paper we consider a flat FLRW universe with bulk viscous Zel’dovich fluid as the cosmic component. Considering the bulk viscosity as characterized by a constant bulk viscous coefficient, we analyze the evolution of the Hubble parameter. Type Ia Supernovae data is used for constraining the model and for extracting the constant bulk viscous parameter and present the Hubble parameter. We also present the analysis of the scale factor, equation of state, and deceleration parameter. The model predicts the later time acceleration and is also compatible with the age of the universe as given by the oldest globular clusters. Study of the phase-space behavior of the model shows that a universe dominated by bulk viscous Zel’dovich fluid is stable. But the inclusion of a radiation component in addition to the Zel’dovich fluid makes the model unstable. Hence, even though the bulk viscous Zel’dovich fluid dominated universe is a feasible one, the model as such fails to predict a prior radiation dominated phase.


Introduction
Observational data on Type-Ia supernovae [1][2][3] and the CMB [4,5] has confirmed with sufficient accuracy that nearly 70 percent of energy of the universe is in an exotic form called dark energy, which is responsible for the current acceleration of the universe. The remaining part of the cosmic components consist of nearly 23 to 24 percent of weakly interacting dark matter [6][7][8][9][10] and a few percent forms the luminous matter and radiation. Even with the overwhelming evidence for the existence of these cosmic components the current observational data does not rule out the possible existence of other form of exotic fluid components. One of the examples for such fluids is the dark radiation which can exist in the early or later stage during the evolution of the universe [11]. a e-mail: rajagopal-physics-61@cusat.ac.in b e-mail: titus@cusat.ac.in Recently, considerable attention has been paid to the study of another exotic fluid, the Zel'dovich fluid or stiff fluid, first studied by Zel'dovich [12]. The Zel'dovich fluid is a perfect fluid in which the speed of sound is equal to the speed of light, so that the equation of state becomes ω z = p z /ρ z = 1, the highest value a fluid can have in consistency with causality. The Zel'dovich fluid or stiff fluid behavior in the cosmological context has been considered by many workers. In dealing with the self-interaction in dark matter the authors in Ref. [13] have shown that the self-interaction field behaves like a stiff fluid. The existence of a Zel'dovich fluid was confirmed in the Horava-Lifshitz gravity-based cosmological models when the so-called detailed balancing conditions [14,15] were relaxed [16,17], and the stiff fluid in such a situation had been studied in Refs. [18][19][20]. The relevance of the existence of the Zel'dovich fluid in the early universe was discussed in Ref. [21]. In certain inhomogeneous cosmological models a stiff fluid appears as an exact nonsingular solution [22,23]. In the standard evolution of the Friedmann universe the density of the Zel'dovich fluid is found to be decreasing faster than radiation and matter. Consequently its effect on the early universe would be significant. One of the important phenomena that took place in the early universe is the primordial nucleosynthesis which might be influenced by the presence of a stiff fluid. In Ref. [24] the authors have found a limit on the density of the stiff fluid from the constraints on the abundances of the light elements.
In an expanding universe there can be deviations from local thermodynamic equilibrium. Consequently there can arise a bulk viscosity in the cosmic fluid which will restore the equilibrium [25]. This bulk viscosity modifies the effective pressure of the fluid in order to facilitate reestablishment of the equilibrium situation. As soon as the equilibrium is reached the bulk viscous pressure vanishes [25,26]. In the context of inflation in the early universe it has been shown that an imperfect fluid with bulk viscosity can cause the early accelerated expansion [27]. The later time viscous universe was studied in Ref. [28]. A considerable number of studies by including bulk viscosity in the dark matter setting were carried out by many authors in the context of the late acceleration of the universe [29][30][31].
Recently considerable interest has been focused on the study of a viscous Zel'dovich fluid in an expanding universe. In Ref. [32] the authors studied the evolution of a viscous Zel'dovich fluid in a flat universe and found that it can have a considerable effect even in the late universe. This work shows that the bulk viscous Zel'dovich fluid can cause the recent acceleration of the universe. To account for the recent acceleration of the universe, various dark energy models with non-exotic components have been studied in the recent literature. One of the actively studied proposals in this regard is bulk viscous dark matter [29,31]. But there are arguments that this model fails in predicting many properties, including the age of the universe [31]. It is in this context, when there is no consensus on any particular dark energy model, that the authors of Ref. [32] considered the bulk viscous Zel'dovich fluid as a possible candidate for causing the recent acceleration of the universe. In Ref. [33] it was shown that, with a bulk viscous Zel'dovich fluid as cosmic component, the age of the universe also can be predicted reasonably well.
In Ref. [32], it was shown that the non-viscous Zel'dovich fluid density will evolve as ρ z ∼ a −6 , in a flat universe with expansion scale factor a. In such a case the Zel'dovich fluid has a notable effect only on the early evolution of universe compared to ordinary dark matter, whose density evolves as a −3 . But the evolution of the bulk viscous Zel'dovich fluid is different from the non-viscous one and it can cause the late acceleration of the universe. The authors of [32] considered this model and proved that a bulk viscous Zel'dovich fluid does produce late acceleration of the universe. However, they did not try to constrain the model with cosmological observational data to arrive at a realistic picture. The question whether the model would successfully allow for a prior radiation dominated phase is also to be considered. In the present work we analyze the evolution of a flat universe with a bulk viscous Zel'dovich fluid and compare the predictions of the model with the latest cosmological data on type Ia supernovae. We evaluate the model parameters including the transport coefficient of bulk viscosity and study the evolution, particularly in the later stage, of the universe.
The paper is organized as follows. In Sect. 2, we consider a flat universe with a bulk viscous Zel'dovich fluid as the only component and derive the Hubble parameter. Then we constrain the model using type Ia supernova data to extract the constant bulk viscous parameter and the present value of the Hubble parameter. We also include the evolution of the equation of state parameter and deceleration in this section. In Sect. 3, we present our analysis of the phase-space structure of the model by including radiation also as a cosmic component to see whether the model naturally allows a prior radiation dominated phase. This is followed by our conclusions in Sect. 4.

The bulk viscous Zel'dovich fluid model
In this section we consider a flat universe with a bulk viscous Zel'dovich fluid as the only component. The main feature of the Zel'dovich fluid is that sound velocity in the fluid is equal to that of light. The equation of state [12] is given by A similar equation of state was studied with reference to some special case by Masso and Rota [34]. In including the viscosity of the fluid we will follow Eckart's formulation [35], which deals with viscous dissipative processes occurring in a thermodynamical system when it deviates from local equilibrium. An equivalent formulation was developed by Landau and Lifshitz [37]. However, it was noted that the equilibria in Eckart's frame are unstable [38] and signals were propagated through the fluid at superluminal velocities [39]. These drawbacks were rectified in a more general formalism by Israel et al. [40,41], from which Eckart's theory follows as a first order limit. But many authors are still using Eckart's theory because of its simplicity. For example, Eckart's formalism was used in some models on the late acceleration of the universe caused by the bulk viscous dark matter [42][43][44][45][46][47]. In the mean time Hiscock and Salomonson [48] have shown that Eckart's formalism can lead to accelerated expansion in the FLRW model. However, later studies have shown that an inflationary solution is possible with the Israel-Stewart formalism also [27]. In the Eckart formalism, the first order deviations from equilibrium can be expressed as an additional non-adiabatic contribution T μν to the energy momentum tensor [35], where with ζ the coefficient of viscosity (positive due to the second law of thermodynamics [36]), u μ is the four velocity of an observer who measures the energy density and pressure, g μν is the metric tensor, and u γ ;γ = 3H for an expansing universe. In the present work we too follow Eckart's approach. The above results show that the effective pressure of the bulk viscous Zel'dovich fluid can be expressed as where ζ is the coefficient of viscosity and H is the Hubble parameter.
We consider a flat Friedman universe with FLRW metric given as where t is the cosmic time, a(t) is the scale factor of expansion, r, θ, φ are the comoving coordinates. This, when combined with the Einstein field equations, gives the dynamical equations and the conservation equatioṅ Here we follow the standard convention 8π G = 1. When these equations are combined with Eq. (4) for the effective pressure, we have Solving these equations after changing the variable from t to x = log a we get where h = H/H 0 , H 0 is the present Hubble parameter and ζ = ζ /H 0 is the dimensionless viscous parameter. Forζ = 0 the Hubble parameter becomes H ∼ H 0 a −3 . The density of the Zel'dovich fluid will then evolve as ρ ∼ a −6 and the scale factor will evolve as a ∼ (H 0 t) 1/3 and hence the universe would be eternally decelerating and the effect of the Zel'dovich fluid will be relevant to the early epoch of the universe [24]. On the other hand, the presence of the viscosity will modify the density in such a way that In the early universe a term proportional to a −6 would dominate and affect the processes in the early universe like primordial nucleosynthesis [11] etc. But a term proportional to a −6 is drastically decreasing as the universe expands, while the second term, proportional to a −3 , is a matter-like term, which in fact guarantees the effect of the bulk viscous Zel'dovich fluid in the later universe. As is obvious from Eq. (4), the viscous term causes a negative pressure, which in turn causes an accelerated phase of expansion in the later stage of the evolution of the universe. In the extreme limit of a(t) → ∞ the density behaves as ρ z → H 2 0 12ζ 2 . As a result the admissible values ofζ are very important in this model, which is to be evaluated by the observational constraints. Since we are interested in the effect of a Zel'dovich fluid in the late evolution of the universe, the bulk viscous coefficient can be determined by contrasting it with supernova data.

Extraction of the model parameters using Type Ia supernovae data
The best fit values of the model parameterζ and the Hubble parameter H 0 can be extracted using type Ia supernova observational data. We have used Union data which consists of 307 data points [49] in the redshift range 0.01 < z < 1.6. The distance modulus of the supernova at a redshift z is where m is the apparent magnitude, M is the absolute magnitude, and d L (z) is the luminosity distance of a supernova in a flat FLRW universe. The expression for the luminosity distance is where h(z ,ζ ) is identical with the normalized Hubble parameter given in Eq. (9). The distance moduli of supernovae at various redshifts are calculated using Eq. (12) and are compared with the corresponding observational data. We then construct the statistical χ 2 function where μ i is the theoretical value of the distance modulus for a given redshift obtained in the present model, μ k is the observed distance modulus of the kth supernova corresponding to the same redshift as of the theoretical one, σ k is the variance of the measurement, and n is the number of data points. The best estimates of the parameters (ζ , H 0 ) are obtained by minimizing the χ 2 function. The minimum of the χ 2 gives measures of the goodness-of-fit of the model apart from giving the best estimates of the model parameters. The confidence regions in Fig. 1  H 0 = 70.20 ± 0.58. In a previous paper [32] we have worked out a single component Zel'dovich fluid model of the universe having dark energy incorporated in terms of linear bulk viscosity, and there we had shown that when 4 <ζ < 6, a transition from deceleration to acceleration took place in the past.

Evolution of cosmic parameters
The behavior of the scale factor in the Zel'dovich fluid dominated universe can be obtained from Eq. (9) as At sufficiently early time the scale factor can be approximated as a(t) ∼ (1 + 3H 0 (t − t 0 )) 1 3 , which implies a decelerated phase, while at later times the scale factor behaves as a(t) ∼ exp(ζ 6 H 0 (t − t 0 )) showing that the universe shows an accelerated evolution at a later time. Figure 2 shows the evolution of scale factors at various choices ofζ . From the figure it is seen that the behavior of the scale factor is different forζ > 6. If one finds the expression for the age of the universe it is of the form Forζ > 6 the age is not defined and consequently the universe does not have a big bang. But for the casesζ < 6 the universe does have a big bang. For the best estimates of the parameter the age of the universe is found to be around 10-12 Gy.
The equation of state of the bulk viscous fluid can be obtained using the standard relation, where x = ln a. Substituting for h in Eq. (16) using Eq. (9), we have Figure 3 shows the variation of the equation of state parameter against the redshift at various choices ofζ . In the extreme future (z → −1), ω z → −1 and it hence corresponds to a de Sitter universe. Otherwise ω z shows a strong dependence on the bulk viscous parameter. For small viscosity, the equation of state parameter remains 1 but reduces to −1 in the distant future. Forζ < 6, likeζ = 5.25, the best estimated value, ω z , is positive for z > 1. But in the subsequent evolution, it reduces to negative values and finally stabilizes at −1 as z → −1. This means that for z < 1 the bulk viscous Zel'dovich fluid mimics the quintessence nature forζ > 6,. For instance whenζ = 6.5 as in the figure, when ω z ≤ −1 always, it corresponds to a phantom nature. We have also evaluated the deceleration parameter. The basic equation of the deceleration parameter is Substituting for H = H 0 h from Eq. (9) we have As seen in Fig. 4, the deceleration parameter q remains −1 for all possibleζ in the distant future. For a small value of ζ the deceleration remains at 2 until a distant future when it drops to −1. For the best estimated value ofζ = 5.25, the switchover from deceleration to acceleration takes place at about z = 0.52, which is closely in agreement with the observational constraints. The same will proceed with ever increasing acceleration to an asymptotic value of −1 in a distant future. For values ofζ > 6.0, such as the one indicated forζ = 6.5 as in Fig. 4, q is always negative and is increasing as the universe expands, saturating to −1 at z → −1 in the future. So the final state of the universe in this model is a de Sitter universe for any positive valueζ .

Phase-space perspective
A convenient method to understand the global picture of the model is to study the equivalent phase space. A phase-space analysis of the model would indicate whether the model is compatible with a realistic evolution of the universe, i.e. the possible existence of the different stages of the universe like a radiation dominated early phase, and then a matter dominated phase followed by the late accelerating phase. But it is generally difficult to solve the cosmological field equations with more than one cosmic component. For doing the phase-space analysis, one has first to identify the phase-space variables and be able to write down the cosmological equations as a system of autonomous differential equations. The critical points of these autonomous differential equations can then be correlated to the cosmological solutions. The stability of such critical points can be determined by examining the system obtained by linearizing about the critical point. If the critical points were a global attractor, then the trajectories of the autonomous system constructed near a critical point will always be attracted toward it, independent of the initial condition, and it will be stable one.

Analysis of Zel'dovich fluid in two dimensional phase space
In this section we analyze the phase-space behavior of a flat universe dominated with a bulk viscous Zel'dovich fluid. The behavior of the system in the two dimensional phase space with h and z as the coordinates is examined. The corresponding universe is the one with only a bulk viscous Zel'dovich fluid as the only component. The coupled differential equations arė where h = H/H 0 and z = ρ z /3H 2 and the over-dot refers to a derivative with respect to time. By settingḣ = 0 anḋ z = 0, we obtain the following three critical points or roots: The first root, (0,1), has a Hubble parameter zero and hence corresponds to a static universe, while the second root depends on the instantaneous value of z , which follows a trajectory as the Zel'dovich fluid density changes, and hence it is not a fixed point. The third root (h c , zc ) = (0.87667, 1), having a positive definite value of the Hubble parameter and mass parameter of the Zel'dovich fluid, corresponds to an expanding universe dominated by Zel'dovich fluid. If the system is stable in the neighborhood of a critical point, the linear perturbation in its neighborhood in phase space decays with time. The perturbations around the critical points must satisfy the following matrix equation: Here and η are perturbations in h and z , respectively, in the neighborhood about a given critical point. The suffix 0 denotes the value evaluated at the critical point, (h c , zc ). The corresponding Jacobian is If the eigenvalues of the Jacobian matrix are all negative, then the critical point is stable, otherwise the critical point is generally unstable. If the eigenvalues are positive, then the critical point is an unstable node, and if there are both positive and negative eigenvalues, then the critical point is a saddle point.
The eigenvalues corresponding to the first critical point (h c , zc ) = (0, 1) are found to be −368.2 and 184.2. As they are of opposite signs the critical point is a saddle point and hence is unstable. Depending on the initial conditions the nearby trajectories around this point may approach the saddle point, but are repelled by it, finally approaching a possible stable attractor in the future. As in Fig. 5 the trajectories are turning away from the equilibrium point as and when they approach it and finally converge on the critical point shown on the right side of the plot.
The second critical point ( 0.87667 z , z ) is not an isolated point, but varies with z . As per the relationship between h and z , it represents a rectangular hyperbola with the axes h = 0 and z = 0 as asymptotes. The eigenvalues are found to be (−184.1, −161.394 ). Since both eigenvalues are negative and real, the neighboring trajectories will converge on to the hyperbola and hence the critical point is a stable one. The hyperbola along which the z dependent critical point moves has the coordinate axes h = 0 and z = 0 as the asymptotes, h = − z as the directrix, and (0.9363, 0.9363) as the focus.
The third critical point is (h c , zc ) = (0.87667, 1). It is observed, as in Fig. 5, that this critical point is a global attractor, and physically it corresponds to a flat expanding universe dominated by a bulk viscous Zel'dovich fluid. We  Fig. 6 shows that the field directions are invariably tilted, though slightly, toward a fixed critical point as they approach. However, on low resolution, it seems to be a straight line toward the isolated critical point. This is evident from the continuous plot in the phase-space structure as shown in Fig. 7. So, in reality the isolated critical point exists and the '0' eigenvalue leads to a line segment as the best fit close to the critical point. This is clear from the fact that the straight line does not arise from the original procedure of settingḣ = 0 and˙ z = 0 without the linear approximation and rather results in an isolated point. So, as we said earlier, depending on the initial conditions the trajectories emanating from the surroundings of the saddle critical point are repelled away from it and they finally approach the stable critical point, (0.87667, 1). From the equations of the deceleration parameter (Eq. (18)), it can easily be seen that this critical point corresponds to q < 0 and implies an accelerating phase. This, by and large, implies the stability of the universe dominated by the bulk viscous Zel'dovich fluid, which accelerates the expansion in the later stage. In the realistic case, the universe has a radiation dominated phase followed by a matter dominated phase and a subsequent late accelerated epoch. In order for the present model to imply a realistic evolution of the universe, it predicts a prior radiation dominated phase. Thus, aware of this, we include, besides the Zel'dovich fluid, the conventional radiation also to study the phase-space structure. On including the radiation component also, the first Friedmann equation becomes where ρ γ is the radiation density. The conservation equation for the radiation component, by assuming a pressure p γ = ρ γ /3, iṡ The phase-space variables are h, z , and γ among which the third parameter is γ = ρ γ /3H 2 . The dynamical equations for these parameters are represented by the coupled differential equations, The critical points are obtained by settinġ and they are out of which the first critical point is not an isolated one, since in this case h is inversely proportional to the instantaneous value of z . The second critical point, (0, 1, 0), has the Hubble parameter h = 0, hence corresponds to static universe, and the third one corresponds to an expanding universe dominated by a bulk viscous Zel'dovich fluid. It is to be noted that there is no critical point corresponding to a radiation dominated phase. The stability of the equilibrium points in the case of these three critical points is obtained (this time in the 3D phase-space case), once again by looking at the behavior of phase-space trajectories close to them and generated due to different initial conditions. The coupled differential equations in the linear limit in matrix representation, in the neighborhood of the equilibrium points, are where˙ ,η, andν are first order perturbation terms of P(h, z , γ ) =ḣ, Q(h, z , γ ) =˙ z , and R(h, z , γ ) =˙ γ , respectively, (t), η(t), and ν(t) being the first order linear perturbation terms of h, z , and γ , respectively. The square matrix term in Eq. (32) is the Jacobian evaluated at the critical point. We then decouple the differential equations (32) by means of the secular equation and in the process the eigenvalues corresponding to the equilibrium point ( 0.87667 zc , zc , 0) are obtained: and the phase-space trajectories corresponding to various initial conditions are repelled from, rather than being attracted to the second critical point. So none of the critical points in this case corresponds to a radiation dominated phase and also even the existing critical points are not stable. In fact the third critical point, which corresponds to a Zel'dovich fluid dominated one, is unstable, and it can be concluded that the inclusion of the radiation component may lead to a complete breakdown of the model. The bulk viscous coefficient is taken as a constant in the present study. Since it is a transport coefficient which may depend on the velocity of the fluid component also, such a velocity dependent bulk viscous coefficient may be checked for consistency of a prior radiation dominated phase and this problem we reserve for a future work.

Conclusion
In this paper we have considered a flat universe consisting of a bulk viscous Zel'dovich fluid. The viscosity parameter was incorporated as per Eckart's formalism. We have evaluated the evolution of the Hubble parameter. The model was constrained with SNe Ia data to evaluate the bulk viscous coefficient asζ = 5.25 ± 0.14 and the present value of Hubble parameter H 0 = 70.20 ± 0.58. The behavior of the resulting scale factor shows that the model predicts a late acceleration in the expansion of the universe. Hence the bulk viscous Zel'dovich fluid can mimic the role of the conventional dark energy.
We also studied the model to analyze the stability of the solutions corresponding to various scenarios using the phasespace analysis method. We first analyzed the two dimensional phase-space behavior, where the contribution due to radiation is neglected, and found that there is a past unstable saddle critical point corresponding to a static universe. The phasespace trajectories originating form the vicinity of this saddle like point are repelled away from it and move toward the stable critical point corresponding to an expanding universe dominated by Zel'dovich fluid.
In the second instance we considered a three dimensional phase-space case by incorporating the radiation component too. But obtaining the analytical solution in a two-component universe, especially when one of the fluids is viscous, is practically difficult or impossible, so we have made a suitable phase-space analysis. In this case no critical points are found corresponding to a prior radiation dominated phase; moreover, none of the existing critical points is stable. Hence the present model of the universe with a bulk viscous Zel'dovich fluid, in which the bulk viscosity is characterized by a constant coefficient, first of all, fails to predict a prior radiation dominated phase and, second, the very inclusion of radiation makes the model unstable.