Modified Equations of State for Dark Energy and Observational Limitations

Cosmological models with variable and modified equations of state for dark energy are confronted with observational data, including Type Ia supernovae, Hubble parameter data $H(z)$ from different sources, and observational manifestations of cosmic microwave background radiation (CMB). We consider scenarios generalizing the $\Lambda$CDM, $w$CDM, and Chevallier--Polarski--Linder (CPL) models with nonzero curvature and compare their predictions. The most successful model with the dark energy equation of state $w= w_0+ w_1(1-a^2)/2$ was studied in detail. These models are interesting in possibly alleviating the Hubble constant $H_0$ tension, but they achieved a modest success in this direction with the considered observational data.


I. INTRODUCTION
In modern cosmology based on recent observational data, our Universe includes dominating fractions of dark energy and dark matter, whereas all kinds of visible matter fill about 4% in total energy balance nowadays. The latest estimations of Planck collaboration [1,2] predict about 70% fraction of dark energy, if we apply the standard ΛCDM model, where dark energy may be represented as the cosmological constant Λ or as a matter with density ρ x and pressure p x = −ρ x . Almost all remaining part of matter in this model is cold dark matter with close to zero pressure. Because of the last property, it is convenient to consider cold dark matter together with visible baryonic matter, where the unified density is ρ m = ρ c + ρ b . One should also add the radiation component including relativistic species (neutrinos) with p r = 1 3 ρ r , which was sufficient before and during the recombination era, but is is almost negligible now.
The ΛCDM model successfully describes numerous observations, including Type Ia supernovae (SNe Ia) data, estimates of the Hubble parameter H(z), manifestations of baryonic acoustic oscillations (BAO), cosmic microwave background radiation (CMB), and other data [1,2]. However this model does not explain the nature of dark energy, the small observable value of the phenomenological constant Λ and the approximate equality ρ x and ρ m now (although these densities evolve differently).
The Planck estimations [1,2] are based on the ΛCDM model and the Planck satellite measurements of the CMB anisotropy and power spectra related to the early Universe at redshifts z ∼ 1000, whereas the SH0ES method uses local distance ladder measurements of Cepheids in our Galaxy [3] and in nearest galaxies, in particular, in the Large Magellanic Cloud [5], that implies z close to 0 (the late Universe). This H 0 tension has not diminished during the last years, and now it exceeds 4σ.
In this paper, we consider cosmological scenarios with modified equation of state (EoS) for dark energy; they generalize the ΛCDM model and its simplest extensions: the wCDM model with EoS p x = wρ x , w = const (1) and the models with variable EoS where w depends on the scale factor a. The class (2) includes the following well-known dark energy equations of state: the linear model [21] Chevallier-Polarski-Linder (CPL) parametrization [22,23] and (their generalization) Barboza-Alcaniz-Zhu-Silva (BAZS) EoS [24] Here, the scale factor a is normalized, so a(t 0 ) = 1 at the present time t 0 ; a is connected with redshift z: a = (1 + z) −1 . Obviously, BAZS parametrization (5) transforms into CPL EoS (4) at β = 1, into the linear EoS (3) if β = −1 and into the logarithmic EoS w = w 0 + w 1 log(1 + z), if β → 0.
For all mentioned equations of state, one can integrate the continuity equation for noninteracting dark energy:ρ In this paper, we explore the above scenarios and suggest the following generalization of BAZS EoS: We test these models, confronting them with the following observational data: Type Ia supernovae data (SNe Ia) from the Pantheon sample survey [25], data extracted from Planck 2018 observations [2,26] of cosmic microwave background radiation (CMB) ,and estimations of the Hubble parameter H(z) =ȧ/a for different redshifts z (the dot means d dt ). We use H(z) data from two sources: (a) from cosmic chronometers, that is, measured from differential ages of galaxies, and (b) estimates of H(z) obtained from line-of-sight baryonic acoustic oscillations (BAO) data. This paper is organized as follows. In the next section, we describe H(z), SNe Ia, and CMB observational data analyzed here. Section III is devoted to dynamics and free model parameters for scenarios (1)- (7). In Section IV we analyze the results of our calculations for these models, and estimate values of model parameters including the Hubble constant H 0 , and in Section V we discuss the results and their possible applications for alleviating the Hubble constant tension problem.

II. OBSERVATIONAL DATA
Observational data should be described by the considered cosmological models. For each model, we calculate the best fit for its free parameters from the abovementioned data sources: (a) Type Ia supernovae (SNe Ia) data from Pantheon sample [25], (b) CMB data from Planck 2018 [2,26], and (c) estimates of the Hubble parameter H(z) from cosmic chronometers and line-of-sight BAO data.
The Pantheon sample database [25] for SNe Ia contains N SN = 1048 data points of distance moduli µ obs i at redshifts z i in the range 0 < z i < 1.92. We compare them with theoretical values by minimizing the χ 2 function: Here, θ j are free model parameters, C SN is the covariance matrix [25], and the distance moduli µ th are expressed via the luminosity distance D L depending on the spacial curvature fraction Ω k and the Hubble parameter H(z): Here, m B is a supernova apparent magnitude, and M B is its absolute magnitude. The distance moduli µ th are not Hubble-free and depend on the Hubble constant H 0 = H(t 0 ) (via the summand −5 log 10 H 0 ). On the other hand, µ obs i is essentially connected with the absolute magnitude M B and calculated with corrections coming from deviations of lightcurve shape, SN Ia color, and mass of a host galaxy [25,27]. In the Pantheon sample, these corrections and the connected pair (H 0 , M B ) were considered as nuisance parameters, in particular: "Using only SNe, there is no constraint on H 0 since H 0 and M B are degenerate" [25]. We cannot divide uncertainties in the Hubble constant H 0 and possible uncertainties in M B [28][29][30].
Due to these reasons, we have to consider H 0 in Equation (8) as a nuisance parameter, its estimations cannot be obtained from χ 2 SN , and we minimize this function over H 0 [16][17][18][19][20]. However, SHe Ia data in χ 2 SN is important for fitting other model parameters. Unlike SNe Ia Pantheon data, the CMB observations are related to the photon-decoupling epoch near z * ≃ 1090. We use the following parameters extracted from Planck 2018 CMB observations [2,19,26]: and their estimations for the non-flat ΛCDM + Ω k model [26]: Considering the flat case (Ω k = 0) of these models, we use the flat wCDM data [26]. The comoving sound horizon r s at z * is calculated as the integral with the fitting formula from Refs. [26,31] for the value z * . The resulting χ 2 function is where we minimize over the normalized baryon fraction ω 0 b to diminish the effective number N p of free model parameters. The covariance matrix C CMB = C ij σ i σ j and other details are described in papers [16,26].

CC Data
We describe these uncertainties as corrections ∆C H to the diagonal covariance matrix C d H = diag{σ −2 i } (from Table I) taking into account their diagonal terms in the form [41] ( Here, η(z) is a mean percentage bias depending on redshift z. We consider "the best-case scenario" from the paper [41] for η(z) and include these contributions of stellar population synthesis and metallicity omitting the non-diagonal terms ∆(C H ) ij (they are negligible for metallicity [42]).

III. MODELS
We explore all considered models in a homogeneous isotropic universe with the Friedmann-Lemaître-Robertson-Walker metric where k is the sign of spatial curvature. In this case, the Einstein equations are reduced to the system of the Friedmann equation (14) and the continuity equationρ Here, the total density ρ includes densities of the abovementioned cold pressureless matter (dark matter unified with baryonic matter), radiation, and dark energy: We suppose here that dark energy and the mentioned components do not interact in the form [12][13][14][15] and independently satisfy the continuity Equation (6) or (15). We integrate this equation for cold and relativistic matter: (the index "0" corresponds to the present time t 0 ) and substitute these relations into the Friedmann equation (14) that can be rewritten as or Here, The dark energy fraction Ω x (a) results from the continuity Equation (6) In particular, for BAZS parametrization (5) w = w 0 + w 1 (1 − a β )/β [24], the expression (20) is One should substitute it into Equation (18).
In the above-mentioned particular cases, the BAZS formula (21) takes the form for the linear model (3) and for CPL EoS (4) w = w 0 + w 1 (1 − a) (if β = 1). In the case w 1 = 0 (and w 0 ≡ w), both models (22) and (23) transform into the wCDM model with Its particular case at w = −1 is the ΛCDM model where Ω x = Ω 0 x = Ω Λ = const. For the generalization (7) (with the factor γ) of BAZS parametrization (5), the expression (20) takes the form If γ → 0 it transforms into Equation (21), and in the case β → 0 we have For all considered models, the dark energy fraction Ω 0 x and other Ω j satisfy the equality Ω 0 m + Ω 0 r + Ω k + Ω 0 x = 1, resulting from Equation (18) or (19). Further, if we fix the ratio [17][18][19]58] to diminish the number N p of free model parameters. We will work with the following five free parameters in the linear (22) and CPL (23) models: Here, we consider ω b in χ 2 CMB (12) as a nuisance parameter. In the ΛCDM model the number of parameters N p = 3 (Ω 0 m , Ω k , H 0 ), in the wCDM model (24) N p = 4 with w 0 ≡ w, in the generalized model (7), (25) we have N p = 7 free parameters with additional β and γ to the set (28).
In the next sections, we compare predictions of these models with the observational data from Section II.

IV. RESULTS
We evaluate how the considered models fit the observations taking into account the χ 2 functions for SNe Ia (8), CMB (12), and H(z) data (13) in the form For the Hubble parameter data, we separately use (a) only 32 data points with cosmic chronometers (CC) H(z) data and (b) the full set CC + H BAO data (see Table I). When we compare two models with different number N p of free parameters, we can expect that the model with larger N p achieves more success in minimizing χ 2 tot . However, some models are not successful in this sense. In particular, the linear model (3), (22) with N p = 5 parameters (28), with our set of observational data, yields the same minimal value min χ 2 tot ≃ 1092.34 (for CC + H BAO ) as the wCDM model (24) with N p = 4. The reason is the following: the best fitted value w 1 for the linear model (3) is very close to zero (w 1 ≃ −0.0006). In this case, the linear model works as the wCDM model (24) and yields the same min χ 2 tot . Such a behavior of the model (3) is shown in Figure 1, where 1σ and 2σ filled contour plots in the left panels correspond to the full set of H(z) data (here and below, "all" denotes CC + H BAO ). The 1σ, 2σ contours for the wCDM model (red lines) in the Ω 0 m − Ω k plane behave similarly and closely; positions of χ 2 tot minima point (shown as the star and the circle) practically coincide. Here, the contours are drawn for χ 2 functions minimized over all other parameters, in particular, for the linear model (3) in the top-left panel: A similar picture also takes place for only CC H(z) data, where min χ 2 tot ≃ 1074.73 for both models. Equality of these minima for only CC and all H(z) data is illustrated with one-parameter distributions χ 2 tot (Ω 0 m ) in the top-right panel. In one-parameter distributions, we also minimize over all other parameters.
In the Ω 0 m − w 1 plane (the bottom-left panel of Fig. 1) we see that for both models, minima of χ 2 tot are achieved near w 1 = 0. It is also shown in the bottom-right panel, where the likelihood functions are depicted for the linear model. For the same observational data we can observe in Figure 2 more successful behavior of the Chevallier-Polarski-Linder (CPL) model (23) [22,23]. Here, 1σ and 2σ contours are drawn for χ 2 tot (θ i , θ j ) of the type (30) in all planes of two parameters (θ i − θ j ) in notation of Figure 1. We consider four cases: for two models (CPL and wCDM) we calculate χ 2 tot for CC and all H(z) data, positions of all χ 2 tot minima points are shown. Naturally, in the panels with w 1 , only the CPL model is presented.
The likelihood functions L(θ i ) of the type (31) are shown in Figure 2 for all five model parameters (28). They are used for estimating the best fits and 1σ errors for these parameters, summarized below in Table II. The CPL model achieves lower values of min χ 2 tot in comparison with the wCDM and ΛCDM models. It can be seen in Table II and in the top-right panel of Figure 2, where one-parameter distributions χ 2 tot (H 0 ) of these models are compared. The graphs χ 2 tot (H 0 ) and the correspondent likelihoods L(H 0 ) in the top-left panel demonstrate that the best fitted values H 0 are very close for the wCDM and CPL models and differ more essentially when we compare CC and all H(z) data.
The best fits of Ω 0 m depend stronger on the chosen model and vary from Ω 0 m ≃ 0.206 for CPL, all H(z) to Ω 0 m ≃ 0.289 for wCDM, CC. The best fits of Ω k behave similarly, but with the maximal estimate for CPL, all H(z). The best CPL fits for w 1 are close to 0.66 for both variants of H(z) data. This value is far from zero; in other words, the CPL model with the considered observational data behaves differently to the wCDM model. The CPL model achieves good results, but it is not the most successful scenario in the frameworks of the generalized model (7), (25). However, this generalized model has N p = 7 free parameters, including β and γ in addition to five parameters (28). This large number is the serious disadvantage of the generalized model (25) if we keep in mind informational criteria [18][19][20] and difficulties in calculations.
Following these reasons, we calculated the χ 2 tot function (29) for the generalized model (18), (25) with all H(z) data searching its minimum in the β − γ plane. The results are presented in Figure 3 as one-parameter distributions χ 2 tot (β) and χ 2 tot (γ) (minimized over all other parameters). We see in Figure 3 that the absolute minimum for the generalized model (25) min χ 2 tot ≃ 1084.0 is achieved near the point β = 2, γ = 0. One may conclude that the Barboza-Alcaniz-Zhu-Silva (BAZS) model [24] (corresponding to γ = 0) with Ω x (a) (21) and β = 2 appeared to be the most successful for the considered observational data Ia + CMB + CC + H BAO . When we fix β = 2, the BAZS model with w = w 0 + w 1 (1 − a 2 )/2 (32) will have five free model parameters (28). We investigate this model (denoted below as "β = 2") in detail; the results are presented in Table II and in Figure 4. Figure 4 illustrates the BAZS model with β = 2 (32): 1σ and 2σ contour plots are shown for all H(z) data (filled contours) and for only CC data. They are compared with the corresponding 1σ contours of the CPL model (shown also in Figure 2). One can compare the related one-parameter distributions  approximately coincide for all considered models; they depend on a chosen dataset: only CC data or all H(z) data.

Figure 4 and
The best fitted values of w 0 are close for the β = 2, CPL, and wCDM models in the CC case and slightly differ for all H(z) data. The optimal values of w 1 are larger for the β = 2 model (w 1 = 1.33 +0.42 −1.56 for CC and 1.12 +0.40 −1.09 for all H data) if we compare with the CPL model; however, one should take into account the factor 1 2 in EoS (32): w = w 0 + w 1 (1 − a 2 )/2. To compare models with different number N p of free model parameters, we use here the Akaike information criterion [18,19,59]: This criterion emphasizes the advantage of models with small number of N p . It can be seen in Table III for the mentioned models, where, for only CC H(z) data, the minimal Akaike values (33) are achieved for the ΛCDM (AIC = 1079.81) and the wCDM (AIC = 1080.97) models. However, for all H(z) data, the β = 2 model (32) with N p = 5 parameters appeared to be more successful than ΛCDM not only in min χ 2 tot ≃ 1084.03, but also with Akaike information (33): AIC = 1094.03. The lowest AIC = 1093.19 is achieved here for the wCDM model.
On the right side of the Table III, we present the predictions of the same models on their flat case k = 0 (Ω k = 0). In this case, the values min χ 2 tot and AIC appeared to be larger, and the β = 2 model (32) loses its advantage in comparison with CPL scenario for all H(z) data.

V. DISCUSSION
We considered different cosmological models with variable equations of state (EoS) for dark energy of the type (2) p x = w(a) ρ x , more precisely, models with EoS (7), generalizing the ΛCDM, wCDM, Chevallier-Polarski-Linder (CPL), and Barboza-Alcaniz-Zhu-Silva (BAZS) models [22][23][24]. These scenarios with nonzero spatial curvature and with Ω k = 0 were confronted with observational data described in Section II and including SNe Ia, CMB data, and two classes of the Hubble parameter estimates H(z): from cosmic chronometers (CC) and from line-of-sight baryonic acoustic oscillations (H BAO ) data. The results of our calculations for different models are presented in Tables II and III, including  minima of χ 2 tot , Akaike information criterion (33), and the best fitted values with 1σ estimates of model parameters. We also investigated the linear model (3), (22) with N p = 5 parameters (28), however it appeared to be unsuccessful because it achieved the best fitted value of χ 2 tot at w 1 very close to zero (see Figure 1). In other words, the linear model (22) with the considered observational data was reduced to the wCDM model with only N p = 4 parameters, but both models have the same min χ 2 tot . Unlike the linear model (22) CPL scenario (4), (23) with the same N p = 5 parameters (28) appeared to be more successful, in particular, for all H(z) data, the CPL model yields min χ 2 tot ≃ 1084.28 in comparison with 1085.19 for the wCDM and 1089.42 for the ΛCDM models. In this case, the best fitted value w 1 ≃ 0.62 is far from zero, so CPL is not reduced to the wCDM model.
We should remember that a large number N p of free model parameters is a drawback of any model, and when we use, in Table III, the Akaike information criterion (33), the wCDM model with AIC = 1093.19 will have the advantage over CPL with AIC = 1094.28 (and more essential advantage over the ΛCDM model with 1095.42). If we consider the generalized model (7), (25) with N p = 7 and additional parameters β and γ, the Akaike expression (33) becomes too large and the model looks worse in comparison with others.
However, our analysis and Figure 3 showed that the minimum min χ 2 tot ≃ 1084.03 of the generalized model (25) may be achieved if we fix β = 2, γ = 0; the resulting model "β = 2" (32) has the same N p = 5 parameters (28) and absolutely minimal min χ 2 tot for all H(z) data (its AIC = 1094.03 is behind only the wCDM AIC). The behavior of one-parameter distributions χ 2 tot (H 0 ) and their minima for all models with all H(z) data are shown in Figure 5.
As mentioned above, if we use only CC H(z) data, we observe smaller differences between minima of χ 2 tot for the considered models in Table II. In this case, the Akaike criterion (33) gives advantage to the ΛCDM model with minimal N p . Table II and Figures 2 and 4 demonstrate that the success of the CPL scenario is achieved at lower (best fitted) values of Ω 0 m and at larger values of Ω k , if we compare these results with the wCDM model. This tendency is strengthened for the β = 2 model (32). It looks natural, because the CPL model (23) is the particular case of the BAZS model (21) with β = 1. One can see in Tables II and III and Figures 2, 4, and 5 that the best fitted values of the Hubble parameter H 0 are very close for models wCDM, CPL, and β = 2 in general and spatially flat cases. However, the predictions of these models diverge with that of the ΛCDM model for all H(z) data. The last result is illustrated in Figure 5; it becomes more clear if we look at Figure 2, keeping in mind that the ΛCDM model is the particular case of the wCDM model when w = −1. One should add the observed difference between H 0 estimates of all models when we compare CC and all H(z) data.
In the left panels of Figure 5 we draw one-parameter distributions χ 2 tot (H 0 ) and likelihood functions (31) L(H 0 ) for the considered models with all H(z) data to clarify their best fits of H 0 . For the considered four models we demonstrate here the results for the flat cases k = 0 (Ω k = 0) as dotted lines.
One can see that in the flat case k = 0, the mentioned models yield appreciably larger minima min χ 2 tot , but the best fitted values of H 0 change unessentially for the wCDM, CPL, and β = 2 models. For the flat ΛCDM model, the graph χ 2 tot (H 0 ) is inscribed between the correspondent graphs of the general ΛCDM and flat wCDM models; hence, the H 0 estimate for the flat ΛCDM model lies between two estimates of the mentioned scenarios.
All H 0 predictions of the wCDM, CPL, and β = 2 models with all H(z) data fit the Planck 2018 estimation of the Hubble constant [2] (see Figure 5), but they are far from the SH0ES 2021 estimation [4]. All H 0 predictions with CC H(z) data are very close for the four considered models and their flat variants: they are larger but overlap the Planck 2018 value; however, they cannot describe the tension with the SH0ES data.
We can see from Table III that the best values of AIC (33) for all H(z) data demonstrate the wCDM and β = 2 (32) models, but for CC H(z) data, the best AIC has the flat ΛCDM model.