1 Introduction

Ever more precise observations of various cosmic phenomena [1,2,3,4,5,6] reveal a present state of universe which cannot be understood only in terms of General Relativity and forms of matter known from local physics, such as radiation or baryonic matter. Available observations point to the presently accelerated cosmic expansion, whereas the dynamics at the level of galaxies and clusters of galaxies, among other places, reveals additional gravitational interaction which could be explained by the presence of large quantities of, yet not directly observed, dark matter. The mechanism behind the accelerated cosmic expansion is usually attributed to a cosmic component with the negative pressure, called dark energy (DE). It has been shown that the concept of dark energy can be realized in many different ways such as cosmological constant, dynamical cosmological term, quintessence, phantom energy, k-essence or interacting dark energy [7,8,9,10,11,12]. Numerous alternatives have been proposed to both dark matter and dark energy, frequently as a modification of gravitational interaction at scales from galactic to cosmic  [13,14,15]. Yet, even if the effects such as accelerated cosmic expansion or galactic rotation curve dynamics do not originate from cosmic components, concepts of dark matter and dark energy (including their unifications) remain very useful effective concepts. Present observational data reveal a large tension in the value of \(H_0\) inferred from low and high redshift measurements assuming the benchmark \(\Lambda \)CDM model (for a recent review see  [16]). Some of proposed solutions to this puzzle are nontrivial dynamics of dark energy  [6] and dark matter-dark energy interaction  [17].

Presently, various models of dark energy have been proposed that available observational data cannot efficiently discriminate. Physically very distinct DE models can produce very similar global DE evolution and, correspondingly, very similar history of global cosmic expansion. Without a preferred DE model, the fits to observational data have to be performed for a large number of dynamically near-degenerate models. In such a situation a number of researchers have adopted a phenomenological approach of modeling the equation of state (EoS) parameter as a function of the scale factor, \(w=w(a)\) (or equivalently of the cosmic redshift, \(w=w(z)\))  [18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36,37,38,39,40]. This approach simplifies the analysis of DE dynamics and allows the analysis of physically interesting w(a) functions. Although such parametrizations may constitute a phenomenological approach of their own to the modeling of dark energy, their main purpose is the simplification of the fits to the observational data. In this way a single simple w(a) parametrization may represent the dynamical behavior of a large number of DE models. Yet, it is important to know to which extent the choice of some parametrization limits its representation by some specific physical model. In particular, it would be interesting to know if some DE models cannot be represented by some w(a) parametrization. In this paper we particularly focus on barotropic fluid models of dark energy and investigate which w(a) parametrizations are compatible with barotropic fluid DE. Here barotropic fluid is understood as a fluid for which the fluid pressure is a function of fluid energy density only. In determining the compatibility we do not employ the comparison of particular parametrizations with the observational data, but rely on fundamental physical constraints on the fluid DE speed of sound.

If we assume that the dark energy specified by some particular w(a) (or equivalently w(z)) parametrization is physically a barotropic fluid, an explicit expression for the barotropic speed of sound squared can be obtained from w(a). Inserting the definition of speed of sound squared for barotropic fluids, \(c_s^2 = \frac{d \, p}{d \, \rho }\) and the Equation of State (EoS) parameter \(w=\frac{p}{\rho }\) into the continuity equation

$$\begin{aligned} d \rho + 3 (\rho +p) \frac{d a}{a} =0, \end{aligned}$$
(1)

yields the dynamical equation for the EoS parameter

$$\begin{aligned} a \frac{d w}{ d a} = -3 (1+w)(c_s^2-w). \end{aligned}$$
(2)

This equation can be easily rearranged to obtain the expression for \(c_s^2\) in terms of w and \(a \frac{d w}{d a}\):

$$\begin{aligned} c_s^2=w - \frac{1}{3(1+w)} a \frac{d w}{d a}. \end{aligned}$$
(3)

If \(a \frac{d w}{d a}\) can be expressed as a function of w, then the speed of sound squared can also be expressed as a function of w, i.e. \(c_s^2=c_s^2(w)\). This line of modeling has been successfully applied to the description of cosmological constant boundary crossing  [41] and dark energy-dark matter unification [42, 43].

The barotropic fluid speed of sound squared is physically constrained to be nonnegative and not larger than speed of light squared, \(c^2\). As we work in system of units where \(c=1\), these requirements translate to \(0 \le c_s^2 \le 1\). For known w(a) one can obtain \(c_s^2(a)\) from (3) and model parameters for which \(0 \le c_s^2(a) \le 1\) is satisfied for the entire past cosmic expansion, i.e. for the entire \([0,a_0]\) interval. In this way we can select w(a) parametrizations which are suitable for the description of barotropic fluid dark energy as those for which the condition on speed of sound squared is satisfied at least for some model parameters. The allowed region of model parameters is further analyzed if some of its portion corresponds to presently accelerating component (corresponding to \(\rho +3 p <0\)). This program, though physically simple, turns out to be quite restrictive for a large number of DE parametrization models.

The paper is organized as follows. The first section brings the introduction and the presentation of the main idea. In the second section we present analytical approach to determination of allowed model parameters and apply it to a one-parameter model and elaborate general methods useful in analytical treatment. In the third section we present numerical approaches to determination of allowed parameter values and apply them to a large number of parametrizations available in the literature. In the following section we discuss the obtained results and finish the paper with conclusions. In the “Appendix” we bring the analytical solution for the Chevallier–Polarski–Linder (CPL) model  [20, 21].

2 Analytical results

The feasibility of constraining the model parameters analytically crucially depends on the form of the w(a) function and constraints can be obtained analytically only in specific cases. Even in cases where the said constraints can be obtained using analytical techniques, the very procedure can be quite involved and the obtained results are not very transparent and informative. Still, analytically tractable cases can be very useful for the verification of more generally applicable numerical approaches and they can provide additional insights that numerical approaches do not provide. As an illustration we describe the analytical procedure for obtaining parameter constraints for a one-parameter model [25] and a more general approach suitable for two-parameter models such as the CPL model  [20, 21].

2.1 \(w=w_0 \frac{a}{a_0}\) model

Starting from the parametrization  [25]

$$\begin{aligned} w=w_0 \frac{a}{a_0}, \end{aligned}$$
(4)

from (3) one obtains

$$\begin{aligned} c_s^2=w_0 \frac{a}{a_0} \frac{2 + 3 w_0 \frac{a}{a_0} }{3(1 + w_0 \frac{a}{a_0})} = w \frac{2+3w}{3(1+w)} = w \left( 1-\frac{1}{3(1+w)} \right) . \end{aligned}$$
(5)

The first and second derivative of \(c_s^2\) with respect to w are:

$$\begin{aligned} \frac{d c_s^2}{d w}= & {} 1 - \frac{1}{3(1+w)^2}, \end{aligned}$$
(6)
$$\begin{aligned} \frac{d^2 c_s^2}{d w^2}= & {} \frac{2}{3} \frac{1}{(1+w)^3}. \end{aligned}$$
(7)

Stationary points of \(c_s^2\) are at \(w_1=-1-\frac{1}{\sqrt{3}}\) and \(w_2=-1+\frac{1}{\sqrt{3}}\) where at \(w_1\) the \(c_s^2\) has a maximum and at \(w_2\) it has a minimum. Both of these stationary points correspond to negative values of w. At \(w=-1\) there is a singularity in \(c_s^2\) corresponding to the crossing of the cosmological constant boundary.

During the cosmic expansion the EoS parameter (4) does not change its sign, i.e. as a increases from 0 to \(a_0\), w changes from 0 to \(w_0\) (increases for positive \(w_0\) and decreases for negative \(w_0\)).

For negative \(w_0\) (negative vaules of w), expression (5) reveals that \(c_s^2\) is negative for \(-\frac{2}{3}<w<0\). Therefore, if \(w_0 > -\frac{2}{3}\), \(c_s^2\) is negative during the entire cosmic past, whereas if \(w_0 < -\frac{2}{3}\), \(c_s^2\) is negative for a from 0 to some finite \(a_* < a_0\). In both cases the condition \(c_s^2 \ge 0\) is violated in the cosmic past and \(w_0 < 0\) does not correspond to viable fluid model of dark energy.

For positive \(w_0\), there are no stationary points in the interval of w between 0 and \(w_0\). The expression (6) reveals that \(\frac{d c_s^2}{d w}\) is positive in the interval \((0,w_0)\) and \(c_s^2\) is a growing function of w. The expression (5) shows that \(c_s^2\) is positive in the entire considered interval. Therefore, to fulfill the requirement \( 0 \le c_s^2 \le 1\) at the entire interval, it suffices to require that \(c_s^2(w_0) \le 1\). Straightforward calculation shows that this is satisfied for

$$\begin{aligned} 0 \le w_0 \le \frac{1+\sqrt{37}}{6}. \end{aligned}$$
(8)

As this allowed paramter range corresponds to \(w \ge 0\), the parametrization from  [25] is clearly unsuitable as a model of barotropic fluid dark energy.

2.2 General analytical approach

From the condition \(0 \le c_s^2 \le 1\) it is possible to obtain general analytical constraints between \(a \frac{d w}{d a}\) and w [44]. In particular, for \(1+w >0\) the condition translates to

$$\begin{aligned} -(1-w)(1+w) \le \frac{a}{3} \frac{d w }{d a} \le w(1+w), \end{aligned}$$
(9)

whereas for \(1+w < 0\), the condition results in

$$\begin{aligned} w(1+w) \le \frac{a}{3} \frac{d w }{d a} \le -(1-w)(1+w). \end{aligned}$$
(10)
Table 1 Overview of studied models, their w(a) (or w(z)) parametrizations and allowed regions of parameters
Fig. 1
figure 1

The allowed parameter regions for two-parameter models  [20, 21] in plot (a),  [23] (model 1) in plot (b),  [23] (model 2) in plot (c), [37] (model 1) in plot (d),  [28] in plot (e) and [40] in plot (f). In all plots the symbol \(w_0\) on the axis denotes the present value of the w(a) function, whereas the other symbols refer to parameters in the corresponding w(a) parametrizations

Fig. 2
figure 2

The allowed parameter regions for two-parameter models [30] in plot (a),  [32] in plot (b),  [33] (model 1) in plot (c),  [33] (model 2) in plot (d),  [36] in plot (e) and three-parameter model  [27] for \(a_*=a_0/3\) in plot (f). In all plots the symbol \(w_0\) on the axis denotes the present value of the w(a) function, whereas the other symbols refer to parameters in the corresponding w(a) parametrizations

Fig. 3
figure 3

The allowed parameter regions for three-parameter models  [29] (model 1) for \(n=3\) in plot (a),  [37] (model 2) and  [29] for \(n=3\) in plot (b) and  [35] for \(w_0=-0.6\) in plot (c) and four-parameter models  [22] for \(p=1\), \(w_0 = -0.6\) and \((\frac{a_s}{a_0})^p = \frac{w_b}{w_a} \frac{w_a - w_0}{w_0 - w_b}\) in plot (d),   [24] for \(p=1\), \(w_0 = -0.6\) and \((\frac{a_c}{a_0})^p = \frac{w_a w_b - w_0}{w_0 - w_a}\) in plot (e) and  [39] \(w_0=-0.6\) and \(w_b=-1\) in plot (f). In all plots the symbol \(w_0\) on the axis denotes the present value of the w(a) function, whereas the other symbols refer to parameters in the corresponding w(a) parametrizations

Fig. 4
figure 4

The allowed paremeter regions for conditions \(0 \le c_s^2 \le 1\) (dots) and \(0 \le c_s^2\) (dots + crosses) for models: [20, 21] in plot (a),  [27] for \(a_* = a_0 / 3\) in plot (b),  [30] in plot (c) and  [37] (model 1) for \(n = 3\) in plot (d). In all plots the symbol \(w_0\) on the axis denotes the present value of the w(a) function, whereas the other symbols refer to parameters in the corresponding w(a) parametrizations

Furthermore, it is straightforward to show that the condition \(0 \le c_s^2 \le 1\) is equivalent to

$$\begin{aligned} c_s^2 (c_s^2-1) \le 0, \end{aligned}$$
(11)

which can also be presented in the form

$$\begin{aligned} \frac{f(a) g(a)}{9 (1+w)^2} \le 0. \end{aligned}$$
(12)

Here, bearing in mind that \(w=w(a)\),

$$\begin{aligned} f(a)= 3 w (1+w) - a \frac{d w}{d a}, \end{aligned}$$
(13)

and

$$\begin{aligned} g(a) = 3 (w-1) (1+w) - a \frac{d w}{d a}. \end{aligned}$$
(14)

From (12) one can determine the regions of allowed model parameters as those for which \(f(a) g(a) \le 0\) for all \( a \in [0, a_0]\). In the parametrizations for which \( a \frac{d w}{d a}\) can be expressed as a function of w, this condition then reads as \(f(w) g(w) \le 0\). For a majority of parametrizations the determination of the regions of allowed model parameters cannot be pursued analitically. In the “Appendix” we systematically apply this approach for the CPL model \(w=w_0+w_1(1-\frac{a}{a_0})\), introduced in  [20, 21]. Although in this case both f and g are (only) quadratic functions of w, the analytical calculations require examination of a number of various cases.

3 Numerical results

The compatibility of a particular parametrization with the barotropic fluid DE can in general be established only numerically. For the studied parametrizations we determine the allowed region of the model parametric space using two approaches and present the results in Table 1. For models which have a nonvanishing parameter region corresponding to \(0 \le c_s^2 \le 1\) for the entire cosmic past these regions are depicted in Figs. 1, 2 and 3.

The graphs in Figs. 1, 2 and 3 were made combining two methods:

  1. 1.

    Shaded areas: An analytical solution was rearranged to put one parameter on each axis (x and y) and solution space for discrete values of a was graphed (0, \(a_0\) and one or two points in the middle; usually, but not always, \(a/a_0=0.5\), depending which value created a better illustration of the effect). The intersection of these three (or four) areas approximates a solution for the whole range of \(a \in [0, a_0 ]\).Footnote 1

  2. 2.

    Dots: A numerical solution was computed with two parameters laid on x and y axes. The condition \(0 \le c_s^2 \le 1\) was tested at a set of scale factor a values, where the values in the set were chosen for each model to cover the entire cosmic past, but also to produce the best coverage of the allowed parameter region. In most cases, one hundred values of a were calculated using \(a = e^{i/10} a_0\) with i taking integer values from -100 to 0. This produced an array of values for a more dense near \(a = 0\) and more spread out near \(a = a_0\). A point was placed on the graph for all values of the two parameters where \(c_s^2 (a) \in [0, 1] \) was true for all values of \(a \in [0, a_0]\).

Following the initial test for \(c_s^2 (a) \in [0, 1] \), this condition was loosened to only \(c_s^2 (a) \ge 0\) without the \(c_s^2 (a) \le 1\) condition. As expected, the \((w_0, w_1)\) space expanded with the following restrictions:

  1. 1.

    If the favourable area under the condition \(c_s^2 (a) \in [0, 1] \) existed only for values \(w_0 > 0\), the expansion happened solely in the direction \(w_0 > 0\).

  2. 2.

    If the favourable area existed for values \(w_0 < 0 \), but did not reach \(w_0 = -1\), it expanded towards \(w_0 = -1\), but did not cross to \(w_0 < -1\).

For a number of parametrizations the allowed parameter regions for the condition \(c_s^2 \ge 0\) are presented in Fig. 4.

Plots in Figs. 1 and 2 present the two-parameter models with nonvanishing allowed regions of parameters. In all these plots the parameter at the x axis is \(w_0\), corresponding to the present value of the EoS parameter \(w(a_0)\) (or \(w(z=0)\)). It is interesting to observe which models allow \(w_0 < -1/3\) values (which can in principle serve as accelerating components). In Fig. 1 this condition is satisfied for models in plots (a)  [20, 21], (d) [37] and (e)  [28]. In Fig. 2 presently accelerating component is possible for the model [32] in plot (b) and the model  [27] in plot (f).

For three-parameter and four-parameter models, presented in Fig. 3, all models presented may describe barotropic fluid dark energy. The plots demonstrate that for the selected values of \(w_0=w(a_0)\) there are nonvanishing allowed regions of other model parameters. This fact shows that the capability of these models to represent baryonic DE model is not the result of some contrived combinations of model parameters, but a generic feature of these models.

4 Discussion and conclusions

In total, of the three one-parameter models studied in this paper, two have a nonvanishing parameter space consistent with the requirement \(0 \le c_s^2 \le 1\), but none of the allowed parameters corresponds to the presently accelerating component. For the 16 two-parameter models we analyzed, 11 of them have a nonvanishing parameter space consistent with \(0 \le c_s^2 \le 1\), but only four of these can describe a presently accelerating component. For ten of the three-parameter and four-parameter models, seven of them satisfy the requirement \(0 \le c_s^2 \le 1\) and are compatible with the barotropic fluid dark energy.

One can observe that virtually none of studied models allows values of \(w_0\) very close to \(-1\). This fact may be attributed to the term \(\sim \frac{1}{1+w}\) in the expression (3) for \(c_s^2\) which diverges when \(w \rightarrow -1\), which correspondingly blows up the value of \(c_s^2\). A natural question arising is which of two requirements (\(c_s^2 \ge 0\) or \(c_s^2 \le 1\)) is responsible for such a behavior. Indeed, if the \(c_s \le 1\) requirement is relaxed, the allowed values of \(w_0\) are much closer to \(-1\), as is evident from Fig. 4.

As one could expect, the larger the number of model parameters, it is easier to find their combination for which the model may be represented as a barotropic fluid DE in the entire cosmic past. For the studied single parameter models it is found that none of them can satisfy the requirement on \(c_s^2\) and presently have a sufficiently negative w. For two-parameter models, only four out of 16 models are capable of representing the barotropic fluid dark energy. For three-parameter and four-parameter models seven of the ten studied models fit the requirement of baryonic fluid dark energy.

It is important to notice (possibly even somewhat surprising) that these strong restrictions have been obtained on purely theoretical, but fundamental grounds. The parametrizations w(a) that are found to be able to represent a barotropic fluid DE by the studied requirements on \(c_s^2\) still need to be compared against the available observational data which will further constrain the parametric space obtained in this paper.

The results of this paper indicate that the suitability of a phenomenological parametrization w(a) to describe some physically motivated dark energy model need not come automatically. Internal theoretical features of a chosen physical DE models may constrain, or even eliminate some phenomenological parametrizations. Possible situations in which this kind of argumentation might be applicable, apart from the barotropic fluid DE models, comprise k-essence DE models or effective DE description of modified gravity theories.