A hybrid model of viscous and Chaplygin gas to tackle the Universe acceleration

Motivated by two seminal models proposed to explain the Universe acceleration, this paper is devoted to study a hybrid model which is constructed through a generalized Chaplygin gas with the addition of a bulk viscosity. We call the model a Viscous Generalized Chaplygin Gas (VGCG) and its free parameters are constrained through several cosmological data like the Observational Hubble Parameter, Type Ia Supernovae, Baryon Acoustic Oscillations, Strong Lensing Systems, HII Galaxies and using Joint Bayesian analysis. In addition, we implement a Om-diagnostic to analyze the VGCC dynamics and its difference with the standard cosmological model. The hybrid model shows important differences when compared with the standard cosmological model. Finally, based on our Joint analysis we find that the VGCG could be an interesting candidate to alleviate the well-known Hubble constant tension.


I. INTRODUCTION
The origin of the accelerated expansion of the Universe in recent epochs is believed to be caused by an unknown component of the universe called the Dark Energy (DE). DE is one of the cornerstones in the modern cosmology. Despite that the DE is a strange phenomenon which evolves into an Einstein-D'Sitter structure, it is confirmed by several observations such as the Cosmic Microwave Background Radiation (CMB) [1], the distance modulus of type Ia Supernovae (SnIa) [2][3][4] and recently through the large-scale structure studies [5].
The simplest candidate to explain DE comes from adding a cosmological constant (CC with value Λ) [6] in the Einstein field equations, i.e., a geometrical correction to them. On the other hand, such CC can also be seen as an effective fluid with an equation of state (EoS) w = −1 when it is included through the energy-momentum tensor. Nevertheless, the CC origin as quantum vacuum fluctuations is not well understood [7,8], encouraging the necessity of new approaches to understand its nature. Motivated by the CC problems, the community has been * ahalmada@uaq.mx exploring alternatives to the CC which can be classified in two kinds of solutions: modify the Einstein field equations or include effective fluids to them. In this sense, a zoo of DE models has been appearing in the literature. We can list for instance, models related to modified gravity as f (R) theories [9][10][11], the extra dimensions [12,13], unimodular gravity [14][15][16][17]. Conversely, models associated to effective fluids as DE are phenomenological dark energy which is related with emergent DE models [18][19][20], quintessence and phantom DE [21], dissipative fluids [22][23][24], Chaplygin gas (CG) [25,26], Generalized Chaplygin Gas (GCG) [27][28][29], and like-Chaplygin gas [30]. For a review of these models see for instance [21].
An interesting property of Chaplygin fluids is that they are able to describe Dark Matter (DM) and DE as an unique perfect fluid through an EoS of the form p c = −Aρ −α c , where 0 ≤ α ≤ 1 and A > 0, where p c and ρ c are pressure and density respectively. This kind of models have been analyzed for inflationary epochs through the formalism of Hamilton-Jacobi [31], and later were studied at late epochs [27][28][29][30]. In spite of the advantages, this kind of models experience several problems such as oscillations and exponential blowup in the DM power spectrum, in disagreement with the observations [32], which in turn affects the formation of structures. Additionally, several authors [33] have shown that the trajectories in the context of dynamical systems are unstable, inducing a fine tuning in the initial condition.
On the other hand, non perfect fluids are another approach to tackle the Universe acceleration. In fact, they also allow the description of DM and DE in the context of unified fluids through the inclusion of kinds of viscosity, bulk and shear. However, the bulk viscosity is the best candidate to face the DE problem because satisfies the cosmological principle and it could be studied through two main approaches: Eckart formalism [34] and Israel-Steward (IS) theory [35]. The main difference between them is that the Eckart formalism is a non causal theory and, as a consequence, it is simpler than the IS. Indeed, only one bulk viscosity model of the form ξ ∼ ρ s (where ρ is the energy density considered typically as dust matter) has been studied [36][37][38][39][40] in the IS theory, while several models have been explored [41][42][43][44][45][46][47][48] in the Eckart one.
In this paper, we present a hybrid model which consists of a fusion of viscous effects and the features of Chaplygin gas, hereafter Viscous Generalized Chaplygin Gas (VGCG), to tackle the late accelerated expansion of the Universe. Naturally, VGCG is also under the unified dark fluids approach as GCG and viscous fluids. This kind of hybrid models was originally proposed by [49] and are able to avoid causality problems that arise when only dissipative fluids are considered. This model has been studied by several authors [50][51][52][53][54][55]. For instance, [54] studied the spherical top-hat collapse both in Einstein gravity and in loop quantum cosmology. Additionally, cosmological bounds of the free parameters were established by [52], finding that the viscosity terms (ξ 0 = 7.08 × 10 −4 ) may affect the CMB power spectrum mainly on the height of the acoustic peaks related to the matter density but consistent to those obtained by ΛCDM. Additionally, one of the most important points to propose this kind of hybrid models is that they could alleviate the oscillations that cause the blowup in the DM power spectrum for GCG models [32,56]. The problem lies in the equation for the density fluctuations δ k , in particular in the nonzero sound speed where, under some regime, the fluctuations become violently unstable and grow exponentially. Based on the results regarding CMB power spectrum for VGCG [52], the viscosity terms in the VGCG could also act as a damping term in the DM power spectrum as occurs in the CMB spectrum (ξ 0 = 7.08×10 −4 ), reducing such instabilities. Therefore, this kind of models must be studied in detail, starting from the background cosmology.
This paper is devoted to analyze this hybrid model considering the Eckart formalism to add its viscous effects from the cosmographics and thermodynamics point of view. Firstly, we constrain the parameter space of the VGCG by performing a Bayesian analysis and using the following samples: observational Hubble parameter (OHD) [57], type Ia supernovaes (SnIa) [4], baryon acoustic oscillations (BAO) [58], strong lensing systems (SLS) [59], and HII Galaxies (HIIG) [60][61][62][63][64]. Furthermore, we also perform a Bayesian joint analysis combin-ing the mentioned samples to study the cosmography of the model, {q, j} phase-space and a test related to the well-known Om-diagnostic. Finally, we discuss the near equilibrium condition imposed by the Thermodynamics to any dissipative fluid.
The manuscript is structured as follow: Sec. II presents the GCG model in presence of viscosity and give general expressions to reconstruct and analyze the cosmography parameters. Sec. III gives a summary of the cosmological data and presents the methodology to establish bounds on the phase-space of the model parameters. Sec. IV discusses the results and finally, Sec. V makes a summary and shows the conclusions.

II. VISCOUS GENERALIZED CHAPLYGIN GAS REVISITED
In a spatially flat Universe described by the Friedman-Lemaitre-Robertson-Walker (FLRW) metric, a generalized Chaplygin gas in presence of viscosity is studied. The model describes the evolution of baryons (b) and relativistic species (r) as perfect fluids with EoS w b = p b /ρ b = 0 and w r = p r /ρ r = 1/3 respectively, where ρ b,r is their energy density. Additionally, it is considered a non perfect fluid following the GCG EoS p c = −Aρ −α c , where p c and ρ c are the pressure and energy density respectively and A and α are appropriate constants. Based on Eckart formalism, the viscous effects are introduced in the energy-momentum tensor T µν = ρ t U µ U ν + p ef f (g µν + U µ U ν ) as an effective pressure, p ef f , where U µ is the four-velocity. In this context, the total effective pressure is the sum of the barotropic pressure of the fluids and the viscous pressure, where Π = −3ξH, ξ is the bulk viscosity, p i (i = b, r, c) is the pressure of the fluid and H ≡ȧ/a is the Hubble parameter. Then the Friedmann and continuity equations are where κ 2 ≡ 8πG and dots are derivatives with respect to cosmic time. By integrating the first two continuity equations, we have ρ i = ρ i0 a −3(wi+1) where ρ i0 is the initial energy density, for i = b, r. To integrate the one for the VGCG we will assume that ξ = ξ 0 ρ 1/2 c , being ξ 0 an appropriate constant and H 2 ≈ κ 2 ρ c /3 (a good approximation for the late stage of the Universe), thus we obtain Hence, the dimensionless Hubble parameter E(z) ≡ H(z)/H 0 in terms of the redshift z is given by where we have used the relation a = (1 + z) −1 and which allows us to recover the GCG perfect fluid case. It is useful to explore the cosmography parameters, mainly the deceleration and jerk parameters, to distinguish the viscous effects included in the GCG through the {q, j}-phase space, and also to compare VGCG with GCG and ΛCDM.
The deceleration parameter can be estimated by q(z) = −1 + (1 + z)E(z) −1 (dE(z)/dz), then using (7) we obtain, and the jerk parameter takes the form where q is given by Eq. (8). Finally, it is also useful to analyse the viscous effects through a test related to the well-known Om-diagnostic proposed in [65], illustrated in a E 2 versus (1 + z) 3 panel, that gives a clear evidence of whether the cosmological models behave as phantom and/or quintessence and their deviations with respect to ΛCDM. In this context, ΛCDM is represented by a trajectory with constant slope and a upper (lower) trajectory for the other cosmological models over the ΛCDM one for z > 0 suggests a quintessence (phanthom) behaviour.

III. DATA AND METHODOLOGY
A Bayesian Markov Chain Monte Carlo (MCMC) analysis is performed to constrain the phase-space parameter Θ of the viscous generalized Chaplygin gas using OHD, SnIa, SLS, BAO, HIIG data and a joint analysis. Then, under the environment of the emcee Python package [66], we use the Gelman-Rubin criterion [67] to achieve the convergence of the chains, and then we set 4000 chains with 250 steps each one to establish bounds on the model parameters. The phase-space priors are: Gaussian for Ω b0 h 2 = 0.02242 ± 0.00014 [1] and flat distributions for the rest of the parameter in the ranges h : [0.6, 0.8], B s : [0, 2], α : [0, 2], and ξ 0 : [0, 1]. For the joint analysis, the figure-of-merit to be optimized is given by where the χ 2 sub indexes correspond to the names of the samples. The rest of the section is devoted to discuss the details of the samples used to constrain Θ at background level.

A. Observational Hubble Parameter
The observational Hubble parameter data (OHD) represents the most direct way to constrain the parameter space. The largest OHD sample is compiled by [57] and contains a total of 51 points covering the range 0 < z < 2.36. From them, only 31 points are cosmological independent measurements obtained using differential age (DA) tools and the rest comes from BAO measurements. Because these measurements are considered to be uncorrelated, the chi square function can be expressed as where H th (z i ) is the theoretical estimate using (7) and H obs (z i ) ± σ i obs is the observational Hubble parameter with its uncertainty at the redshift z i .

B. Type Ia Supernovae
Authors in Ref. [4] provide 1048 luminosity modulus measurements, known as Pantheon sample, from Type Ia Supernovae which cover a region 0.01 < z < 2.3. Due to in this sample the measurements are correlated, it is convenient to build the chi square function as where e = ∆1 T · Cov −1 P · ∆1, and ∆μ is the vector of residuals between the theoretical distance modulus and the observed one, ∆1 = (1, 1, . . . , 1) T , Cov P is the covariance matrix formed by adding the systematic and statistic uncertainties, i.e. Cov P = Cov P,sys + Cov P,stat . The super-index T on the above expressions denotes the transpose of the vectors.
The theoretical distance modulus is estimated by where M is a nuisance parameter which has been marginalized in (12). The dimensionless luminosity distance, denoted as d L (z), is computed by where c is the speed of light.

C. Baryon Acoustic Oscillations
Another way to establish a constraint of the model parameter is through the standard rules known as Baryon Acoustic Oscillations (BAO). These are primordial signatures produced by the interaction between baryons and photons in a hot plasma in the pre-recombination epoch. Authors in [58] collected 15 transversal BAO scale measurements which come from luminous red galaxies distributed in the redshift range 0.110 < z < 2.225. As these points are considered uncorrelated, the chi square function is built as where θ i BAO ± σ θ i BAO is the BAO angular scale and its uncertainty at 68% measured at z i . On the other hand, the theoretical BAO angular scale, denoted as θ th , is estimated by where D A = d L (z)/(1 + z) 2 is the angular diameter distance at z and d L (z) was defined in (15). The parameter r drag is defined by the sound horizon at baryon drag epoch. For this work, we set the r drag = 147.21 ± 0.23 obtained by Planck collaboration [1].

D. Strong lensing systems
A compilation of 204 strong lensing systems (SLS) in the redshift 0.0625 < z l < 0.958 for the lens and 0.196 < z s < 3.595 for the source is provided by [59]. The chi square function for SLS takes the form where the observable to confront is D obs = c 2 θ E /4πσ 2 , where θ E is the Einstein radius of the lens obtained by assuming the gravitational lens potential is modeled by a Singular Isothermal Sphere (SIS) defined by In the above expression, σ SIS is the velocity dispersion of the lens galaxy, D s is the angular diameter distance to the source, and D ls is the angular diameter distance from the lens to the source. It is interesting to remark that as this SLS data assumes a lens model for θ E and σ SIS comes from spectroscopy, the sample is independent of h, and as consequence the parameter constraints do not depend on h. The uncertainty of D obs is estimated by where δθ E and δσ are the uncertainties of the Einstein radius and the observed velocity dispersion, respectively. The theoretical counterpart is estimated by the ratio where D ls is the angular diameter distance from the lens to the source given by and D s = d L (z)/(1+z) 2 is the angular diameter distance to the source which is obtained using (15).

E. HII Galaxies
Authors [60-64, and references therein] argued that the correlation between the measured luminosity L and the inferred velocity dispersion σ of the ionized gas (e.g. Hβ, Hα, [OIII] emission lines) in extreme starburst galaxies (i.e. containing a population of O and/or B stars) may be used as a cosmological tracer to constrain cosmological model parameters. A compilation of 153 HII galaxies (HIIG), containing apparent magnitude, emission line luminosity and velocity dispersion is provided by [64,68]. Then, the chi square function is estimated by where In the above expressions, µ i obs ± σ i obs is the observed distance modulus with its uncertainty at redshift z i . The theoretical estimate at the redshift z is obtained by using (15) and where µ 0 is a nuisance parameter which has been marginalized. Table I presents the cosmological constraints of VGCG and GCG for OHD, SnIa, SLS, BAO, HIIG samples and the joint analysis respectively. Each best-fit parameter value includes its uncertainty at 68% confidence level (CL) and is consistent, within 1σ, with those reported in the literature [52]. Figure 1 shows the 1D marginalized posterior distributions for each data and joint analysis and also the 2D phase space distribution at 68% (1σ), 99.7% (3σ) CL. According to the χ 2 value, both models are in good agreement with the data 1 . However, as the VGCG includes an extra parameter over GCG, it is better to confront them statistically using the corrected Akaike information criterion (AICc) [69][70][71] and the Bayesian information criterion (BIC) [72] defined as AICc = χ 2 min + 2k + (2k 2 + 2k)/(N − k − 1) and BIC = χ 2 min + k log(N ) respectively, where k is the number of free parameters and N is the number of data points. The criteria establishes that the model with lower values of AICc and BIC is preferred by data. For AICc, if a difference between a given model and the best one, ∆AICc, is ∆AICc < 4, both models are supported by the data equally. If 4 < ∆AICc < 10 the data still support the given model but less than the preferred one. A value of ∆AICc > 10 indicates that the data does not support the given model. For BIC, the difference between a candidate model and the best model ∆BIC is interpreted as the evidence against the candidate model being the best model. A yield of ∆BIC < 2 indicates there is no evidence against the candidate model. A value in the range 2 < ∆BIC < 6 suggests that there is modest evidence against the candidate model. A strong evidence against the candidate model is given when 6 < ∆BIC < 10, and a stronger evidence against is if ∆BIC > 10. In summary, based on our ∆AICc results given in Table I, the data show a similar preference for both models, VGCG and GCG. In contrast, ∆BIC results give a strong evidence against VGCG for SLS and SnIa data and no evidence against for the rest of the samples, including the combined sample (Joint). These results are expected considering that BIC penalizes free parameters more strongly than AICc and that the penalization function is proportional to the number of points in the dataset.

A. Cosmography parameters
Based on the joint analysis, the parameters of deceleration (q) and jerk (j) at z = 0 are estimated, obtaining respectively values of q 0 = −0.523 +0.021 VGCG where the uncertainties are at 1σ. These q 0 values are in agreement with the ΛCDM value [73] within 1.6σ and 1.1σ respectively. Regarding the jerk parameter, our results are consistent within 1σ with respect to the expected ΛCDM value. Additionally, it is interesting to compare our q 0 (j 0 ) values of VGCG with other viscous models. For instance, we find a good agreement with respect to models which consider a constant viscosity within 1.5σ (0.6σ) and polynomial viscosity within 1σ (1.2σ) (see Fig. 4 of [73]). It is worth mention that although the VGCG model differs from the mentioned  models in the EoS (w = −1) and the viscosity coefficient (constant and polynomial), we find consistent results in the current values of q and j. An alternative way to find differences between cosmological models is through the {q, r}-panel where r = j 2 . Figure 3 shows such {q, r}-space reconstruction for both models and ΛCDM in the range −1 < z < 150. The arrows represent the direction of the universe evolution, being the evolution from a deceleration phase to an accelerated epoch. The black square markers over the trajectories represent the current states and the bands around 2 Typically this space is shown by using the notation r instead of j.
the VGCG trajectory are its uncertainties up to 3σ. It is interesting to observe that for both models, their states q, r converge to ΛCDM state in the future (z → −1) in contrast of the observed behaviour by considering other viscous models as those presented in [73]. Additionally, we observe the path of GCG is close to the one by ΛCDM, while for the VGCG we find a deviation of more than 3σ in the deceleration epochs from ΛCDM model. B. Near equilibrium condition Figure 4 displays E(z) 2 with respect to (1 + z) 3 for VGCG, GCG and ΛCDM. This panel shows that VGCG and GCG behave as quintessence in the past (0 < z < 2.5) and as phantom in the future (−1 < z < 0). Furthermore, we observe agreement up to 3σ between VGCG and GCG, and a deviation of more than 3σ between VGCG and ΛCDM in z 2.3. Because the VGCG model produces an accelerated expansion stage using two mechanism (the viscous pressure and the EoS), this model does not need to satisfy the near equilibrium condition required in the Eckart formalism. In this context, a natural question is how far the VGCG is from this condition given by |Π/p T | 1, where p T is the total pressure of the fluids. Figure 5 shows the behavior of |Π/p T | for VGCG with its uncertainties up to 3σ (blue bands) and the case of ΛCDM in presence of viscosity (VΛCDM) when ξ ∼ ρ 1/2 [74]. Although, we find a good agreement between both models, within 3σ in the region −0.5 < z < 2.5, we observe that VGCG (VΛCDM) presents a redshift region z 0.2 (z 0.9) with |Π/p T | > 1 which are values lower than those found in [73], who find a redshift region z 4 (z 13) by considering a polynomial (constant) bulk viscosity using OHD+SnIa+SLS data.

V. CONCLUSIONS AND OUTLOOKS
This paper was dedicated to revisit a hybrid model that contains a non-perfect Chaplygin gas, called VGCG, to produce the accelerated expansion stage. VGCG was studied mainly in two parts at the background level. The first part was devoted to perform a Bayesian statistics analysis to confront the model with the most recent cosmological data (OHD, SnIa, BAO, SLS, HII Galaxies) and give updated values of the VGCG parameters. We found consistent constraints of the phase space parameters within 1σ. Based on our Joint analysis results shown in Table I, the second part was dedicated to explore the evolution of the cosmographic parameters and the near equilibrium condition required in the Eckart formalism. In particular, we analysed the {q, r}space which is useful to find differences between cosmological models. Thus, we observed a deviation of VGCG with respect to the ΛCDM and GCG of more than 3σ CL in the deceleration phase (see Figure 3). Additionally, we reported the current cosmographic values and l 0 = −0.110 +1.960 −3.701 for the VGCG model corresponding to deceleration, jerk, snap, and lent parameters respectively. We obtained consistent values with those values reported in [73] for alternative viscous models. On the other hand, we found a quintessence behaviour for VGCG and GCG models in the past z > 0 and a phantom one for the future as shown in Figure 4, with the difference that GCG presents divergence at z = −1. On the other hand, we reconstructed the |Π/p T | to explore how far the VGCG is from the near equilibrium condition. We found consistent results up to 3σ with the behaviour by considering a ΛCDM model with viscosity presented in [74]. In this context, we observed that VGCG presents a tighter phase around current epochs which |Π/p T | < 1 with respect to those regions for the viscous ΛCDM reported in [73,74].
Finally, based on our results presented in Table I, it is noteworthy that the VGCG suggests a central value on h between the value obtained by Planck [1] and the value by Riess et. al. [75], contributing to diminish the tension between these measurements. In fact, we remark that although we are using data based on local measurements such as OHD that favor the local value (Riess value), we obtain a central value on H 0 close to the Planck one. Hence, VGCG is a promising candidate to lessen the discrepancy generated in H 0 . Further studies of these models is needed.