Viscous Cosmology for Early- and Late-Time Universe

TL;DR

This paper reviews how incorporating bulk viscosity into cosmological models impacts early universe inflation, late-time acceleration, and potential singularities, offering a unified perspective on cosmic evolution.

Contribution

It provides a comprehensive analysis of viscous effects on both early and late universe dynamics, including inflation, acceleration, and singularity avoidance, using modern observational data.

Findings

Viscosity influences inflationary observables.

Viscous effects can unify inflation and late-time acceleration.

Viscosity may induce cosmological bounces and alter the quintessence-phantom divide.

Abstract

From a hydrodynamicist's point of view the inclusion of viscosity concepts in the macroscopic theory of the cosmic fluid would appear most natural, as an ideal fluid is after all an abstraction (excluding special cases such as superconductivity). Making use of modern observational results for the Hubble parameter plus standard Friedmann formalism, we may extrapolate the description of the universe back in time up to the inflationary era, or we may go to the opposite extreme and analyze the probable ultimate fate of the universe. In this review we discuss a variety of topics in cosmology when it is enlarged in order to contain a bulk viscosity. Various forms of this viscosity, when expressed in terms of the fluid density or the Hubble parameter, are discussed. Furthermore, we consider homogeneous as well as inhomogeneous equations of state. We investigate viscous cosmology in the early universe, examining the viscosity effects on the various inflationary observables. Additionally, we study viscous cosmology in the late universe, containing current acceleration and the possible future singularities, and we investigate how one may even unify inflationary and late-time acceleration. Finally, we analyze the viscosity-induced crossing through the quintessence-phantom divide, we examine the realization of viscosity-driven cosmological bounces, and we briefly discuss how the Cardy-Verlinde formula is affected by viscosity.

Tables

Table 1. Table 1 : The equation-of-state parameter of the reconstructed viscous inflationary models of ( 2.88 ), ( 2.89 ), ( 2.92 ), and ( 2.93 ), along with the conditions for a graceful exit. The parameter C 1 subscript 𝐶 1 C_{1} is always positive, while | C 2 ​ ρ | ≫ 1 much-greater-than subscript 𝐶 2 𝜌 1 \left|C_{2}\rho\right|\gg 1 and C 2 < 0 subscript 𝐶 2 0 C_{2}<0 for Models (A) and (B), and | C 2 ​ ρ | ≪ 1 much-less-than subscript 𝐶 2 𝜌 1 \left|C_{2}\rho\right|\ll 1 for Models (C) and (D). From [ 42 ] .

Case	Model	EoS	Conditions for graceful exit
(i)	(a)	$p=-\rho +\left[2C_2/\left(3C_1\right)\right]\rho -\left[3C_2^2/\left(C_1{\kappa }^2\right)\right]H^4$	No condition
(i)	(b)	$p=-\rho -\left[C_2^2/\left(3C_1\right)\right]{\rho }^2+\left[2C_2/\left(C_1\kappa \right)\right]H^2$	No condition
(ii)	(c)	$p=-\rho -\left[1/\left(3C_1\right)\right]+\left[2C_2/\left(C_1\kappa \right)\right]H^2$	or $C_2>\left(1/36\right){\left(\sqrt{\kappa }/H_{\inf }\right)}^2$
(ii)	(d)	$p=-\rho +\left[2C_2/\left(3C_1\right)\right]\rho -\left[1/\left(3C_1\right)\right]$	or $C_2>\left(1/18\right){\left(\sqrt{\kappa }/H_{\inf }\right)}^2$

Table 2. Table 2 : The present-day values of the deceleration parameter q 𝑞 q , of the dark energy equation-of-state parameter w D ​ E subscript 𝑤 𝐷 𝐸 w_{DE} and its derivative w D ​ E ′ subscript superscript 𝑤 ′ 𝐷 𝐸 w^{\prime}_{DE} , the statefinder parameters ( r , s ) 𝑟 𝑠 (r,s) and the value of the transition redshift z t ​ r subscript 𝑧 𝑡 𝑟 z_{tr} , for the non-interacting model, for several values of the viscosity parameter χ 𝜒 \chi in ( 3.118 ), and with α = 1.15 𝛼 1.15 \alpha=1.15 , α 0 = 0.15 subscript 𝛼 0 0.15 \alpha_{0}=0.15 , α 1 = 0.2 subscript 𝛼 1 0.2 \alpha_{1}=0.2 , α 2 = 0.25 subscript 𝛼 2 0.25 \alpha_{2}=0.25 . We have set H 0 = 0.7 subscript 𝐻 0 0.7 H_{0}=0.7 and Ω D ​ M = 0.27 subscript Ω 𝐷 𝑀 0.27 \Omega_{DM}=0.27 . From [ 123 ] .

$\chi$	$q$	$(w_{DE}^{\prime },w_{DE})$	$(r,s)$	$z_{tr}$
$-0.25$	$-0.766$	$(0.433,-0.952)$	$(1.356,-0.094)$	$1.21$
$-0.1$	$-0.666$	$(0.487,-0.954)$	$(1.797,-0.228)$	$0.82$
$0.0$	$-0.599$	$(0.529,-0.953)$	$\left(2.111-0.337\right)$	$0.67$
$0.1$	$-0.533$	$(0.572,-0.952)$	$(2.441,-0.464)$	$0.53$
$0.25$	$-0.433$	$(0.635,-0.952)$	$(2.965,-0.701)$	$0.4$

Table 3. Table 3 : The present-day values of the deceleration parameter q 𝑞 q , of the dark energy equation-of-state parameter w D ​ E subscript 𝑤 𝐷 𝐸 w_{DE} and its derivative w D ​ E ′ subscript superscript 𝑤 ′ 𝐷 𝐸 w^{\prime}_{DE} , the statefinder parameters ( r , s ) 𝑟 𝑠 (r,s) and the value of the transition redshift z t ​ r subscript 𝑧 𝑡 𝑟 z_{tr} , for the non-interacting model, for several values of the parameter α 𝛼 \alpha in ( 3.118 ), and with χ = − 0.1 𝜒 0.1 \chi=-0.1 , α 0 = 0.15 subscript 𝛼 0 0.15 \alpha_{0}=0.15 , α 1 = 0.2 subscript 𝛼 1 0.2 \alpha_{1}=0.2 , α 2 = 0.25 subscript 𝛼 2 0.25 \alpha_{2}=0.25 . We have set H 0 = 0.7 subscript 𝐻 0 0.7 H_{0}=0.7 and Ω D ​ M = 0.27 subscript Ω 𝐷 𝑀 0.27 \Omega_{DM}=0.27 . From [ 123 ] .

$\alpha$	$q$	$(w_{DE}^{\prime },w_{DE})$	$(r,s)$	$z_{tr}$
$0.75$	$-0.508$	$(0.587,-0.952)$	$(1.894,-0.296)$	$-$
$0.85$	$-0.554$	$(0.559,-0.952)$	$(1.977,-0.309)$	$-$
$0.95$	$-0.595$	$(0.534,-0.952)$	$\left(1.967-0.295\right)$	$1.2$
$1.15$	$-0.666$	$(0.487,-0.952)$	$(1.797,-0.227)$	$0.82$
$1.2$	$-0.682$	$(0.477,-0.952)$	$(1.736,-0.208)$	$0.82$

Table 4. Table 4 : The present-day values of the deceleration parameter q 𝑞 q , of the dark energy equation-of-state parameter w D ​ E subscript 𝑤 𝐷 𝐸 w_{DE} and its derivative w D ​ E ′ subscript superscript 𝑤 ′ 𝐷 𝐸 w^{\prime}_{DE} , the statefinder parameters ( r , s ) 𝑟 𝑠 (r,s) and the value of the transition redshift z t ​ r subscript 𝑧 𝑡 𝑟 z_{tr} , for the interacting model, for several values of the interaction parameter β 𝛽 \beta of ( 3.121 ), and with χ = 0.1 𝜒 0.1 \chi=0.1 , α = 1.15 𝛼 1.15 \alpha=1.15 , α 0 = 0.15 subscript 𝛼 0 0.15 \alpha_{0}=0.15 , α 1 = 0.2 subscript 𝛼 1 0.2 \alpha_{1}=0.2 , α 2 = 0.25 subscript 𝛼 2 0.25 \alpha_{2}=0.25 . We have set H 0 = 0.7 subscript 𝐻 0 0.7 H_{0}=0.7 and Ω D ​ M = 0.27 subscript Ω 𝐷 𝑀 0.27 \Omega_{DM}=0.27 . From [ 123 ] .

$b$	$q$	$(w_{DE}^{\prime },w_{DE})$	$(r,s)$	$z_{tr}$
$0.0$	$-0.533$	$(0.572,-0.952)$	$(2.441,-0.465)$	$0.54$
$0.01$	$-0.548$	$(0.572,-0.966)$	$(2.342,-0.427)$	$0.58$
$0.03$	$-0.578$	$(0.572,-0.993)$	$(2.152,-0.356)$	$0.65$
$0.05$	$-0.608$	$(0.568,-1.021)$	$(1.972,-0.292)$	$0.72$
$0.07$	$-0.634$	$(0.561,-1.048)$	$(1.800,-0.234)$	$0.82$

Equations

$$h_{\mu \nu }=g_{\mu \nu }+U_{\mu }{U}_{\nu },$$

$${\omega }_{\mu \nu }=h_{\mu }^{\alpha }{h}_{\nu }^{\beta }{U}_{(\alpha ;\beta )}=\frac{1}{2}(U_{\mu ;\alpha }{h}_{\nu }^{\alpha }-U_{\nu ;\alpha }{h}_{\mu }^{\alpha }).$$

$${\theta }_{\mu \nu }=h_{\mu }^{\alpha }{h}_{\nu }^{\beta }{U}_{(\alpha ;\beta )}=\frac{1}{2}(U_{\mu ;\alpha }{h}_{\nu }^{\alpha }+U_{\nu ;\alpha }{h}_{\mu }^{\alpha }),$$

$${\sigma }_{\mu \nu }={\theta }_{\mu \nu }-\frac{1}{3}{h}_{\mu \nu }\theta ,$$

$$U_{\mu ;\nu }={\omega }_{\mu \nu }+{\sigma }_{\mu \nu }+\frac{1}{3}{h}_{\mu \nu }\theta -A_{\mu }{U}_{\nu },$$

$$d{s}^2=-d{t}^2+a^2(t)\left[\frac{d{r}^2}{1-k{r}^2}+r^2(d{\theta }^2+{\sin }^2\theta d{\phi }^2)\right],$$

$$U_{\mu ;\nu }=H{h}_{\mu \nu },$$

$${\omega }_{\mu \nu }={\sigma }_{\mu \nu }=0,A_{\mu }=0,$$

$$\theta =3H.$$

$$Q^{\mu }=-K{h}^{\mu \nu }(T_{,\nu }+T{A}_{\nu }),$$

$$T_{\mu \nu }=\rho U_{\mu }{U}_{\nu }+(p-3H\zeta )h_{\mu \nu }-2\eta {\sigma }_{\mu \nu }+Q_{\mu }{U}_{\nu }+Q_{\nu }{U}_{\mu },$$

$$H^2+\frac{k}{a^2}=\frac{\kappa \rho }{3}$$

$$2\dot{H}+3H^2=-\kappa p,$$

$$\dot{H}=-(\kappa /2)(\rho +p)$$

$$p=w(\rho )\rho -B(a(t),H,\dot{H}...),$$

$$B(a(t),H,\dot{H}...)=3H\zeta (H),$$

$$\frac{dS}{dt}=\frac{2\eta }{T}({\theta }_{ik}-\frac{1}{3}{\delta }_{ik}\nabla \cdot \mathbf{u}{)}^2+\frac{\zeta }{T}{(\nabla \cdot \mathbf{u})}^2+\frac{K}{T^2}(\nabla T{)}^2,$$

$${\theta }_{ik}\to {\theta }_{\mu \nu },{\delta }_{ik}\to h_{\mu \nu },\nabla \cdot \mathbf{u}\to 3H,-K{T}_{,k}\to Q_{\mu },$$

$${S^{\mu }}_{;\mu }=\frac{2\eta }{T}{\sigma }_{\mu \nu }{\sigma }^{\mu \nu }+\frac{9\zeta }{T}{H}^2+\frac{1}{K{T}^2}{Q}_{\mu }{Q}^{\mu },$$

$$S^{\mu }=n{k}_B\sigma U^{\mu }+\frac{1}{T}{Q}^{\mu }.$$

$$H^2=\frac{\kappa }{3}\rho =\frac{\kappa }{3}\left(\frac{1}{2}{\dot{\varphi }}^2+V\right),$$

$$\frac{\ddot{a}}{a}=-\frac{\kappa }{6}(\rho +3p),$$

$$\rho =\frac{1}{2}{\dot{\varphi }}^2+V,$$

$$p=\frac{1}{2}{\dot{\varphi }}^2-V,$$

$$p=w\rho ,$$

$$w=\frac{\frac{1}{2}{\dot{\varphi }}^2-V}{\frac{1}{2}{\dot{\varphi }}^2+V}.$$

$$\ddot{\varphi }+3H\dot{\varphi }=-V^{\prime },$$

$$\dot{\rho }+3H(\rho +p)=0.$$

$$N=\ln \left(\frac{a_f}{a}\right).$$

$$\epsilon$$

$$\eta$$

$$\xi$$

Full text

1 Introduction

The introduction of viscosity coefficients in cosmology has itself a long history, although the physical importance of these phenomenological parameters has traditionally been assumed to be weak or at least subdominant. In connection with the very early universe, the influence from viscosity is assumed to be the largest at the time of neutrino decoupling (end of the lepton era), when the temperature was about $10^{10}~{}$K. Misner [1] was probably the first to introduce the viscosity from the standpoint of particle physics; see also Zel’dovich and Novikov [2]. Nevertheless, on a phenomenological level, the viscosity concept was actually introduced much earlier, with the first such work being that of Eckart [3].

When considering deviations from thermal equilibrium to the first order in the cosmic fluid, one should recognize that there are in principle two different viscosity coefficients, namely the bulk viscosity $\zeta$ and the shear viscosity $\eta$. In view of the commonly accepted spatial isotropy of the universe, one usually omits the shear viscosity. This is motivated by the WMAP [4] and Planck observations [5], and is moreover supported by theoretical calculations which show that in a large class of homogeneous and anisotropic universes isotropization is quickly established. Eckart’s theory, as most other theories, is maintained at first-order level. In principle, a difficulty with this kind of theory is that one becomes confronted with a non-causal behavior. In order to prevent this one has to go to the second order approximation, away from thermal equilibrium.

The interest in viscosity theories in cosmology has increased in recent years, for various reasons, perhaps especially from a fundamental viewpoint. It is well known among hydrodynamicists that the ideal (nonviscous) theory is after all only an approximation to the real world. For reviews on both causal and non-causal theories, the reader may consult Grøn [6] (surveying the literature up to 1990), and later treatises by Maartens [7, 8], and Brevik and Grøn [9].

The purpose of the present review is to explore how several parts of cosmological theory become affected when a bulk viscosity is brought into the formalism. After highlighting the basic formalism in the remaining of the present section, in Section 2 we consider the very early (inflationary) universe. We briefly present the conventional inflation theory, covering “cold”, “warm” and “intermediate” inflation, and we extract various inflationary observables. Thereafter we investigate the viscous counterparts in different models, depending on the form of bulk viscosity as well as on the equation of state.

In Section 3 we turn to the late universe, including the characteristic singularities in the far future, related also to the phantom region in which the equation-of-state parameter is less than $-1$. The different types of future singularities are classified, and we explore the consequences of letting the equation of state to be inhomogeneous. A special case is the unification of inflation with dark energy in the presence of viscosity, a topic which is dealt with most conveniently when one introduces a scalar field. Additionally, we discuss holographic dark energy in the presence of a viscous fluid.

In Section 4 we discuss various special topics, amongst them the possibility for the viscous fluid to slide from the quintessence region into the phantom region and then into a future singularity, if the magnitude of the present bulk viscosity is large enough. Comparison with estimated values of the bulk viscosity derived from observations, indicate that this may actually be a realistic scenario. In the same section we also discuss the viscous Big Rip realization, and finally we see how the Cardy-Verlinde formula becomes generalized when viscosity is accounted for, since the thermodynamic (emergent) approach to gravity has become increasingly popular.

Finally, in Section 5 we summarize the obtained results and we discuss on the advantages of viscous cosmology.

1.1 Basic formalism

We begin by an outline of the general relativistic theory, setting, as usual, $k_{B}$ and $c$ equal to one. The formalism below is taken from Ref. [10]. We adopt the Minkowski metric in the form $(-+++)$, and we use Latin indices to denote the spatial coordinates from 1 to 3, and Greek indices to denote spacetime ones, acquiring values from 0 to 3. $U^{\mu}=(U^{0},U^{i})$ denotes the four-velocity of the cosmic fluid, and we have $U^{0}=1,U^{i}=0$ in a local comoving frame.

With $g_{\mu\nu}$ being a general metric tensor we introduce the projection tensor

[TABLE]

and the rotation tensor

[TABLE]

The expansion tensor is

[TABLE]

and has the trace $\theta\equiv\theta_{\mu}^{\mu}={U^{\mu}}_{;\mu}$. The third tensor that we shall introduce is the shear tensor, namely

[TABLE]

which satisfies $\sigma_{\mu}^{\mu}=0$. Finally, it is often useful to make use of the three tensors above in the following decomposition of the covariant derivative of the fluid velocity:

[TABLE]

where $A_{\mu}$ stands for the four-acceleration, namely $A_{\mu}={\dot{U}}_{\mu}=U^{\nu}U_{\mu;\nu}$.

The above formalism is for a general geometry. In the following we will focus on Friedmann-Robertson-Walker (FRW) geometry, which is of main interest in cosmology, whose line element is

[TABLE]

where $a(t)$ is the scale factor and $k=1,0,-1$ the spatial curvature parameter. In this case the coordinates $x^{\mu}$ are numerated as $(t,r,\theta,\varphi)$. In these coordinates the covariant derivatives of the velocity acquire the simple form

[TABLE]

with $H=\dot{a}/a$ the Hubble parameter. The rotation tensor, the shear tensor, and the four-acceleration all vanish, i.e

[TABLE]

and the relation between scalar expansion and Hubble parameter is simply

[TABLE]

As a next step we consider the fluid’s energy-momentum tensor $T_{\mu\nu}$ in the case where viscosity as well as heat conduction are taken into account. If $K$ is the thermal conductivity, considered in its nonrelativistic framework, then for the spacelike heat flux density four-vector we have the expression

[TABLE]

with $T$ the temperature. The last term in this expression is of relativistic origin. The coordinates used in (1.1) are comoving, with freely moving reference particles having vanishing four-acceleration. Thus, one obtains the usual expression $Q_{\hat{i}}=-KT_{,\hat{i}}$ for the heat flux density through a surface orthogonal to the unit vector ${\bf{e}}_{\hat{i}}$. Hence, assembling everything, in an FRW metric we can now introduce the energy-momentum tensor as

[TABLE]

with $\rho$ and $p$ the fluid’s energy density and pressure respectively, and where $\zeta$ is the bulk viscosity and $\eta$ the shear viscosity.

Taking all the above into consideration, we conclude that for a universe governed by General Relativity in the presence of a viscous fluid, in FRW geometry the two Friedmann equations read as:

[TABLE]

[TABLE]

with $\kappa$ the gravitational constant. Note that these equations give

[TABLE]

for a flat universe. We mention that the energy density and pressure can acquire a quite general form. For instance, a quite general parametrization of an inhomogeneous viscous fluid in FRW geometry is [11, 12, 13]

[TABLE]

where $w(\rho)$ can depend on the energy density, and the bulk viscosity $B(a(t),H,\dot{H}...)$ can be a function of the scale factor, and of the Hubble function and its derivatives. A usual subclass of the above general equation of state is to assume that

[TABLE]

with $\zeta(H)>0$ the bulk viscosity, which can be further simplified to the subclass with $\zeta(H)=\zeta=const.$.

Let us now focus on thermodynamics, and especially on the production of entropy. The simplest way of extracting the relativistic formulae is to generalize the known formalism from nonrelativistic thermodynamics. We use $\sigma$ to denote the dimensionless entropy per particle, where for definiteness as “particle” we mean a baryon. The nonrelativistic entropy density thus becomes $nk_{B}\sigma$, where $n$ is the baryon number density. Making use of the relationship [14]

[TABLE]

where $\bf u$ denotes the nonrelativistic velocity and $\nabla$ the three-dimensional Laplace operator, we can generalize to a relativistic formalism simply by imposing the effective substitutions

[TABLE]

whereby we obtain the desired equation

[TABLE]

in which $S^{\mu}$ denotes the entropy current four-vector

[TABLE]

More detailed derivations of these results can be found, for instance, in Refs. [15] and [16].

In summary, viscous cosmology is governed by the Friedmann equations (1.12) and (1.13), along with various considerations of the fluid’s equation of state. Hence, these relations will be the starting point of the discussion of this manuscript. In the following sections we investigate viscous cosmology in detail.

2 Inflation

We start the investigation of viscous cosmology by focusing on early times, and in particular on the inflationary realization. Inflation is considered to be a crucial part of the universe cosmological history, since it can offer a solution to the flatness, horizon and monopole problems [17, 18, 19]. In order to obtain the inflationary phase one needs to consider a suitable mechanism, which is either a scalar field in the framework of General Relativity [20, 21, 22], or a degree of freedom arising from gravitational modification [23, 24]. In this section we will see how inflation can be driven by a viscous fluid.

2.1 Inflation: The basics

Before proceeding to the investigation of viscous inflation, let us briefly describe the basic inflationary formulation and the relation to various observables. For convenience we review the scenarios of cold and warm inflation separately.

• •

Cold Inflation

We first start with the standard inflation realization, also called as “cold” inflation, in which a scalar field $\phi$ plays the role of the inflaton field. The Friedmann equations are

[TABLE]

[TABLE]

where $\rho$ and $p$ are respectively the energy density and pressure of the inflaton field, and $V=V(\phi)$ is the corresponding potential. In (2.1) we have used the fact that the scalar field can be viewed as a perfect fluid with

[TABLE]

and hence its equation-of-state (EoS) parameter reads

[TABLE]

The fluid interpolates between an invariant vacuum energy with $w=-1$ for a constant inflaton field, and a stiff (Zel’dovich) fluid with $w=1$ and $V=0$.

The scalar-field equation of motion takes the simple form

[TABLE]

where $V^{\prime}=dV/d\phi$, which can be re-written as a continuity equation

[TABLE]

Finally, we can define the quantity $N$, i.e. the number of e-folds in the slow-roll era, as the logarithm of the ratio between the final value $a_{f}$ of the scale factor during inflation and the initial value $a(N)=a$, namely

[TABLE]

In inflationary theory it proves very convenient to define the so-called slow roll parameters. One set of such parameters is defined via derivatives of the potential with respect to the inflaton field. These “potential” slow roll parameters, conventionally called $\varepsilon$, $\eta$, $\xi$, are defined as [25]

[TABLE]

Since these should be small during the slow-roll period, the potential $V(\phi)$ must have a flat region.

One may also define the slow roll parameters in a different way, by taking the derivatives of the Hubble parameter with respect to the e-folding number (such an approach has a more general application, since it can be also used in inflationary realizations that are driven from modified gravity, where a field and a potential are absent) [25]. In particular, these horizon-flow [26, 27, 28] parameters $\epsilon_{n}$ (with $n$ a positive integer), are defined as

[TABLE]

with $\epsilon_{0}\equiv H_{ini}/H$ and $N$ the e-folding number, and $H_{ini}$ the Hubble parameter at the beginning of inflation (inflation ends when $\epsilon_{1}=1$). Thus, the first three of them are calculated as

[TABLE]

We now briefly review the formalism that is used to describe the temperature fluctuations in the Cosmic Microwave Background (CMB) radiation. The power spectra of scalar and tensor fluctuations are written as [29]

[TABLE]

with

[TABLE]

Here $k$ is the wave number of the perturbation, and $k_{*}$ is a reference scale usually chosen as the wave number at horizon crossing (the pivot scale). Often one chooses $k=\dot{a}=aH$, with $a$ the scale factor. The quantities $A_{s}$ and $A_{T}$ are amplitudes at the pivot scale, while $n_{s}$ and $n_{T}$ are called the spectral indices of scalar and tensor fluctuations. Moreover, $-\delta_{ns}=n_{s}-1$ and $n_{T}$ are called the tilts of the power spectrum, since they describe deviations from the scale invariant spectrum where $\delta_{ns}=n_{T}=0$. The factors $\alpha_{s}$ and $\alpha_{T}$ are called running spectral indices and are defined by

[TABLE]

Finally, the tensor-to-scalar ratio $r$ is defined as

[TABLE]

Analysis of the observations from the Planck satellite give the result $n_{s}=0.968(6)\pm 0.006$ [30, 5]. Furthermore, the observations give $\alpha_{s}=-0.003\pm 0.007$. The tilt of the curvature fluctuations is $\delta_{ns}=0.032$. The combined BICEP2/Planck and LIGO data give $n_{T}={-0.76}^{+1.37}_{-0.52}$ [31], while the BICEP/Planck data alone constrain the tensor tilt to be $n_{T}={0.66}^{+1.83}_{-1.44}$.

From the above equations we derive

[TABLE]

where the quantities are evaluated at the horizon crossing ($k=k_{*}$), and as we mentioned $k=aH$. Hence, we can finally extract the expressions that relate the inflationary observables, namely the tensor-to-scalar ratio, the scalar spectral index, the running of the scalar spectral index, and the tensor spectral index, with the potential-related slow-roll parameters (2.10)-(2.12), which read as [25]:

[TABLE]

Hence, a consistency relation between $r$ and $n_{T}$ follows from Eqs. (2.13), (2.16) and (2.19), namely $n_{T}=-\frac{r}{8}$. The preferred BICEP2/Planck value of $r=0.05$ then gives $n_{T}=-0.006$.

Lastly, when the horizon flow slow-roll parameters are used, the inflationary observables read as [25]

[TABLE]

Definitely, in cases where both the potential slow-roll parameters and the horizon flow slow-roll parameters can be used, the final expressions for the observables coincide.

• •

Warm Inflation

Let us now briefly review the scenario of “warm” inflation. Usually, one is concerned with cold inflationary models described above, for which dissipation arising from the decay of inflaton energy to radiation is omitted. Nevertheless, this contrasts the characteristic feature of the so-called warm inflation, where dissipation is included as an important factor, and inflaton energy dissipates into heat [32, 33, 34, 35]. This implies in turn that the inflationary period lasts longer than it does in the cold case. Additionally, no reheating at the end of the inflationary era is needed, and the transition to radiation era becomes a smooth one.

The main characteristic for the warm inflationary models is that the inflaton field energy $\rho_{\phi}$ is considered to depend on the temperature $T$ [36], in a same way as the radiation density $\rho_{r}$ depends on $T$. The first Friedmann equation writes as

[TABLE]

and the continuity equations for the two fluid components read

[TABLE]

where $\Gamma$ is a dissipation coefficient describing the transfer of dark energy into radiation and it is in general time dependent. In warm inflation the inflaton energy is the dominating component, $\rho_{\phi}\gg\rho_{r}$, and $H$, $\phi$ and $\Gamma$ vary slowly such that $\ddot{\phi}\ll H\dot{\phi}$, $\dot{\rho}_{r}\ll 4H\rho_{r}$ and $\dot{\rho}_{r}\ll\Gamma\dot{\phi}^{2}$. In the slow roll epoch, the radiation is produced by dark energy dissipation. Thus

[TABLE]

Defining the so-called dissipative ratio by

[TABLE]

we see that in the warm inflation era Eq. (2.30) yields

[TABLE]

During warm inflation $T>H$ (in geometric units), and it turns out that the tensor-to-scalar ratio is modified in comparison to the cold inflation case, namely [37]

[TABLE]

and thus this ratio is suppressed by a factor $(T/H)(1+Q)^{5/2}$ compared to the cold inflationary case.

The slow roll parameters in the present models are calculated at the beginning $t=t_{i}$ of the slow roll epoch. From the definition equation (2.8) we acquire

[TABLE]

Comparing with (2.10) we see that the first slow-roll parameter of the warm inflation scenario is modified with the factor $1+Q$ relative to the corresponding cold inflation parameter. Furthermore, manipulation of the above equations then yields for the parameter $\eta$

[TABLE]

For convenience we introduce the quantity $\beta=\Gamma^{\prime}V^{\prime}/(\kappa\Gamma V)$, and therefore this quantity appears in the expression for the relative rate of change of the radiation energy density, namely

[TABLE]

Introducing also

[TABLE]

one can find that [38]

[TABLE]

When warm inflation is strong, $Q\gg 1$ , $\omega\gg 1$, and thus

[TABLE]

whereas when it is weak, $Q\ll 1$, and therefore

[TABLE]

Finally, the cold inflationary case corresponds to the limit $Q\to 0$ and $T\ll H$, and then $\omega\to 0$ and

[TABLE]

Visinelli found the following expression for tensor-to-scalar ratio in warm inflation [38]:

[TABLE]

Hence, in the limit of strong dissipative warm inflation we have

[TABLE]

while in the limit of cold inflation we re-obtain the standard result (2.20), namely $r\to 16\varepsilon$. Thus, the warm inflation models with $Q\gg 1$ and $\omega\gg 1$ yield a very small tensor-to-scalar ratio.

• •

Intermediate inflation

Intermediate inflation scenario, introduced by Barrow in 1990 [39] (see also [40, 41]), also uses a scalar field. We consider the scale factor to take the form

[TABLE]

with $0<\alpha<1$, and where $A$ is a positive dimensionless constant, while $a_{p}$ refers to the Planck time ($\hat{t}=t/\sqrt{\kappa},\ t_{p}=\sqrt{\kappa}$). The reason that these models are called intermediate, is that the expansion is faster than the corresponding one in power-law inflation, and slower than an exponential inflation (the latter corresponding to $\alpha=1$). It follows from (2.46) that

[TABLE]

and since $\dot{H}<0$ for $\alpha<1$, the Hubble parameter decreases with time. Inserting these equations into Eqs. (2.1) and (2.2) we obtain

[TABLE]

Since $\rho+p=\dot{\phi}^{2}$ we obtain by integration, using the initial condition $\phi(0)=0$, that

[TABLE]

while since $V=\frac{1}{2}(\rho-p)$ we acquire

[TABLE]

Hence, eliminating $t$ between (2.49) and (2.50) we can express the potential as a function of the inflaton field:

[TABLE]

For this class of models the spectral parameters are most easily calculated from the Hubble slow roll parameters

[TABLE]

The optical parameters $\delta_{ns},\,n_{r}$ and $r$ can be expressed in terms of the Hubble slow roll parameters to lowest order as

[TABLE]

This gives

[TABLE]

The slow roll parameter $\varepsilon_{H}$ can be expressed in terms of the inflaton field as

[TABLE]

In the intermediate inflation, the e-folding number becomes

[TABLE]

where $\hat{t}_{i}$ and $\hat{t}_{f}$ are the initial and final point of time of the inflationary era, respectively. In these models the beginning of the inflationary era is defined by the condition $\varepsilon_{H}(\hat{t}_{i})=1$, giving

[TABLE]

Hence, the inflationary era ends at a point of time

[TABLE]

The slow roll parameters are evaluated at this point of time, giving

[TABLE]

Inserting the above expressions into (2.53) we can thus write

[TABLE]

Note that the curvature spectrum is scale independent, corresponding to $n_{s}=1$, for $\alpha=2/3$. Furthermore, $n_{s}<1$ requires $\alpha<2/3$. Note that the expression for $n_{s}$ corrects an error of Ref. [40]. For these models the $r,\,\delta_{ns}$ relation becomes

[TABLE]

The constant $\alpha$ can be expressed in terms of $N$ and $\delta_{ns}$ as

[TABLE]

With the Planck values $\delta_{ns}=0.032$ and $N=60$ we get $\alpha=0.4$ giving $r=0.38$. This value of $r$ is larger than permitted by Planck observations. However, the more general models with non-canonical inflaton fields studied in Refs. [40] and [41], contain an adjustable parameter in the expressions for the observables, leading to agreement with observational data. Below we shall consider warm intermediate inflation models, which lead naturally to a suppression of the curvature perturbation, resulting to a small value of $r$.

2.2 Viscous Inflation

Having described the basics of inflation, in this subsection we will see how inflation can be realized in the framework of viscous cosmology, that is if instead of a scalar field inflation is driven by a viscous fluid [42]. We start from the two Friedmann equations (1.12) and (1.13), namely

[TABLE]

[TABLE]

Concerning the viscosity of the fluid we consider a subclass of (1.15) and parametrize the equation of state as

[TABLE]

with $A$,$\beta$ constants, and $\zeta(H)$ the bulk viscosity considered with a dynamical nature in general, i.e. being a function of the Hubble parameter. As a specific example we consider

[TABLE]

with $\bar{\zeta}$, $\gamma$ parameters.

From the Friedmann equation (2.63) and for an expanding flat universe ($H>0,k=0$), we acquire

[TABLE]

Therefore, $\zeta(H)$ can be expressed in terms of $\rho$, i.e $\zeta(H)=\zeta(H(\rho))$. Thus, comparing the general expression for the EoS of a fluid, namely

[TABLE]

with (2.65) and (2.66), we deduce that

[TABLE]

We mention that $f(\rho)$ is expressed as a series of powers in $\rho$ due to the imposed assumption that $\zeta(H)$ is a power of $H$. Hence, this allows us to find analytical solutions and examine the behavior of various inflationary observables.

Since in fluid inflation we do not have a potential, it proves convenient to use the Hubble slow-roll parameters. Inserting the Hubble function from (2.67) into (2.10)-(2.12) and then into the inflationary observables (2.24)-(2.27), after some algebra one can express the tilt, the tensor-to-scalar ratio and the running spectral index as [42]

[TABLE]

where we have also used that $f(\rho)/\rho(N)=w(N)+1$. In these expressions all quantities may be considered as functions of the e-folding number $N$. Hence, if we choose $f(\rho)/\rho(N)=4.35\times 10^{-3}$, we obtain $w=-0.996$, and thus $(n_{\mathrm{s}},r,\alpha_{s})=(0.974,0.104,-1.70\times 10^{-4})$. These results are consistent with the Planck data, namely $n_{\mathrm{s}}=0.968\pm 0.006\,(68\%\,\mathrm{CL})$, $r<0.11\,(95\%\,\mathrm{CL})$, and $\alpha_{s}=-0.003\pm 0.007\,(68\%\,\mathrm{CL})$, [5, 43].

Let us now use the required scalar spectral index in order to reconstruct the EoS of the fluid through a corresponding effective potential, following [44, 42]. In order to achieve this, we first express the Friedmann equations using derivatives in terms of the e-folding number $N$ as

[TABLE]

and similarly for the slow-roll parameters (2.8a)-(2.8c), namely

[TABLE]

Hence, one can use these quantities in order to reconstruct the equation-of-state of the corresponding fluid. In particular, having the $\delta_{ns}(N)$ as a function of $N$, using (2.2) we can solve for $V(N)$, which will be the effective potential in an equivalent scalar-field description. Then the Hubble function is related to $V(N)$ through (2.72), and thus we obtain $H=H(N)$. Finally, using (2.73) we can reconstruct $f(\rho)$ through (2.68).

Let us give a specific example of the above method, in the case where [42]

[TABLE]

which is valid in Starobinsky inflation [45], and it can be satisfied in chaotic inflation [20], in new Higgs inflation [46, 47], and in models of $\alpha$-attractors [48, 49], too. Combining (2.75) and (2.2) gives

[TABLE]

where $C_{1}(>0)$ and $C_{2}$ are constants. Hence, using (2.76) and (2.2) we acquire

[TABLE]

and thus

[TABLE]

If $\delta_{ns}=0.032,\,r=0.05$ one gets $C_{2}/C_{1}\approx 0.025$. Finally, from (2.74) we find that

[TABLE]

Thus, inserting a reasonable value $N=60$ we obtain $\alpha_{s}=-5.56\times 10^{-4}$, in agreement with Planck analysis.

In the case of a fluid model one uses the equation-of-state parameter from (2.68) instead of the scalar potential. Hence, one can have $\left(3/\kappa\right)\left(H(N)\right)^{2}=\rho(N)\approx V(N)$, since the last approximation arises from the slow-roll condition that the kinetic energy is negligible comparing to the potential one. Therefore, using (2.76) we obtain

[TABLE]

with $\left(C_{1}/N\right)+C_{2}>0$. Additionally, inserting $\rho\approx V$ into (2.76) results to

[TABLE]

Thus, inserting (2.80) into (2.72) and (2.73) gives

[TABLE]

Finally, comparing (2.68) with (2.82) leads to

[TABLE]

where we have also used (2.72) and (2.81).

We now focus on fluid inflationary models with $n_{s}$ and $r$ in agreement with observations. From (2.68) and (2.69) we obtain

[TABLE]

Therefore, we suitably choose the model parameters $A$, $\bar{\zeta}$, $\beta$, and $\gamma$, in order for relation (2.75) to be satisfied. For convenience we will focus in the regimes $\left|C_{2}\rho\right|\gg 1$ and $\left|C_{2}\rho\right|\ll 1$ separately following [42].

• •

Case I: $\left|C_{2}\rho\right|\gg 1$

In this case expression (2.83) leads to

[TABLE]

with $C_{2}<0$ in order to have a positive $N$ from (2.81). From (2.84) and (2.85) we acquire

[TABLE]

where we have also used that $\left(-C_{2}\right)/C_{1}\approx 1/N$. For instance, if $\left|C_{2}\rho\right|=\mathcal{O}(10)$, $\left(-C_{2}\right)/C_{1}\approx 1/N$, and $N\gtrsim 60$, relation (2.86) leads to $w\approx-1$, and hence the de Sitter inflation can be realized, with a scale-factor of the form

[TABLE]

It should be noted that for $\left(-C_{2}\right)/C_{1}<1/N$, relation (2.77) for $N\gtrsim 73$ provides a tensor-to-scalar ratio $r>1$, in disagreement with observations.

Comparing (2.85) and (2.69) we deduce that we obtain equivalence for two combinations of parameters:

[TABLE]

and

[TABLE]

Hence, the corresponding fluid equation of state can be reconstructed.

• •

Case (II): $\left|C_{2}\rho\right|\ll 1$

In this case expression (2.83) leads to

[TABLE]

Thus, (2.81) with $\left|C_{2}\rho\right|\ll 1$ gives $C_{1}\rho\approx N\gg 1$ and thus $\left|C_{2}\right|/C_{1}\ll 1$. Hence, (2.84) and (2.85), give

[TABLE]

where we have used that $C_{1}\rho\approx N$. Similarly to the previous subcase, (2.91) with $1/N\ll 1$ and $\left|C_{2}\right|/C_{1}\ll 1$, leads to $w\approx-1$, i.e to the realization of the de Sitter inflation, with a scale factor given by (2.87). Moreover, for $C_{2}>0$ and $C_{2}/C_{1}\lesssim 1/N$, and for $N\gtrsim 60$, relation (2.77) gives $r<0.11$ in agreement with Planck results. On the other hand, for $C_{2}<0$ and $\left|C_{2}\right|/C_{1}<1/N$, we need to have $N\gtrsim 73$ in order to get $r<0.11$, similarly to the previous Case (I). Finally, comparing (2.90) and (2.69) we deduce that we obtain equivalence for two combinations of parameters:

[TABLE]

and

[TABLE]

Having analyzed the basic features of inflationary realization from a viscous fluid, let us examine the crucial issue of obtaining a graceful exit and the subsequent entrance to the reheating stage [42]. In particular, we will investigate the instability of the de Sitter solution characterized by $H=H_{\mathrm{inf}}=const.$ under perturbations. One starts by perturbing the Hubble function as [50]

[TABLE]

where $\left|\delta(t)\right|\ll 1$. Thus, the second Friedmann equation writes as a differential equation in terms of the cosmic time $t$, namely

[TABLE]

Without loss of generality we choose

[TABLE]

with $\lambda$ a constant, and therefore a positive $\lambda$ would correspond to an unstable de Sitter solution. This instability implies that the universe can exit from inflation. On the other hand, a stable inflationary solution is just an eternal inflation.

Inserting (2.94) and (2.96) into (2.95), and keeping terms up to first order in $\delta(t)$, we obtain

[TABLE]

with

[TABLE]

The solutions of (2.97) read as

[TABLE]

and therefore if $\mathcal{Q}>0$ we may obtain $\lambda=\lambda_{+}>0$, which implies the realization of a successful inflationary exit.

Let us now check whether the four fluid models described in (2.88), (2.89), (2.92), and (2.93) above, can give rise to a graceful exit, i.e whether they can give a positive $\mathcal{Q}$ in (2.98). Substituting the corresponding values of $A$, $\bar{\zeta}$, $\beta$, and $\gamma$ into (2.99), we obtain the expressions of $\mathcal{Q}$ as [42]:

[TABLE]

[TABLE]

[TABLE]

[TABLE]

Hence, Models (A) and (B) have always $\mathcal{Q}>0$. On the other hand, Models (C) and (D) have $\mathcal{Q}>0$ for f $C_{2}<0$, while for $C_{2}>0$ they have $\mathcal{Q}>0$ if

[TABLE]

In summary, we can see that the models of viscous fluid inflation can have a graceful exit without any tuning. In Table 1 we summarize the obtained results. From the corresponding equation-of-state parameters, and comparing with (2.65), we can immediately see the term inspired by the bulk viscosity. Finally, as we described in detail above, in these models the inflationary observables are in agreement with observations. In particular, the spectral index from (2.75) is $n_{s}=0.967$ for $N=60$. The running of the spectral index is given by $\alpha_{s}=-2/N^{2}$ in (2.79), leading to $\alpha_{s}=-5.56\times 10^{-4}$.

We close this subsection by studying the singular inflation in the above viscous fluid model. The finite-time singularities are classified into four types [51], and hence one can see that Type IV singularity can be applied in singular inflation since there are no divergences in the scale factor and in the the effective (i.e. total) energy density and pressure. In particular, in Type IV singularity, as $t\to t_{\mathrm{s}}$, with $t_{\mathrm{s}}$ the singularity time, we have $a\to a_{\mathrm{s}}$, $\rho\to 0$ and $\left|p\right|\to 0$. Here, $a_{\mathrm{s}}$ is the value of $a$ at $t=t_{\mathrm{s}}$. Nevertheless, the higher derivatives of the Hubble function diverge.

Let us consider the above viscous fluid inflationary realization, assuming that

[TABLE]

with $\bar{H}$, $q$, and $\bar{a}$ the model parameters. From the two Friedmann equations (2.63),(2.64) we straightforwardly acquire

[TABLE]

Therefore, a Type IV singularity appears at $t=t_{\mathrm{s}}$, since (2.107) and (2.108) imply that as $t\to t_{\mathrm{s}}$ the quantities $a$, $\rho$, and $p$ asymptotically approach finite values, while from (2.106) we deduce that higher derivatives of $H$ diverge. From (2.106) and (2.108) we find the following equation-of-state parameter of the cosmic fluid:

[TABLE]

with

[TABLE]

In the case where $H_{\mathrm{inf}}/\sqrt{\kappa\rho/3}=H_{\mathrm{inf}}/H\ll 1$ we have

[TABLE]

Thus, one can clearly see from (2.111) that the function $f(\rho)$ includes a linear combination of two powers of $\rho$, as in (2.69) and (2.85). Hence, indeed this scenario can be realized by the viscous fluid models reconstructed above.

From (2.111), using (2.108), we find

[TABLE]

and therefore for $\bar{H}/H^{q+1}\ll 1$ we get $f(\rho)/\rho\ll 1$. Thus, $n_{s}$, $r$, and $\alpha_{s}$ can be approximately given by (2.70) and be in agreement with observations, which act as an additional advantage of singular inflation.

We now examine the limit $\bar{\zeta}=0$ in (2.65), in which the fluid equation of state in (2.65) becomes $p=-\rho+A\rho^{\beta}$. In this limit from (2.69) we deduce that $f(\rho)=A\rho^{\beta}$, i.e $f(\rho)$ has only one power of $\rho$. However, from (2.111) and (2.112) we see that $f(\rho)$ consists of two $\rho$ powers. Thus, $f(\rho)$ can be given by (2.111) and (2.112) only if the singular inflation is realized. Hence, for a non-viscous fluid, i.e. for a fluid without the $\zeta(H)$-term in (2.65), singular inflation cannot be realized. From this feature we can see the importance of the viscous term, and its significant effect on the dynamics of the early universe. This important issue will be studied in more detail in the following subsection.

In summary, in the present subsection we studied the realization of inflation in a fluid framework, whose equation-of-state parameter has an additional term corresponding to bulk viscosity. Firstly, we saw that the obtained inflationary observables, namely $n_{s}$, $r$ and $\alpha_{s}$, are in agreement with Planck data. Secondly, we presented a reconstruction procedure of the fluid’s equation of state, when a specific $n_{s}$ is given, while the tensor-to-scalar ratio is still in agreement with observations. Thirdly, we analyzed the stability of the inflationary, de Sitter phase, showing that a graceful exit and the pass to the subsequent thermal history of the universe is obtained without fine tuning. Finally, we investigated the realization of singular inflation, corresponding to Type IV singularity, in the present viscous fluid model. Hence, viscous fluid inflation can be a candidate for the description of early universe.

2.3 Viscous warm and intermediate inflation

In this subsection we show how the viscous cold inflationary models considered above can be generalized to the warm case. These kind of models are most likely more physical than the idealized cold ones, since they take into account the presence of massive particles produced from the decaying inflaton field. Moreover, an important advantage of warm scenarios is that they give rise to a much smaller tensor-to-scalar ratio than the cold models, and hence are easier to be in agreement with the Planck data. The presence of massive particles provides a natural way to explain why the cosmic fluid can be associated with a bulk viscosity.

We abstain from using the simple equation of state $p=(1/3)\rho$ holding for radiation, and we assume instead the more general form $p=w\rho$, where $w$ is constant. For convenience, one can introduce the form $p=(\gamma-1)\rho$ with $\gamma=1+w$. The effective pressure becomes $p_{\rm eff}=p+p_{\zeta}$, where

[TABLE]

is the viscous part of the pressure and $\zeta$ the bulk viscosity. In this case Eq. (2.30) generalizes to [52]

[TABLE]

The usual condition about quasi-stationarity implies $\dot{\rho}\ll 3H(\gamma\rho-3\zeta H)$ and $\dot{\rho}\ll\Gamma{\dot{\phi}}^{2}$.

We will henceforth follow the formalism of [52] for the strong dissipative case, namely for $Q\gg 1$ (see also [53]). As it was mentioned in subsection 2.1 above, in intermediate inflation the scale factor and the Hubble parameters are given by (2.46) and (2.47). We will base the analysis on the basic assumptions

[TABLE]

where the proportionality of $\zeta$ to $\rho$ is a frequently used assumption (a similar analysis can be performed for the case where $\Gamma$ and $\zeta$ are assumed constants [52, 53], however we will not go into further details and focus on the general case). From Eqs. (2.31) and (2.33) we then have $Q=\sqrt{\kappa}H$. Focusing on the strong dissipative case $Q\gg 1$, manipulation of the equations gives the following expression for the inflaton field as a function of time:

[TABLE]

This equation, predicting $\phi(t)$ to increase with time, is seen to be different from the corresponding Eq. (2.49) for cold intermediate inflation.

Taking into account the expression (2.47) for $H$ we can express the potential as a function of time:

[TABLE]

which can alternatively be represented as a function of $\phi$ as

[TABLE]

Since

[TABLE]

we see that it is necessary for the constant $\zeta_{1}$ in (2.115) to satisfy the condition $\zeta_{1}<\gamma/3H$ in order to make $\rho$ positive. The density varies with time as

[TABLE]

while when considered as a function of the inflaton field it reads

[TABLE]

Additionally, the number of e-folds becomes in this case

[TABLE]

where $\phi_{f}$ is the inflaton field at the end of the slow-roll epoch. Finally, the slow-roll parameters in the strong dissipative epoch ($Q\gg 1$) become

[TABLE]

giving in turn for the spectral parameter $\delta_{ns}$

[TABLE]

Hence, the Harrison-Zel’dovich spectrum (independent of scale) corresponds to $\alpha=2/3$.

2.4 Singular inflation from fluids with generalized equation of state

In the end of subsection 2.2 we presented a brief discussion on the possibility to realize singular inflation in the framework of viscous cosmology. Since this is an important issue, in this subsection we investigate it in detail following [54], considering more general viscous equation of states. We consider an inhomogeneous viscous equation-of-state parameter of the form

[TABLE]

which is a subclass of the general ansatz (1.15). Thus, when the function $G(H)$ becomes zero we re-obtain the homogeneous case. An even more general equation of state would be to consider

[TABLE]

In the following we desire to investigate the realization of type IV singularity in inflation driven by a fluid with the above EoS’s.

As we mentioned earlier, a type IV singularity occurs at $t\to t_{s}$, if the scale factor and the effective energy density and pressure remain finite, but the higher derivatives of the Hubble function diverge. A general form of the Hubble function which can describe a Type IV singularity reads as

[TABLE]

with $f_{1}(t)$, $f_{2}(t)$ being arbitrary differentiable functions. Hence the type IV singularity occurs when $\alpha>1$, and without loss of generality we can consider it to take the form

[TABLE]

with $n$, $m$ positive integers.

Let us start from a simple example of type IV singularity realization, namely we consider $f_{1}(t)=0$ and $f_{2}(t)=f_{0}$, with $f_{0}$ a positive parameter. In this case the two Friedmann equations, namely $\rho=\frac{3}{\kappa}H^{2}$ and $p=-\frac{1}{\kappa}\left(3H^{2}+2\dot{H}\right)$, become

[TABLE]

and hence eliminating $t_{s}-t$ we get the result

[TABLE]

Hence, a viscous fluid with this equation of state can generate the Hubble function (2.127) and hence the type IV singularity. Defining $\tilde{\alpha}\equiv\frac{\alpha-1}{2\alpha}$, a type IV singularity will occur if $0<\tilde{\alpha}<\frac{1}{2}$ (or equivalently, $\alpha>1$).

Observing the equation-of-state parameter in (2.131) we deduce that it can be seen either as a homogeneous one, of the form (2.125) with $G(H)=0$ and

[TABLE]

or as an inhomogeneous one, of the form (2.125) with $f(\rho)=0$ and

[TABLE]

(since $\rho=\frac{3}{\kappa}H^{2}$).

Let us now consider a more general Hubble function inside the class (2.127), namely

[TABLE]

where $c_{0}$, $f_{0}$ are constants, and $\alpha$,$\beta>1$. Thus, two type IV singularities appear at $t=t_{1}$ and $t=t_{2}$. We choose $t_{1}$ to correspond to the inflation end and $t_{2}$ to lie at late times. In order to simplify the expressions, we focus our analysis in the vicinity of the type IV singularity. In this region, inserting (2.134) into the two Friedmann equations leads to

[TABLE]

and therefore the equation of state reads

[TABLE]

Interestingly enough, we observe that the late-time type IV singularity is related to the early-time type IV singularity and the corresponding equation-of-state paramater. In the same lines, the early-time singularity is related to the effective equation of state that gives rise to the late time one.

One can proceed in similar lines, and study the scenario where

[TABLE]

In this case, in the vicinity of the early-time singularity at $t_{1}$ we obtain [54]

[TABLE]

where these relations are again determined by the late-time singularity.

Finally, one can study the scenario where

[TABLE]

which is reproduced by

[TABLE]

At both type IV singularities at $t_{1}$ and $t_{2}$, the effective energy density and pressure become

[TABLE]

and thus the corresponding equation-of-state parameter becomes $-1$.

Let us now proceed to the calculation of the slow-roll parameters, which as usual are used for the calculation of the various inflationary observables, since the effect of the type IV singularity can be significant. The starting point is that in flat FRW geometry one can express the various quantities as a function of the number of $e$-foldings $N$, namely [54]

[TABLE]

where $H^{\prime}(N)=\mathrm{d}H/\mathrm{d}N$. Assuming that the equation of state is given by the general ansatz:

[TABLE]

then (2.144) gives

[TABLE]

Since the usual conservation equation is valid, namely

[TABLE]

with $\rho^{\prime}(N)=\mathrm{d}\tilde{f}(\rho(N))/\mathrm{d}N$, using (2.146) we find

[TABLE]

Finally, inserting (2.148) into (2.145) we acquire

[TABLE]

with $\tilde{f}^{\prime}(\rho(N)\equiv\mathrm{d}\tilde{f}(\rho)/\mathrm{d}\rho$.

Now, for a given $H(t)$, the slow-roll parameters $\epsilon$, $\eta$ and $\xi$ write as [54]

[TABLE]

Hence, if $H(t)$ is given by (2.127), and if $\alpha>1$, i.e when a type IV singularity is obtained at $t\sim t_{s}$ , the slow-roll parameters at the vicinity of the singularity become

[TABLE]

Therefore, we deduce that if $f_{1}(t)$ is a smooth function then the slow-roll parameter $\epsilon$ diverges when $2>\alpha>1$, whereas it remains regular for $\alpha>2$. Moreover, $\eta$ diverges for $3>\alpha>1$. Finally, $\xi^{2}$ diverges when $2>\alpha>1$ and when $4>\alpha>2$. In summary, when $\alpha>4$ all slow-roll parameters are non-singular near the Type IV singularity.

In order to provide a more concrete example, we consider the simplified case where $H(t)=f_{0}\left(t-t_{s}\right)^{\alpha}$. Thus, the slow-roll parameters become [54]

[TABLE]

[TABLE]

[TABLE]

Hence, we immediately observe that the slow-roll parameters exhibit singularities at $t=t_{s}$, as mentioned above. However, such singularities in the slow-roll parameters can be viewed as rather unwanted features.

We close this subsection by making a comparison with observations. In order to achieve this it proves convenient to express the various quantities in terms of the number of $e$-foldings $N$. In particular, for a given $H(N)$ the slow-roll parameters read as [54]

[TABLE]

[TABLE]

[TABLE]

Therefore, we can use (2.143) and (2.144) in order to calculate the usual inflationary observables, namely the spectral index $n_{s}$, the tensor-to-scalar ratio $r$ and the running spectral index $a_{s}$ as [54]

[TABLE]

where

[TABLE]

and

[TABLE]

In order to obtain a qualitative picture of the above observables we consider an equation of state of the form $p=-\rho+f(\rho)$, with $f(\rho)=A\rho^{\alpha}$. Since according to (2.143) and (2.144) the scale factor reads as

[TABLE]

we can express $\rho$ as a function of $N$ as

[TABLE]

Additionally, since a Type IV singularity is obtained when $0<\alpha<\frac{1}{2}$, we can choose $\frac{f(\rho)}{\rho}\ll 1$. Inserting these into (2.170)-(2.172) we finally acquire

[TABLE]

Now, the 2015 Planck results [5] provide the following values:

[TABLE]

Hence, if in the scenario at hand we choose $(N,\alpha)=(60,1/20)$, we obtain

[TABLE]

Therefore, concerning the spectral index we acquire a good agreement, however the values of the tensor-to-scalar ration and of the running spectral index are not inside the observational bounds. Nevertheless, we can obtain satisfactory agreement in more sophisticated models instead of the simple example $f(\rho)=A\rho^{\alpha}$.

3 Late-time acceleration

According to the concordance model of cosmology the universe is currently accelerating, while it entered this era after being in a long matter-dominated epoch. This behavior, similarly to the early accelerated era of inflation, cannot be reproduced within the standard framework of general relativity and flat $\Lambda$CDM-model with dust and vacuum energy, and therefore extra degrees of freedom should be introduced. One can attribute these extra degrees of freedom to new, exotic forms of matter, such as the inflaton field at early times and/or the dark energy concept at late times (for reviews see [55, 56]). Alternatively, one can consider the extra degrees of freedom to have a gravitational origin, i.e. to arise from a gravitational modification that possesses general relativity as a particular limit (see [57, 24, 58, 59, 60] and references therein). In this section we will show how the late-time acceleration can be driven by the fluid viscosity [15, 10, 63, 61, 62, 12, 74, 64, 65, 66, 67, 68, 69, 73, 70, 42, 71, 72].

3.1 Late-time viscous cosmology

We start our investigation by studying the basic scenario of late-time viscous cosmology, presenting the main properties of the viscous cosmic fluid, following [9]. As usual we assume a homogeneous and isotropic FRW universe with geodesic fluid flow, and thus the two Friedmann equations are given by (1.12),(1.13).

Let us very briefly discuss the non-viscous case. According to the standard model the total energy density and pressure are

[TABLE]

where $\rho_{\Lambda}=\Lambda/8\pi G$ is the Lorentz invariant vacuum energy density and $p_{\Lambda}=-\Lambda/8\pi G$ is the vacuum pressure corresponding to a positive tensile stress. With the critical energy density $\rho_{c}$, the matter density parameter $\Omega_{M}$, and the Einstein gravitational constant $\kappa$ defined as

[TABLE]

we obtain for the scale factor [75]

[TABLE]

with $t_{\Lambda}=\frac{2}{3H_{0}\sqrt{\Omega_{\Lambda 0}}},\quad K_{s}=\frac{1-\Omega_{\Lambda 0}}{\Omega_{\Lambda 0}}$, where the subscript zero refers to the present time $t=t_{0}$ (as usual we impose $a(t_{0})=1$). The present age of the universe is

[TABLE]

which leads to $t_{\Lambda}=11.4\times 10^{9}$ years if we insert that $t_{0}=13.7\times 10^{9}$ years and $\Omega_{\Lambda 0}=0.7$. In terms of these quantities for the Hubble parameter we obtain

[TABLE]

whereas the deceleration parameter becomes

[TABLE]

Inserting Eq. (3.4) the present value for for the deceleration parameter of the $\Lambda$CDM- universe is

[TABLE]

With $\Omega_{\Lambda 0}=0.7$ we obtain $\hat{q}_{0}=-0.55$.

It is of interest to determine the time $t=t_{1}$ when deceleration turns into acceleration. The condition for this is $q(t_{1})=0$, and leads to $t_{1}=t_{\Lambda}\rm{arctanh}\frac{1}{\sqrt{3}}$, with corresponding redshift

[TABLE]

that is $t_{1}=7.4\times 10^{9}~{}$years and $z_{1}=0.67$. Finally, let $t_{e}$ be the time of emission of a signal that arrives at the time $t_{0}$. Considering time in units of Gyr and inserting $t_{0}=13.7$ and $\Omega_{\Lambda 0}=0.7$, we acquire the useful expression

[TABLE]

After this brief introduction we now proceed to the investigation of the viscous case, that is we switch on the viscosity in the equation-of-state parameter of the cosmic fluid. For convenience we assume a flat geometry. Without loss of generality we consider the simplest ansatz (1.16), and thus the two Friedmann equations (1.12),(1.13) become

[TABLE]

[TABLE]

while the conservation equation reads

[TABLE]

with

[TABLE]

where in its simplest version $w$ is a constant. Finally, similarly to $\Omega_{M}=\rho/\rho_{c}$, it proves convenient to introduce the density parameters

[TABLE]

where the critical density follows from $3H^{2}=\kappa\rho_{c}$. Thus, we can express the current deceleration parameter as

[TABLE]

If the cosmic fluid is cold, i.e. with $w=0$, as is often assumed, we obtain

[TABLE]

In principle, this equation enables one to estimate the viscosity parameter $\Omega_{\zeta 0}$ if one has at hand accurate measured values of $q_{0}$ and $\Omega_{\Lambda 0}$. It follows from Eqs. (3.7) and (3.16) that

[TABLE]

Hence $\Omega_{\zeta 0}$ is proportional to the deviation of the measured deceleration parameter from the standard $\Lambda$CDM-value as given in Eq. (3.7). This means that one needs to measure the deceleration parameter very accurately in order to obtain information about the viscosity coefficient from its relation to the deceleration parameter. One has so far not been able to determine $\Omega_{\zeta 0}$ in this way.

However, we can indicate its present status. Ten years ago D. Rapetti et al. [76] gave kinematical constraints on the deceleration parameter using type Ia supernovae- and X-ray cluster gas mass fraction measurements, obtaining $q_{0}=-0.81\pm 0.14$ at the 1$\sigma$ confidence level. Inserting $q_{0}>-0.95$ in Eq. (3.17) we obtain $\Omega_{\zeta 0}<0.27$. However, some years later Giostri et al. [77] used SN Ia and BAO/CMB measurements and found $-0.42<q_{0}<-0.20$ with one light curve fitted, and $-0.66<q_{0}<-0.36$ with another. Note that if measurements give $q_{0}<-0.55$ then $\Omega_{\zeta 0}<0$, which is unphysical.

We mention though an interesting study of Mathews et al. [78], in which the production of viscosity was associated with the decay of dark matter particles into relativistic particles in a recent epoch with redshift $z<1$.

Let us review the simplest viscous model in some detail. It was proposed by Padmanabhan and Chitre already in 1987 [79], and is based upon a dust model for matter, vanishing cosmological constant, and constant viscosity coefficient $\zeta=\zeta_{0}$. Equation (3.11) gives

[TABLE]

which upon integration with $H(t_{0})=H_{0}$ leads to

[TABLE]

Another integration with $a(t_{0})=1$ gives

[TABLE]

This implies that the age of the universe when expressed in terms of the present Hubble parameter $H_{0}$ becomes

[TABLE]

Hence, it is seen that for early times, in which $\Omega_{\zeta 0}H_{0}t\ll 1$, the viscosity can be neglected, and we obtain

[TABLE]

corresponding to the evolution of a dust universe. At late times $\Omega_{\zeta 0}H_{0}t\gg 1$, the expansion becomes exponential with $H=\kappa\zeta_{0},\,a\propto\exp(\kappa\zeta_{0}),\,\rho=3\kappa\zeta_{0}^{2}$, and thus the universe enters into a late inflationary era with accelerated expansion. A drawback of this model is however that the time when the bulk viscosity becomes dominant is predicted to be unrealistically large.

Let us now consider briefly the following model, which has attracted attention, namely the one where viscosity is considered to be [80, 81, 66]

[TABLE]

It is based on the physical idea that the dynamic state of the fluid influences its viscosity. We then obtain

[TABLE]

where

[TABLE]

Integrating this equation with $a(0)=0,\,a(t_{0})=1$ and assuming $\kappa(\zeta_{1}+\zeta_{2})<1$ and $w\geq 0$, which lead to $b>0$ and $4bd+c^{2}>0$, we obtain

[TABLE]

with $\quad\hat{H}^{2}=\frac{bd}{a^{2}}+\frac{c^{2}}{4a^{2}}$. The age of the universe in this model becomes

[TABLE]

and thus viscosity increases the age of the universe. Hence, assuming that $\kappa\zeta_{0}\ll H_{0}$ the increase of the age due to viscosity is roughly $\Omega_{\zeta 0}^{2}\,t_{0}$.

Let us return to the solution (3.26) and apply it to the case where the universe does not contain any matter but only dark energy with $w=-1$. Moreover, we assume a linear viscosity $(\zeta_{1}=\zeta_{2}=0)$ and therefore $b=0$. Cataldo et al. [82] found that in this case

[TABLE]

and thus integration with $a(t_{0})=1$ gives

[TABLE]

[TABLE]

Hence, a universe dominated by viscous dark energy with constant viscosity coefficient expands exponentially faster comparing to the corresponding universe without viscosity.

One may now ask the question how does the introduction of a bulk viscosity confront with the observed acceleration of the universe. There have been several works dealing with this issue, for instance see Refs. [83, 84, 85]. In the model of Avelino and Nucamendi [85] it was considered that $\zeta_{1}=\zeta_{2}=0$, $w=0$, $\Omega_{M}=1$, $\Omega_{\Lambda}=0$, and therefore the scale factor can be written as

[TABLE]

which satisfies the boundary conditions $a(0)=0,a(t_{0})=1$. The age of the universe in this model becomes

[TABLE]

Such a universe model was actually considered earlier, by Brevik and Gorbunova [87, 86] and by Grøn [88], and is also similar to the model of Padmanabhan and Chitre considered above [79]. The Hubble parameter reads as

[TABLE]

and it approaches a de Sitter phase for $t\gg 1/\Omega_{\zeta 0}H_{0}$, with a constant Hubble parameter equal to $\Omega_{\zeta 0}H_{0}$. The deceleration parameter is

[TABLE]

and its value at present is

[TABLE]

Hence

[TABLE]

Assuming that accurate measurements will verify the $\Lambda$CDM model, so that $q_{0}=-0.55$, this equation implies that $\Omega_{\zeta 0}=0.7$. This means that for the universe model to be realistic, there must exist a physical mechanism able to produce a viscosity of this magnitude.

The expansion thus starts from a Big Bang with an infinitely large velocity, but decelerates to a finite value. As usual, when $t=t_{1}$ determined by $q(t_{1})=0$ there is a transition to an accelerated eternal expansion, namely at

[TABLE]

at which time the scale factor is

[TABLE]

and the corresponding redshift is

[TABLE]

Under the assumption that this model contains a mechanism producing viscosity so that $\Omega_{\zeta 0}=0.7$, this equation gives $z_{1}=0.8$. This is larger than the corresponding value in the $\Lambda$CDM model. Hence the transition to accelerated expansion happens earlier if the acceleration of the expansion is driven by viscosity than by dark energy.

We deduce that the bulk viscosity must have been sufficiently large, namely $\Omega_{\zeta 0}>1/3$, in order for this transition to have been realized in the past, i.e. at $a(t_{1})<1$. Finally, note that for this universe model, with spatial curvature $k=0$, the matter density is equal to the critical density, namely

[TABLE]

and thus the matter density approaches a constant value, $\rho\rightarrow(3/\kappa)\Omega_{\zeta 0}^{2}H_{0}^{2}$.

In the aforementioned study of Avelino and Nucamendi [85] supernova data were used in order to estimate the value of $\Omega_{\zeta 0}$, giving the best fit for a universe containing dust with constant viscosity coefficient. The result was that $\Omega_{\zeta 0}=0.64$ had to be several orders of magnitude greater than estimates based upon kinetic gas theory [10]. However, as an unorthodox idea we may mention here the probability for producing larger viscosity via dark matter particles decaying into relativistic products [89]. Additionally, the comparison between the magnitude of bulk viscosity and astronomical observations were also performed in a recent paper by Normann and Brevik [71], using the analyses of various experimentally-based sources [90, 91]. Various ansatzes for the bulk viscosity were analyzed: (i) $\zeta=$constant, (ii) $\zeta\propto\sqrt{\rho}$, and (iii) $\zeta\propto\rho$. The differences between the predictions of the options were found to be small. As a simple estimate based upon this analysis, we suggest that

[TABLE]

can serve as a reasonable mean estimate for the present viscosity. With $H_{0}=67.7\,\rm{km}\,\rm{s}^{-1}\,\rm{Mpc}^{-1}$ this corresponds to $\Omega_{\zeta 0}=0.01$.

The behavior of a viscous universe in its final stages has been discussed in [87, 86] and in [82, 92]. Consider first a universe without viscosity and dark energy, containing only a non-viscous fluid with $p=w\rho$. In this case Eq. (3.24) reduces to

[TABLE]

where $b=\frac{3}{2}(1+w)$. For such a universe there is a Big Rip at

[TABLE]

On the other hand, in Ref. [82] a fluid was considered to have $w<-1$ and constant viscosity coefficient $\xi_{0}$, implying $b<0$ and $d=0$. In this case (3.24) reduces to

[TABLE]

The Hubble parameter, scale factor and density for this universe are respectively extracted to be

[TABLE]

[TABLE]

and

[TABLE]

Thus, in this case there is a Big Rip singularity at

[TABLE]

Similar models, with variable gravitational and cosmological “constants” have been investigated by Singh et al. [93, 94]. Furthermore, one can go beyond isotropic geometry and study viscous fluids in spatially anisotropic spaces, belonging to the Bianchi type-I class. The interested reader might consult, for instance, the discussion in Ref. [9].

3.2 Inhomogeneous equation of state of the universe: phantom era and

singularities

In this subsection we examine the appearance of singularities in viscous cosmology. It is well-known that in FRW geometry, when the equation of state modeling the matter content is a linear equation with an equation of state parameter greater than $-1$, the Big Bang singularity appears at early times, where the energy density of the universe diverges. Moreover, dealing with nonlinear equations of state one can see that other kind of singularities such as Sudden singularity [95, 96, 97] or Big Freeze [51, 98, 99, 100] appear.

In fact, the future singularities are classified as follows [51] (see also [101] for a more detailed classification):

• •

Type I (Big Rip): $t\to t_{s}$, $a\to\infty$, $\rho\to\infty$ and $|p|\to\infty$.

• •

Type II (Sudden): $t\to t_{s}$, $a\to a_{s}$, $\rho\to\rho_{s}$ and $|p|\to\infty$.

• •

Type III (Big Freeze): $t\to t_{s}$, $a\to a_{s}$, $\rho\to\infty$ and $|p|\to\infty$.

• •

Type IV (Generalized Sudden): $t\to t_{s}$, $a\to a_{s}$, $\rho\to 0$, $|p|\to 0$ and derivatives of $H$ diverge.

Similarly to the future ones, one can define the past singularities:

• •

Type I (Big Bang): $t\to t_{s}$, $a\to 0$, $\rho\to\infty$ and $|p|\to\infty$.

• •

Type II (Past Sudden): $t\to t_{s}$, $a\to a_{s}$, $\rho\to\rho_{s}$ and $|p|\to\infty$.

• •

Type III (Big Hottest): $t\to t_{s}$, $a\to a_{s}$, $\rho\to\infty$ and $|p|\to\infty$.

• •

Type IV (Generalized past Sudden): $t\to t_{s}$, $a\to a_{s}$, $\rho\to 0$, $|p|\to 0$ and derivatives of $H$ diverge.

For the simple case of a linear equation of state $p=w\rho$ it is well-known that for a non-phantom fluid ($w>-1)$ one obtains a Big Bang singularity, while for a phantom fluid ($w<-1)$ [102, 103, 104, 105, 106] the singularity is a future Type I (Big Rip). Hence, in order to obtain the other type of singularities one has to consider phantom fluids modeled by non-linear equations of state of the form

[TABLE]

where $f$ is a positive function. The simplest model is obtained taking $f(\rho)=A\rho^{\alpha}$ with $A>0$. In this case from the conservation equation $\dot{\rho}=-3H(\rho+p)$ and the Friedmann equation $H^{2}=\frac{\kappa\rho}{3}$ one obtains the dynamical equation

[TABLE]

whose solution is

[TABLE]

Furthermore, in order to obtain the evolution of the scale factor we will integrate the conservation equation, resulting in

[TABLE]

which using (3.53) leads to

[TABLE]

Once we have calculated these quantities, we have the following different situations (see also [51]):

When $\alpha<0$ we have a past singularity of Type II, since the energy density vanishes for $t_{s}=t_{0}-\frac{2}{\sqrt{3\kappa}A}\frac{\rho_{0}^{\frac{1}{2}-\alpha}}{1-2\alpha}<t_{0}$, implying that the pressure diverges at $t=t_{s}$. 2. 2.

When $\alpha=0$ there are no singularities. The dynamics is defined from $t_{s}=t_{0}-\frac{2}{\sqrt{3\kappa}A}\sqrt{\rho_{0}}$ (where the energy density is zero) up to $t\rightarrow\infty$. 3. 3.

When $0<\alpha<\frac{1}{2}$ there are two different cases:

1. (a)

$\frac{1}{1-2\alpha}$ is not a natural number. One has a past Type IV singularity at $t_{s}=t_{0}-\frac{2}{\sqrt{3\kappa}A}\frac{\rho_{0}^{\frac{1}{2}-\alpha}}{1-2\alpha}$, since higher derivatives of $H$ diverge at $t=t_{s}$. 2. (b)

$\frac{1}{1-2\alpha}$ is a natural number. In that case there are not any singularites and the dynamics is defined from $t_{s}=t_{0}-\frac{2}{\sqrt{3\kappa}A}\frac{\rho_{0}^{\frac{1}{2}-\alpha}}{1-2\alpha}$ to $t\rightarrow\infty$. 4. 4.

When $\alpha=\frac{1}{2}$ there are no singularities in cosmic time. 5. 5.

When $\frac{1}{2}<\alpha<1$, one has future Type I singularities, since in this case $\rho$, $p$ and $a$ diverge at $t_{s}=t_{0}-\frac{2}{\sqrt{3\kappa}A}\frac{\rho_{0}^{\frac{1}{2}-\alpha}}{1-2\alpha}>t_{0}$. 6. 6.

When $\alpha=1$ the equation of state is linear, and thus we obtain a Big Rip singularity. 7. 7.

When $\alpha>1$, the energy density and the pressure diverge but the scale factor remains finite at $t=t_{s}$, implying that we have a future Type III singularity.

The remarkable case appears when $0<\alpha<\frac{1}{2}$ and with $\frac{1}{1-2\alpha}$ being a natural number. In this case, from the Friedmann equation $H^{2}=\frac{\kappa\rho}{3}$ and the solution (3.53) one obtains

[TABLE]

with $n=\frac{1}{1-2\alpha}$. As we have already seen, this solution describes a universe in the expanding phase driven by a phantom fluid, which is defined from $t_{s}=t_{0}-\frac{2}{\sqrt{3\kappa}A}\frac{\rho_{0}^{\frac{1}{2}-\alpha}}{1-2\alpha}$ (where $H=0$) up to $t\rightarrow\infty$. However, solution (3.58) could be extended analytically back in time. There are two different cases: When $n$ is odd, this extended solution describes a universe driven by a phantom field that goes from the contracting to expanding phase, bouncing at time $t_{s}$. On the contrary, when $n$ is even the universe moves always in the expanding phase, and before $t_{s}$ it is driven by a non-phantom field, while after $t_{s}$ the universe enters in a phantom era. We will explain this phenomenon in more detail in the next subsection.

Motivated by the introduction of bulk viscous terms in an ideal fluid one can consider a subclass of the general equation of state of (1.15) of the form

[TABLE]

Then, the conservation equation becomes $\dot{\rho}=3H[f(\rho)-G(H)]$, and using the Friedmann equation (1.12) in the expanding phase leads to

[TABLE]

which reveals that this formalism is equivalent with considering a fluid with an effective equation of state given by

[TABLE]

It is clear that, in general, the equation of state (3.59) does not lead to a universe crossing the phantom barrier. A simple way to obtain transitions from the non-phantom to the phantom regime is to explicitly consider an inhomogeneous equation of state of the form $F(\rho,p,H)=0$, for example [12, 74]

[TABLE]

with $C_{0}$ and $H_{0}$ some positive constants. Inserting this into the square of the equation $\dot{H}=-\frac{\kappa}{2}(\rho+p)$, one obtains the bi-valued dynamical equation

[TABLE]

From this equation, since there are two square roots and the effective equation of state parameter is given by $w_{eff}\equiv-1-\frac{2\dot{H}}{3H^{2}}$, one can see that there are two different dynamics: one which corresponds to the branch with $\dot{H}<0$ describing a universe in a non-phantom regime, and one corresponding to the branch $\dot{H}>0$ describing a universe in the phantom era. In fact, (3.63) can be integrated as

[TABLE]

where we have introduced the notation $t_{\pm}=\pm\frac{4}{3C_{0}H_{0}}$. It is easy to check that $H(t)$ is only defined for times between $t_{-}$ and $t_{+}$, since at $t_{\pm}$ the Hubble function $H$ diverges (we obtain a Big Bang at $t_{-}$ and a Big Rip at $t_{+}$). Moreover, it is a decreasing function for $t\in(t_{-},0)$ and an increasing one for $t\in(0,t_{+})$, implying that at $t=0$ the universe crosses the phantom divide (it passes from the non-phantom to the phantom era).

Another interesting example arises from the equation of state

[TABLE]

where $t_{0},H_{0},H_{1}$ are parameters satisfying $H_{0}>H_{1}>0$. The corresponding bi-valued dynamical equation is

[TABLE]

which has two fixed points, namely $H_{0}$ and $H_{0}-H_{1}$. As we have already explained, when $\dot{H}<0$ (resp. $\dot{H}>0$ ) the universe is in a non-phantom (resp. phantom) era. When the universe is in the branch with $\dot{H}<0$ it moves from $H_{0}$ to $H_{0}-H_{1}$, it reaches $H=H_{0}$ and then it enters in the other branch ($\dot{H}>0$) going from $H_{0}-H_{1}$ to $H_{0}$. In fact, in [12, 74] the authors found the following solution:

[TABLE]

which satisfy all the properties described above.

A final remark is in order: One can indeed consider the more general equation of state given in (1.15), namely of the form $F(\rho,p,H,\dot{H},\ddot{H},\cdots)=0$, containing higher order derivatives of the Hubble parameter. In this case, using the Friedmann equations the equation of state becomes the dynamical equation

[TABLE]

A non-trivial example is the following equation of state [12, 74]:

[TABLE]

where $G_{0}$ and $G_{1}$ are constant. Then, the dynamical equation becomes

[TABLE]

We look for periodic solutions of the form $H(t)=H_{0}\cos(\Omega t)$ depicting an oscillatory universe. Inserting this expression into (3.70) we obtain the algebraic system:

[TABLE]

whose solution is given by

[TABLE]

provided that $G_{0}(1+w)>0$ and $G_{1}(1+w)>0$. On the other hand, when $G_{1}(1+w)<0$, one can look for solutions of the form $H(t)=H_{0}\cosh(\Omega t)$, obtaining

[TABLE]

3.3 Unification of inflation with dark energy in viscous cosmology

In this subsection we analyze how one can describe in a unified way the early (inflationary) and late time acceleration, in the framework of viscous cosmology. The simplest way to unify early inflationary epoch with the current cosmic acceleration is by using scalar fields [107]. Starting with the action

[TABLE]

where $\omega$ and $V$ are functions of the scalar field $\phi$, and focusing on flat FRW geometry, one obtains the following dynamical equation

[TABLE]

The relevant fact, is that given a function $f(\phi)$ the equation (3.75) has always the solution $\phi=t$ and $H=f(t)$, provided that (for details see [12, 74])

[TABLE]

An interesting example is obtained when one considers the function

[TABLE]

where $H_{0}$ and $\phi_{s}$ are the two positive parameters of the model. In this case one has

[TABLE]

whose dynamics is given by

[TABLE]

where we have introduced the notation $t_{s}=\phi_{s}$. Since $H$ diverges at $t=0$ and $t=t_{s}$, the dynamics is defined in $(0,t_{s})$. In fact at $t=0$ one has $a=0$, which means that we obtain a Big Bang singularity, while at $t=t_{s}$ the scale factor diverges, implying that we have a Big Rip singularity. On the other hand, the derivative of the Hubble parameter reads as

[TABLE]

that is the universe lies in the non-phantom regime when $0<t<t_{s}/2$, while it lies in the phantom phase for $t_{s}/2<t<t_{s}$. Hence, we conclude that this model could describe the current cosmic acceleration.

In order to examine the behavior at early times, we note that near $t=0$ one can make the approximation $a=a_{0}\left(\frac{t}{t_{s}}\right)^{H_{0}t_{s}}$, and thus its second derivative at early times is approximately

[TABLE]

From this we deduce that if one chooses $H_{0}t_{s}>1$ the universe will have an early period of acceleration.

Another example is to consider

[TABLE]

with $H_{0}$ and $\nu$ positive parameters. A straightforward calculation leads to

[TABLE]

In this case one obtains a non-singular oscillating universe, whose dynamics is given by

[TABLE]

This solution depicts a universe that bounces at time $t=\frac{n\pi}{\nu}$ where $n$ is an integer, and since $\dot{H}=H_{0}\nu\cos(\nu t)$ one can easily check that the universe lies in the phantom regime when $\frac{\pi}{\nu}\left(-\frac{1}{2}+2n\right)<t<\frac{\pi}{\nu}\left(\frac{1}{2}+2n\right)$, while it is in the non-phantom phase when $\frac{\pi}{\nu}\left(\frac{1}{2}+2n\right)<t<\frac{\pi}{\nu}\left(\frac{3}{2}+2n\right)$.

We proceed by considering viscosity, taking $\Lambda=0$ and $w=1$ in (3.11). Based on the equivalence between bulk viscous and open cosmology, where isentropic particle production is allowed [108], we choose the following viscosity coefficient [109]:

[TABLE]

where $\xi_{0}>0$ is a constant. Hence, the second Friedmann equation (3.11) becomes

[TABLE]

which only has $H=\frac{\xi_{0}}{8}$ as a fixed point. If one considers the dynamics in the domain $\frac{\xi_{0}}{8}\leq H\leq\infty$, it is easy to check that the effective equation of state parameter is greater than $-1$, which implies that the Hubble parameter varies from infinity to $\frac{\xi_{0}}{8}$. Moreover, since

[TABLE]

$w_{eff}\cong-1$ at early ($H\gg\xi_{0}$) and late ($H\cong\frac{\xi_{0}}{8}$) times, from which we deduce that this viscous fluid model unifies inflation with the current cosmic acceleration. Additionally, $w_{eff}$ is positive when $\xi_{0}\frac{\sqrt{2}-1}{2\sqrt{2}}<H<\frac{\sqrt{2}+1}{2\sqrt{2}}$, having the maximum value $w_{eff}=1$ at $H=\frac{\xi_{0}}{4}$. Thus, in summary, in the scenario at hand the universe starts from an inflationary epoch, it evolves through a Zel’dovich fluid ($w_{eff}=1$), radiation- ($w_{eff}=1/3$) and matter- ($w_{eff}=0$) dominated epochs, and finally it enters into late-time acceleration tending towards a de Sitter phase.

The solution of equation (3.90) is

[TABLE]

and the scalar field that induces this dynamics, if one chooses $\omega(\phi)\equiv 1$, has the following Higgs-style potential (for details see [109]):

[TABLE]

We stress here that this unified model for inflation and late-time acceleration leads to inflationary observables, namely the spectral index, its running and the ratio of tensor to scalar perturbations, that match at $2\sigma$ Confidence Level with the observational data provided by Planck 2015 announcements [5] (for a detailed discussion see [109, 110]).

We close this section by considering a very simple quintessential-inflation potential which unifies inflation with late time acceleration, namely [111]

[TABLE]

where $\phi_{E}\equiv-\sqrt{\frac{2}{3\kappa}}\frac{H_{E}}{\sqrt{H_{E}^{2}-\frac{\Lambda}{3}}}$, and with $H_{E}>0$ the parameter of the model. This model leads to the following dynamics

[TABLE]

whose solution has the following expression

[TABLE]

with the corresponding scale factor

[TABLE]

We mention that this dynamics arises also from a universe filled with a fluid with the simple linear equation of state of the form

[TABLE]

where $\rho_{E}=\frac{3H_{E}^{2}}{\kappa}$. Equivalently it can arise from a viscous fluid, since effectively, choosing $w=1$ and the following viscosity coefficient

[TABLE]

and inserting it into (3.11) one obtains the dynamics (3.99). Finally, for this model the effective equation-of-state parameter is given by

[TABLE]

which shows that for $H\gg H_{E}$ one has $w_{eff}(H)\cong-1$ (early quasi - de Sitter period). When $H\cong H_{E}$, the equation-of-state parameter satisfies $w_{eff}(H)\cong 1$ (deflationary period dominated by a Zel’dovich fluid), and lastly for $H\cong\sqrt{\frac{\Lambda}{3}}$ one also acquires $w_{eff}(H)\cong-1$ (late quasi - de Sitter period).

3.4 Generalized holographic dark energy with a viscous fluid

In this subsection we investigate the cosmological scenario of holographic dark energy with the presence of a viscous fluid. Holographic dark energy [112, 113] is a scenario in the direction of incorporating the nature of dark energy using some basic quantum gravitational principles. It is based on black hole thermodynamics [114] and the connection of the ultraviolet cut-of of a quantum field theory, which induces the vacuum energy, with the largest distance of this theory [115]. Determining suitably an IR cut-off $L$, and imposing that the total vacuum energy in the maximum volume cannot be greater than the mass of a black hole of the same size, one obtains the holographic dark energy, namely

[TABLE]

with $\kappa$ the gravitational constant, set to $\kappa=1$ in the following for simplicity, and $c$ a parameter. The holographic dark energy scenario has interesting cosmological applications [74, 116, 117, 118, 119]. Concerning the ultraviolet cut-off one uses the future event horizon $L_{f}$

[TABLE]

However, one can generalize the model using the quadratic Nojiri-Odintsov cut-off $L$ defined as [74, 120, 121]

[TABLE]

with $c$, $\alpha_{0}$, $\alpha_{1}$ and $\alpha_{2}$ constants, or the generalized Nojiri-Odintsov cut-off defined in Refs. [74] and [122].

In this subsection we will consider the scenario in which generalized holographic dark energy interacts with a viscous fluid, following [123]. In particular, we consider a dark matter sector with a viscous equation of state of the form

[TABLE]

where $\rho_{DM}$ and $p_{DM}$ are respectively the dark matter energy density and pressure, and with $\alpha$, $\chi$ and $\beta$ the model parameters. Furthermore, we allow for an interaction between viscous dark matter and holographic dark energy:

[TABLE]

[TABLE]

with $w_{DE}$ and $w_{DM}$ respectively the equation-of-state parameters of the dark energy and dark matter sectors, and where $Q$ is a function that determines the interaction. One can impose the following form for $Q$ [123]

[TABLE]

where $b$ is a constant, although more complicated forms could also be used [124]. Finally, the first Friedmann equation reads as

[TABLE]

where the effective (total) energy density is given by $\rho_{eff}=\rho_{DE}+\rho_{DM}$.

We start our analysis by investigating the non-interacting scenario, that is setting $Q$ (i.e. $b$) to zero. In this case, using (3.115),(3.116),(3.118) and (3.122), the deceleration parameter $q\equiv-1-\dot{H}/H^{2}$ is found to be

[TABLE]

where $\hat{q}_{0}=H\Omega_{de}(\dot{L}_{f}+1)+L_{f}\left(H^{2}+p_{DM}\right)$. Moreover, the evolution of the dark matter density parameter $\Omega_{DM}$ is determined by the differential equation

[TABLE]

with $\hat{L}_{f}=\dot{L}_{f}(\alpha_{1}+2\alpha_{2}L_{f})$ and $\hat{A}_{0}=\Omega_{DE}[H(\dot{L}_{f}+1)+L_{f}p_{DM}]+H\Omega_{DE}^{2}(\dot{L}_{f}+1)$, and where primes denote differentiation with respect to $N=\ln a$. Finally, for the non-interacting case (3.115) and (3.119) lead to

[TABLE]

As one can see, the deceleration parameter $q$ starts from positive values, it decreases, and it becomes negative marking the passage to late-time accelerated phase [123]. The role of the viscosity parameter $\chi$ is significant, since larger positive $\chi$ leads the transition redshift $z_{tr}$ (from deceleration to acceleration) to smaller values and the present deceleration parameter to negative values closer to zero. Hence, the larger the fluid viscosity is the more difficult it is for the universe to exhibit accelerated expansion. Lastly, an important feature is that the dark energy equation of state parameter can exhibit the phantom divide-crossing, as can be seen from (3.125) [74, 123, 125]. In Tables 2 and 3 we present the various calculated values for different choices of the model parameters, where the features described above are obvious.

Let us now study the interacting scenario, i.e. considering a non-zero $Q$ in (3.119),(3.120). In this case, using (3.115),(3.116),(3.118) and (3.122), for the deceleration and matter density parameters we find

[TABLE]

and

[TABLE]

with $A_{1}=\Omega_{DE}\left\{H[(3b-1)HL_{f}-2]-L_{f}p_{DM}\right\}$, while for the dark-energy equation-of-state parameter we obtain

[TABLE]

As we observe, the deceleration parameter $q$ exhibits the transition from deceleration to acceleration, and the role of the positive interaction parameter $b$ in (3.121) is to make $z_{tr}$ larger and the present value of $q$ more negative [123]. This is expected since larger positive $b$ implies larger positive $Q$ in (3.119),(3.120) and thus larger energy transfer to the dark energy sector. Moreover, the role of the viscosity parameter $\chi$ is as in the non-interacting case, i.e the larger the $\chi$ is the more difficult it is for the universe to exhibit accelerated expansion. Lastly, the dark energy equation-of-state parameter $w_{DE}$ can exhibit the phantom divide-crossing, too. In Table 4 we present the various calculated values for different choices of the model parameters, where the features described above are obvious.

In summary, as we saw, one can study the scenario of generalized holographic dark energy in the framework of viscous cosmology, allowing additionally for an interaction term between viscous dark matter and holographic dark energy. As one can show, the role of viscosity is to make the the transition to late-time acceleration more difficult, while the role of interaction has the opposite effect. Lastly, the scenario at hand allows for the phantom-divide crossing, which can be an additional advantage revealing its capabilities.

4 Special topics

In this section we discuss various topics of viscous cosmological theory, focusing on investigations in which the present authors have taken part. As a brief remark to the material covered below we think it is appropriate to underscore the great power of the hydrodynamical formalism when applied to quite different problems in cosmology. The formalism robustness is in general striking. Definitely, in view of the considerably large activity in the field of viscous cosmology, there are many aspects that cannot be discussed here. For instance, instead of assuming a one-component fluid model, one might consider an extension of the model in order to encompass two different fluid components. We may here mention the recent study of Ref. [126], where the cosmic fluid was considered to be constituted of a dark matter component endowed with a constant bulk viscosity, and a non-viscous dark energy component. In other related works [73, 127], viscous coupled-fluid models were investigated when the equation of state was assumed to be inhomogeneous. Furthermore, in [128] the authors studied the important self-reproduction problem of the universe, namely the graceful exit from inflation, where it was shown how inflation without self-reproduction can actually be obtained by imposing restrictions on the value of the thermodynamic parameter in the equation of state. Finally, we mention the very different approach which consists in applying particle physics theory and the relativistic Boltzmann equation in order to derive expressions for the bulk and the shear viscosities, and the corresponding entropy production, in the specific lepton-photon era, where the temperature dropped from $10^{12}~{}$K to $10^{10}~{}$K. Calculations of this kind were recently given in Ref. [129], [130] and [131].

4.1 Estimate for the present bulk viscosity and remarks on the future

universe

A significant amount of research has been spent in order to study the behavior of the cosmic fluid in the far future. In such an examination, and as we discussed in detail in subsection 3.2 above, there may appear various kinds of singularities: the Big Rip [61, 62], the Little Rip [132, 67, 133], the Pseudo-Rip [134], the Quasi-Rip [135], as well as other kinds of soft singularities (for instance the so-called type IV finite time singularities [136]).

In the framework of viscous cosmology, the value of the (effective) bulk viscosity at present time is naturally an important ingredient of such investigation. Recent observations from the Planck satellite have given us a better ground for estimating the bulk viscosity value $\zeta=\zeta_{0}$ at present time $t=t_{0}$. As discussed already in the previous Sections, referring to [71], as well as to several other theoretical and experimental manuscripts, the estimate

[TABLE]

was suggested as a reasonable (logarithmic) mean value. However, the corresponding uncertainty is quite large; there have appeared proposals ranging from about $10^{4}~{}\rm{Pa~{}s}$ to about $10^{7}~{}\rm{Pa~{}s}$, depending on analyses of different sources.

We will follow the discussion of [126], in which two different cosmological models were analyzed: (1) a one-component dark energy model where the bulk viscosity $\zeta$ was associated with the cosmic fluid as a whole, and (2) a two-component model where $\zeta$ was associated with a dark matter component $\rho_{m}$ only, the latter component assumed to be non-viscous. For convenience, we focus on the one-component scenario.

We assume the simple equation of state $p=w\rho$, with $w=const.$, and hence the two viscous Friedmann equations acquire the usual form, namely

[TABLE]

and the energy conservation equation reads as

[TABLE]

Solving this equation in the regime around $w=-1$, i.e expanding as $w=-1+\alpha$ and assuming that $\alpha$ is small, we obtain

[TABLE]

where $t_{0}=0$, and the integration extends into the future. For the bulk viscosity we will consider the form adopted in the literature, namely

[TABLE]

with $\lambda$ a constant. In the following we examine two options for the value of $\lambda$, which are both physically reasonable.

• •

Case (i): $\lambda=1/2~{}(\zeta\propto\sqrt{\rho})$.

In this case, from Eq. (4.4) we obtain

[TABLE]

where for convenience we have introduced the dimensionless quantities

[TABLE]

The point that worths attention here is that even if the fluid is initially in the quintessence region $\alpha>0$ at $t=0$ it will, if $X_{0}>0$, inevitably be driven into a Big Rip singularity $(\rho\rightarrow\infty)$ after a finite time [87, 86, 69, 71]

[TABLE]

If on the other hand the combination of equation-of-state parameter $\alpha$ and viscosity $\zeta_{0}$ is such that $X_{0}<0$, then the cosmic fluid becomes gradually diluted as $\rho\propto 1/t^{2}$ in the far future.

• •

Case (ii): $\lambda=0~{}(\zeta=\zeta_{0}=const.)$.

In this case we obtain the solution

[TABLE]

which implies an energy density of the form

[TABLE]

Hence in the far future $\rho\rightarrow const.$, which implies $H\rightarrow const.$, which is just the de Sitter solution. Let us denote the limiting value of the density by $\rho_{\rm dS}$. Then

[TABLE]

From this expression we deduce that both $\alpha$ and $X_{0}$ are important for the future fate of the cosmic fluid.

Thus, this case may be defined as a pseudo-Rip in accordance with the definition given by Frampton et al. [134], since the limiting value of the density reached after an infinite span of time is finite.

We close this subsection by providing some values for the inflationary observables, in order to compare with the 2015 Planck observations. In particular, from Table 5 of [5] we have $w=-1.019^{+0.075}_{-0.080}$. Thus, $\alpha=1+w$ will be lying within two limits, i.e. between

[TABLE]

As mentioned above, we took $\zeta_{0}=10^{6}~{}$Pa s, i.e. $\Omega_{\zeta 0}=0.01$, to be a reasonable mean value of the present viscosity. Then, according to (4.7) we have

[TABLE]

Hence, we recover the cases 2 and 3 above: the future de Sitter energy density will become lower than $\rho_{0}$.

4.2 Is the bulk viscosity large enough to permit the phantom divide crossing?

This subsection is a continuation of the previous one, and is motivated by the following question: is the value of $\zeta_{0}$, as inferred from the analysis of recent observations, actually large enough to permit the crossing of the phantom divide, i.e. the transition from the quintessence region to the phantom region? To analyze this question we have to consider more carefully the uncertainties in the data found from different sources. We will present some material discussing this point, following the recent work [69].

Assume that the bulk viscosity varies with energy density as $\zeta\propto\sqrt{\rho}$. The condition for phantom divide crossing, as noted above, is that the quantity $X_{0}$ defined in Eq. (4.7) has to be positive. In the analysis of Wang and Meng [90] various assumptions for the bulk viscosity in the early universe were considered, and the corresponding theoretical curves for $H=H(z)$ were compared with a number of observations. The detailed comparison is quite complicated, but for our purpose it is sufficient to note that the preferred value of the magnitude $\Omega_{\zeta 0}$ is (compare also with the discussion in [71]):

[TABLE]

corresponding to

[TABLE]

which is a rather high value. In this context, we may compare with the formula for the bulk viscosity in a photon fluid [15], namely

[TABLE]

where $a_{\rm rad}=\pi^{2}k_{B}^{4}/15\hbar^{3}c^{3}$ is the radiation constant and $\tau_{f}$ the mean free time. If we estimate $\tau_{f}=1/H_{0}$ (the inverse Hubble radius), we obtain $\zeta\sim 10^{4}~{}$Pa s, which is considerably lower. In summary, it seems that one has to allow for a quite wide span in the value of the present bulk viscosity. All suggestions in the literature can be encompassed if we write

[TABLE]

We can now rewrite the condition for phantom divide crossing as

[TABLE]

where we have inserted $H_{0}=67.80$ km s*-1* Mpc${}^{-1}=2.20\times 10^{-18}~{}{\rm s}^{-1}$. As noted above, from the observed data we derive the maximum value of $\alpha$ to be $\alpha_{\rm max}=0.056$. This yields

[TABLE]

Thus, comparison between (4.17) and (4.19) implies that, on the basis of available data, a phantom divide crossing is actually possible even if $\alpha=\alpha_{\rm max}$.

4.3 Bounce universe with a

viscous fluid

In this subsection we investigate the realization of bouncing solutions in the framework of viscous cosmology following [137] (see also [138]). Bouncing cosmological evolutions offer a solution to the initial singularity problem [139]. Such models have been constructed in modified gravity constructions, such as in the Pre-Big-Bang [140] and in the Ekpyrotic [141] scenarios, in $f(R)$ gravity [142, 143, 144], in $f(T)$ gravity [145], in braneworld scenarios [146, 147], in loop quantum cosmology [148, 149] etc. Additionally, non-singular bounces can be obtained using matter forms that violate the null energy condition [150, 151].

In order to be more general, in the following we will allow also for a spatial curvature, and hence the two Friedmann equations write as

[TABLE]

with $k=-1,0,1$ corresponing to open, flat or closed geometry. Additionally, concerning the fluid’s equation of state we will consider a general inhomogeneous viscous one of the form (1.15), namely

[TABLE]

where $w(\rho)$ can depend on the energy density, but the bulk viscosity $B(a(t),H,\dot{H}...)$ is allowed to be a function of the scale factor, and of the Hubble function and its derivatives. Thus, the fluid stress-energy tensor writes as

[TABLE]

with $u_{\mu}=(1,0,0,0)$ the four velocity. Hence, the standard conservation law $\dot{\rho}+3H(\rho+p)=0$ leads to

[TABLE]

We now proceed to the investigation of simple bounce solutions in the above framework, and we discuss the properties of the viscosity of the fluids that drive such solutions. A first example is the bounce with an exponential scale factor of the form

[TABLE]

with $n$ a positive integer and $a_{0}\,,\alpha$ positive parameters. We consider $t_{0}>0$ to be the bounce point, i.e. for $t<t_{0}$ we have a contracting universe and when $t>t_{0}$ expansion takes place. We mention that if $n$ is non-integer then singularities may arise, while the simplest case $n=1/2$ corresponds to just the Sitter solution $H(t)=const.$ (in general for $n=m/2$, with $m$ an odd integer, the bounce is absent). Finally, note that for the ansatz (4.26) we have

[TABLE]

and hence we obtain (early-time) acceleration after the bounce, which is a significant phenomenological advantage.

Inserting the bouncing scale factor (4.26) into (4.20) we acquire

[TABLE]

Since for $k=-1$ the above quantity might become negative, we focus on the $k=0$ and $k=+1$ cases where it is always positive definite. As we observe, in the flat case $\rho$ decreases in the contracting phase, it becomes zero at $t=t_{0}$, and it increases in the expanding regime. On the other hand, for $k=+1$, and when $n>1$, there is a region around the bouncing point where $\rho$ increases in the contracting phase, it reaches the value $\rho=3/(a_{0}\kappa)$ at $t=t_{0}$, and then it decreases (these can be seen by examining the derivatives of (4.28)). However, for $t\gg t_{0}$, the energy density starts to increase. This behavior may have an important effect on the cosmological parameter $\Omega=1+\frac{k}{a^{2}H^{2}}$, which for the bouncing scale factor (4.26) becomes

[TABLE]

and thus it exhibits a decreasing behavior. Such a post-bounce acceleration, with the simultaneous decrease of $\rho$ and of $\Omega$ may be compatible with the inflationary phenomenology, in which at the end of inflation $\Omega$ is very close to $1$. Definitely, in order to stop the aforementioned early-time acceleration we need to add additional fluids that could become dominant and trigger the transition to the matter era.

Let us now analyze what kind of fluids with equation of state given by (4.22) can produce the bouncing solution (4.26). We first consider an inhomogeneous but non-viscous fluid, namely we assume $B(a(t),H,\dot{H}...)=0$. In this case, for the flat geometry, equations (4.20) and (4.22) lead to

[TABLE]

and thus to

[TABLE]

As we mentioned earlier, if the exponent of $\rho$ in (4.31) is negative then we obtain the bounce realization, however if it is positive then we have the appearance of a singularity.

As a second example we switch on viscosity, considering

[TABLE]

with $\zeta(H)>0$ the bulk viscosity, and for simplicity and without loss of generality we consider $w=-1$. In this case, for the flat geometry, equations (4.20) and (4.22) lead to

[TABLE]

Nevertheless, for $k=+1$ the equation of state for the fluid becomes complicated and therefore it is necessary to go beyond (4.32) and consider a viscosity that depends on the scale factor too. Such a case could be

[TABLE]

which then leads to

[TABLE]

As we can see, and as expected, for large scale factors the above relation coincides with (4.33), and therefore we can treat the closed geometry as the flat one.

Since we have analyzed the exponential bounce, we now proceed to the investigation of other bouncing solutions. In particular, we will focus on the power-law bouncing scale factor of the form

[TABLE]

with $n$ a positive integer, $a_{0}\,,\alpha$ positive parameters, and $t_{0}>0$ the bounce point. This relation leads to

[TABLE]

which implies that the post-bounce expansion is accelerated. Inserting (4.37) into the first Friedmann equation (4.20) we acquire

[TABLE]

Since for the open universe case $\rho$ can become negative (in particular, $\rho=-3/(a_{0}\kappa)^{2}$ at $t=t_{0}$), we focus on the $k=0$ and $k=+1$ cases, where $\rho$ is positive definite. Taking the derivative of (4.39) we find

[TABLE]

thus near the bounce point we have

[TABLE]

from which we deduce that the energy density decreases in the contracting phase before the bounce and increases immediately after it. Nevertheless, for $|t|\gg t_{0}$ we have

[TABLE]

which implies that after a suitable amount of time in the expanding phase $\rho$ starts decreasing again. Finally, the cosmological parameter $\Omega=1+\frac{k}{a^{2}H^{2}}$ behaves as

[TABLE]

and therefore it exhibits a decreasing behavior. Hence, similarly to the case of the exponential bounce analyzed earlier, such behaviors could be interesting for the description of the post-bouncing universe and the correct subsequent thermal history, since it will leave a universe with $\Omega$ very close to $1$ and a decreasing $\rho$.

Lastly, let us investigate what kind of fluids with equation of state given by (4.22) can produce the power-law bouncing solution (4.37). Considering an inhomogeneous viscous fluid with equation of state

[TABLE]

and inserting (4.37), we find the bulk viscosity as

[TABLE]

Note that away from the bouncing point, namely when $a(t)\gg a_{0}$, the bulk viscosity becomes

[TABLE]

Hence, if $0<\zeta<2/3$, which corresponds to $n>1/2$, then the bounce can be realized. On the other hand if $2/3<\zeta$ then, as we mentioned earlier, singularities might appear.

In summary, in this subsection we saw that viscous fluids can offer the mechanism to violate the null energy condition, which is the necessary requirement for the bounce realization. Hence, various bouncing solutions can be realized, driven by fluids with suitably reconstructed viscosity. As specific examples we studied the exponential and the power-law bounces, which are also capable of describing the accelerated post-bouncing phase, with the additional establishment of the spatial flatness. These features reveal the capabilities of viscosity.

4.4 Inclusion of isotropic turbulence

In this subsection we discuss turbulence issues in the framework of viscous cosmology. From hydrodynamics point of view the inclusion of turbulence in the theory of the cosmic fluid seems most natural, at least in the final stage of the universe’s evolution when the fluid motion may well turn out to be quite vigorous. The local Reynolds number must then be expected to be very high. On a local scale this brings the shear viscosity concept into consideration, as it has to furnish the transport of eddies over the wave number spectrum until the local Reynolds number becomes of order unity, marking the transfer of kinetic energy into heat. Due to the assumed isotropy in the fluid, we must expect that the type of turbulence is isotropic when looked upon on a large scale. According to standard theory of isotropic turbulence in hydrodynamics we then expect to find a Loitziankii distribution for low wave numbers (energy density varying as $k^{4}$), whereas for higher $k$ we expect an inertial subrange in which the energy distribution is

[TABLE]

where $\alpha$ denotes the Kolmogorov constant and $\epsilon$ is the mean energy dissipation per unit mass and unit time. When $k$ reaches the inverse Kolmogorov length $\eta_{K}$, i.e.

[TABLE]

with $\nu$ the kinematic viscosity, then the dissipative region is reached.

In the following we will consider a dark fluid developing into the future from the present time $t=0$, when turbulence is accounted for. We will perform the analysis in two different ways: either assuming a two-fluid model with one turbulent constituent, or assuming simply a one-component fluid, following [9, 68, 152, 153].

We start by considering a two-component model, where the effective energy is written as a sum of two parts, namely

[TABLE]

with $\rho$ denoting the conventional energy density. Taking $\rho_{\rm turb}$ to be proportional to the scalar expansion $\theta=3H$, and calling the proportionality factor $\tau$, we acquire

[TABLE]

Additionally, the effective pressure $p_{\rm eff}$ is split in a similar way as

[TABLE]

For both components we assume homogeneous equations of state, namely

[TABLE]

The Friedmann equations can thus be written (recall that $\kappa=8\pi G$) as

[TABLE]

[TABLE]

leading to the following governing equation for $H$:

[TABLE]

where we used the standard notation

[TABLE]

Finally, when the energy dissipation is assumed to be

[TABLE]

the energy balance may be written as

[TABLE]

In summary, the input parameters in this model are $\{w,w_{\rm turb},\tau\}$, all of them assumed to be constants. In the following we analyze the cases of two specific choices for $w$ and $w_{\rm turb}$.

• •

The case $w_{\rm turb}=w<-1$

This assumption implies that we equalize the ordinary and turbulent components as far as the EoS is concerned. From Eq. (4.55) we acquire

[TABLE]

Hence, we have a Big Rip singularity after a finite time

[TABLE]

and we obtain correspondingly

[TABLE]

In the vicinity of $t_{s}$, using that $Z=1-t/t_{s}$, we find

[TABLE]

[TABLE]

which reveal the same kind of behavior for $H$ and $a$ as in conventional cosmology, nevertheless the singularity in $\rho$ has become more weak. The physical reason for this is obviously the presence of the factor $\tau$.

It is interesting to see how these solutions compare with our assumed form (4.58) for the energy equation. The left hand side of Eq. (4.58) can be calculated, and we obtain in the limit $t\rightarrow t_{s}$ (details omitted here) the following expression for the present energy dissipation:

[TABLE]

This result could hardly have been seen without calculation; it implies that the specific dissipation $\epsilon_{0}$ is closely related to the EoS parameter $\gamma$ and the parameter $\tau$.

• •

The case $w<-1,~{}w_{\rm turb}>-1$

In general, the turbulent component is accordingly not only a passive component in the fluid. The assumption of the present case, namely $w<-1,~{}w_{\rm turb}>-1$, encompasses the region $-1<w_{\rm turb}<0$, in which the turbulent pressure will be negative as before. However, it also covers the region $w_{\rm turb}>0$, where the turbulent pressure becomes positive as in ordinary hydrodynamics.

The governing equation (4.55) can be solved with respect to $t$ as

[TABLE]

showing that the kind of singularity encountered in this case is of the Little Rip type. As $t\rightarrow\infty$, the Hubble function $H$ approaches the finite value

[TABLE]

Physically, $\gamma_{\rm turb}$ plays the role of softening the evolution towards the future singularity.

We close this subsection by investigating the case of a one-component scenario. In particular, instead of assuming the fluid to consist of two components as above, we can introduce a one-component model in which the fluid starts from $t=0$ as an ordinary viscous non-turbulent fluid, and then after some time, marked as $t=t_{*}$, it enters a turbulent state of motion. This picture is definitely closer to ordinary hydrodynamics.

Let us follow the development of such a fluid, assuming as previously that $w<-1$, in order for the fluid to develop towards a future singularity. After the sudden transition to turbulent motion at $t_{*}$, we have that $w\rightarrow w_{\rm turb}$ and correspondingly $p_{\rm turb}=w_{\rm turb}\,\rho_{\rm turb}$. Similarly to the two-component scenario, we assume $w_{\rm turb}>-1$, and for simplicity we assume that $\zeta$ is a constant.

We can now easily solve the Friedmann equations, requiring the density of the fluid to be continuous at $t=t_{*}$. It is convenient to introduce the “viscosity time”, namely

[TABLE]

Hence, for $0<t<t_{*}$ we obtain [87, 86]:

[TABLE]

[TABLE]

[TABLE]

whereas for $t>t_{*}$ we acquire:

[TABLE]

[TABLE]

[TABLE]

Thus, the density $\rho$ at first increases with time, and then decreases again until it goes to zero as $t^{-2}$ when $t\rightarrow\infty$. Note that in the turbulent region, $p_{*}=w_{\rm turb}\,\rho_{*}$ will even be greater than zero in the case where $w_{\rm turb}>0$.

As a final remark of this subsection, we mention that the presence of turbulence may alternatively be dealt with in terms of a more general equation of state of the form (1.15), admitting inhomogeneity terms too.

4.5 Viscous Little Rip cosmology

As discussed in detail in subsection 3.2 above, it is well known that there exist several theories for singularities in the future universe [51, 101]. Amongst them, the Little Rip scenario proposed by Frampton et al. [132, 133] (for nonviscous fluids) is an elegant solution, which we will consider in more detail in this subsection, generalized to the case of viscous fluids. The essence of the original model, as well as of its viscous counterpart, is that the dark energy is predicted to increase with time in an asymptotic way, and therefore an infinite span of time is required to reach the singularity. This implies that the equation-of-state parameter is always $w<1$, but $w\rightarrow-1$ asymptotically. In the following we will survey the essentials of this theory, as were developed by Brevik et al. [67]. In most cases the appearance of a bulk viscosity turns out to promote the future singularity.

For concreteness we assume an equation of state of the form

[TABLE]

where $A$ is a constant and $\xi(H)$ a viscosity function (not the viscosity itself). This is an inhomogeneous equation of state. Assuming a spatially flat FRW universe the first and second Friedmann equations write as

[TABLE]

while the conservation equation for energy, namely ${T^{0\nu}}_{;\nu}=0$, becomes

[TABLE]

Let us study separately the non-viscous and viscous cases.

• •

(I) Non-viscous case.

For comparison, we start from the non-viscous case $\xi(H)=0$ [132]. Setting the present scale factor $a_{0}$ equal to one, we obtain

[TABLE]

This relation reveals the Little Rip property: the singularity $\rho\rightarrow\infty$ is not reached in a finite time. Additionally, the density $\rho$ can be expressed as a function of the scale factor as

[TABLE]

Using the first Friedmann equation we can also express $a$ as a function of $t$, namely

[TABLE]

• •

(II) Viscous case.

Let us now switch on the viscous term in (4.74). In this case the second Friedmann equation, as well as the energy conservation equation, will change. We shall consider here only the simplifying ansatz where the viscosity function is constant, namely

[TABLE]

This choice is motivated mainly from mathematical reasons. Then, from the governing equations above, it follows that

[TABLE]

Inverting this equation we acquire

[TABLE]

Hence, the state $\rho\rightarrow\infty$ can indeed be reached, however it requires an infinite time interval. This is precisely the Little Rip characteristic, now met under viscous conditions. The term $\xi_{0}/A$ multiplying the exponential tends to promote the singularity, as mentioned. The influence from the last term $\xi_{0}/A$ becomes negligible at large times.

4.6 Viscous cosmology and the Cardy-Verlinde formula

In this subsection we will discuss the connection of viscous cosmology with thermodynamics. The apparent deep connection between general relativity, conformal field theory (CFT), and thermodynamics, has aroused considerable interest for several years. In the following we will consider one specific aspect of this subject, namely to what extent the Cardy-Verlinde entropy formula remains valid if we allow for bulk viscosity in the cosmic fluid. For simplicity we will assume a one-component fluid model, and we assume the bulk viscosity $\zeta$ to be constant. For more details, the reader may consult Refs. [154, 87, 86, 155, 156], and additionally the related Ref. [157].

We start with the Cardy entropy formula for an (1+1) dimensional CFT:

[TABLE]

where $c$ is the central charge and $L_{0}$ the lowest Virasoro generator [158, 159]. Comparing with the first Friedmann equation for a closed universe ($k=+1$) when $\Lambda=0$, namely

[TABLE]

we deduce (as pointed out by Verlinde [160]) that formal agreement is achieved if we choose

[TABLE]

where $E=\rho V$ is the energy in the volume $V$. One noteworthy fact is evident already at this stage: the correspondence is valid also if the fluid possesses viscosity, since there is no explicit appearance of viscosity in the first Friedmann equation. Moreover, the equation of state for the fluid is so far not involved.

In order to highlight the physical importance of the formal substitutions (4.85), let us consider the thermodynamic entropy of the fluid. As is known, there exist several definitions, the Bekenstein entropy, the Bekenstein-Hawking entropy, and the Hubble entropy. We will consider only the last quantity here, called $S_{H}$. Its order of magnitude can be easily estimated by observing that the holographic entropy $A/4G$ ($A$ is the area) of a black hole with the same size as the universe may be written in the form

[TABLE]

since $A\sim H^{-2}$ and hence $V\sim H^{-3}$. Various arguments have been provided to assume the universe’s maximum entropy to be identified with the entropy of a black hole having the same size as the Hubble radius [164, 161, 162, 163]. Nevertheless, more precise arguments of Verlinde [160] lead to the replacement of the factor 4 in the denominator with a factor 2, that is

[TABLE]

Therefore, one can see that this relation coincides with the last relation of (4.85), indicating that the formal substitutions above have a physical basis.

Consider now the Casimir energy $E_{C}$, defined in this context to be

[TABLE]

We may make use of scaling arguments for the extensive part $E_{E}$ and the Casimir part $E_{C}$ that make up the total energy $E$. These arguments finally give (details omitted here) $E(S,V)=E_{C}(S,V)+\frac{1}{2}E_{C}(S,V)$. An essential point is the property of conformal invariance, that the products $E_{E}a$ and $E_{C}a$ are volume independent and depend only on $S$. Hence, we acquire

[TABLE]

where $\alpha,\beta$ are constants. Their product arises from CFT arguments as $\sqrt{\alpha\beta}=3$ for $n=3$ spatial dimensions. From the formulae above we obtain

[TABLE]

which is the Cardy-Verlinde formula. With the substitutions $Ea\rightarrow L_{0}$ and $E_{C}\rightarrow c/12$ it is seen that expressions (4.89) and (4.83) are in agreement, apart from a numerical prefactor. This is caused by our assumption about $n=3$ spatial dimensions instead of the $n=1$ assumption in the Cardy formula.

The above arguments were made for a radiation dominated, conformally invariant, universe. Hence, the question that arises naturally is whether the same arguments apply to a viscous universe too. The subtle point here is the earlier pure entropy dependence of the product $Ea$, which is now lost. To analyze this question we may consider the following equation, holding for a $k=1,\Lambda=0$ universe with EoS $p=\rho/3$, namely

[TABLE]

This is essentially an equation for the rate of change of the quantity $Ea$. Let us compare this relation with the entropy production formula

[TABLE]

where $n$ is the particle number density and $\sigma$ the entropy per particle. As we observe, both time derivatives in (4.92) and (4.91) are proportional to $\zeta$. If $\zeta$ is small we can insert the usual solution for the scale factor of the nonviscous case, namely $a(t)=\sqrt{(8\pi G/3)\rho_{\rm in}a_{\rm in}^{4}}\,\sin\eta$, with $\eta$ the conformal time (“in” denotes the initial time). As the densities $\zeta^{-1}\rho a^{4}$ and $\zeta^{-1}n\sigma$ can then be regarded as functions of $t$ (recall that $\zeta=$constant), we conclude that $\rho a^{4}$ can be regarded as a function of $n\sigma$. This implies in turn that $Ea$ can be regarded as a function of $S$. This property, originally based upon CFT, can thus be carried over to the viscous case too, assuming that the viscosity is small.

At this stage we should pay attention to the following conceptual point. The specific entropy $\sigma$ in (4.92) is a conventional thermodynamic quantity, whereas the identification $S\rightarrow HV/(2G)$ in (4.85) is based on the holographic principle. The latter entropy is identified with the Hubble entropy $S_{H}$, and thus we can set $n\sigma_{H}=H/(2G)$, with $\sigma_{H}$ the specific Hubble entropy. The quantity $\sigma_{H}$ is holography-based, whereas the quantity $\sigma$ is not.

Finally, note that the same kind of arguments can be also applied in the more general situation where the EoS has the form

[TABLE]

with $\gamma$ a constant. For the non-viscous case this analysis was performed by Youm [165], with the result

[TABLE]

Lastly, in this case the application to weak viscosity can also be performed as in [154, 87, 86], and when $\gamma=4/3$ the radiation dominated result is recovered.

5 Conclusions

From a hydrodynamicist’s point of view the inclusion of viscosity concepts in the macroscopic theory of the cosmic fluid seems most natural, as an ideal fluid is after all an abstraction (unless the fluid is superconducting). Modern astronomical and cosmological observations permit us to look back in history, evaluating the Hubble parameter up to a redshift $z$ of about 2. Armed with such observational data, and having at one’s disposal the formalism of FRW cosmology with bulk viscosity included, one would like to extrapolate the description of the universe back in time up to the inflationary era, or go to the opposite extreme and analyze the probable ultimate fate of the universe, which might well be in the form of a Big Rip singularity. In the present review we have undertaken this quite extensive program.

After fixing the notation in subection 1.1, we began in Section 2 with a presentation of the theory of the inflationary epoch, covering cold as well as warm inflation in the presence of bulk viscosity. We investigated in detail the viscosity effects on the various inflationary observables, showing that they can be significant. A point to be noted in this context is that viscous effects may be represented by a generalized and inhomogeneous equation of state.

In Section 3 we turned to viscous theory in the late universe. We considered the phantom era with its characteristic singularities. Additionally, we discussed how one can describe in a unified way the inflationary and late-time acceleration in the framework of viscous cosmology. The simplest way to achieve this task is to introduce scalar fields. Moreover, we investigated the cosmological scenario of holographic dark energy in the presence of a viscous fluid, a subject which is related to black hole thermodynamics.

In the final Section 4 of our review we dealt with specific topics. We classified various options for the ultimate fate of the universe. We gave an analysis of whether the magnitude of bulk viscosity derived from observations is sufficient to drive the cosmic fluid from the quintessence into the phantom region. Numerical estimates indicated that such a transition might well be possible. Furthermore, we investigated viscous bounce cosmology, and we made use of isotropic turbulence theory from hydrodynamics to describe the late cosmic fluid. Moreover, we discussed the Little Rip occurrence in the presence of viscosity. Finally, we examined how viscosity influences the Cardy-Verlinde formula, which is a topic that relates cosmology with thermodynamics, and falls within the emergent gravity program.

We close this work mentioning that it would be both interesting and necessary to apply the cosmography formalism [55, 166, 167, 168] in order to impose constraints on the viscosity parameters. Contrary to standard observational constraints, the advantage of cosmography is that it is model-independent, since one expands the scale factor independently of the solution of the cosmological equations. In particular, one introduces $H=\frac{1}{a}\frac{da}{dt}$, $q=-\frac{1}{a}\frac{d^{2}a}{dt^{2}}\ H^{-2}$, $j=\frac{1}{a}\frac{d^{3}a}{dt^{3}}\ H^{-3}$, $s=\frac{1}{a}\frac{d^{4}a}{dt^{4}}\ H^{-4}$, $l=\frac{1}{a}\frac{d^{5}a}{dt^{5}}\ H^{-5}$, known respectively as Hubble, deceleration, jerk, snap and lerk parameters [55]. One can show that these parameters are related through [166, 167, 168],

[TABLE]

which can easily lead to the distance - redshift relation [55]. Hence, using the Friedmann equation in the case of viscous cosmology, we can relate the above cosmographic quantities with the present value of viscosity parameter. This could be a significant advantage, since the obtained constraint would be model-independent, and thus more robust that the ones discussed in subsection 4.1. Nevertheless, the detailed investigation of this subject lies beyond the scope of this review, and it is left for a future project.

Mostly, this review is based on a theoretical approach. We have however provided information concerning quantities related to observations, giving estimations on the inflationary observables, as well as on the magnitude of the current bulk viscosity itself.

In summary, from the above analysis one can see the important implications and the capabilities of the incorporation of viscosity, which make viscous cosmology a good candidate for the description of Nature.

Acknowledgments

This work is supported by MINECO (Spain), Project FIS2013-44881, FIS2016-76363-P and by CSIC I- LINK1019 Project (S. D. O. and I. B.) and by Ministry of Education and Science of Russia, Project N. 3.1386.2017.