Non-linear causal bulk viscosity in Unified Dark Matter Cosmologies

TL;DR

The paper develops a fully causal, nonlinear Müller–Israel–Stewart bulk-viscosity model for a unified dark matter fluid in a flat FRW universe with radiation. By treating the bulk-viscosity exponent s in ξ ∝ ρ_m^s as a free parameter and recasting the dynamics as an autonomous dynamical system, it reveals three late-time attractor classes: de Sitter for s<1/2, phantom for s>1/2, and a coexistence of de Sitter/quintessence-like behaviour at s=1/2. Numerical evolution confirms theoretical fixed-point predictions and shows that acceleration can arise without imposing large ξ_0, with the non-linear parameters k and v shaping the detailed evolution and the ultimate fate (de Sitter or phantom). The study also highlights meaningful differences from linear Israel–Stewart theory, emphasizing the importance of far-from-equilibrium thermodynamics in cosmology and offering a potential unified mechanism for late-time acceleration within the dark sector.

Abstract

We propose a bulk viscous unified dark matter scenario based on a nonlinear extension of the full causal Israel-Stewart theory. This framework allows the viscous fluid to remain far from equilibrium, an essential feature for a physically consistent description of viscosity-driven accelerated expansion. We adopt the standard parametrization for the bulk viscosity, $ξ= ξ_{0} ρ_{m}^{s}$, treating $s$ as a free parameter (in contrast to most previous works), and study the model in a spatially flat Friedmann-Robertson-Walker background. By reformulating the cosmological equations as an autonomous dynamical system, we obtain both asymptotic analytical solutions and a numerical characterization of the phase space. At early times, the viscous component can mimic a stiff fluid, while at intermediate epochs it behaves like dark matter. With a suitable choice of dynamical variables, the system admits three distinct classes of late-time attractors. Two of them are separated by a basin-boundary saddle point: (i) a generic quintessence solution for $s = 1/2$, which encompasses a de Sitter-like behavior when $ξ_{0}$ satisfies a specific relation involving the nonlinear parameters; (ii) a global exact de Sitter attractor for $s < 1/2$; and (iii) a phantom-like solution that emerges for $s \ge 1/2$. In contrast to the generic $s \ne 1/2$ case, the $s = 1/2$ scenario exhibits a qualitatively different stability structure, allowing de Sitter and phantom attractors to coexist. All solutions respect entropy production, and cosmic acceleration emerges independently of $ξ_{0}$, relaxing the strong bounds $ξ_{0} \sim \mathcal{O}(1)$ required in Eckart-based viscous models.

Figures (7)

• Figure 1: Numerical trajectories in phase space with several initial conditions that capture the physical behavior of the dynamical system for $s\neq 1/2$. This evidences the attractor nature of the fixed point $P_{Ia}$. Initial conditions are taken close to the radiation point $P_{II}^{+}$. The latter is bridged to its negative branch $P_{II}^{-}$, as evidenced. Thus, trajectories starting strictly on $P_{II}^{+}$ will end up on $P_{II}^{-}$ and, subsequently, on the attractor point $P_{Ia}$. The point $P_{III}^{+}$ is a repeller, describing a stiff fluid at early times, while its negative branch, $P_{III}^{-}$, is a phantom fixed point. This cannot be visualized here due to the resolution scale and is instead represented in Figure \ref{['fig:phase_space__sfree_phantom']}. The point $P_{IV}$, corresponding to dark matter with $w_{\rm eff}=0$, is the only one that responds sensitively to changes in the bulk viscosity coefficient $\xi_{0}$. The blue region corresponds to the physical region where acceleration can take place, along with the condition $U > 0$.
• Figure 2: Phase-space trajectories for the dynamical system with $s = 0.6$. Trajectories departing from $P_{II}^{-}$ evolve either toward the de Sitter (saddle) point $P_{Ia}$ (blue and red curves) or directly toward the effective phantom attractor $P_{III}^{-}$, characterized by $w_{\Pi} = -\sqrt{2} \, \gamma \, v \approx -1.00409$ with $v = 0.71$ (magenta and black curves), depending on the initial conditions. The basin-boundary saddle point $P_{II}^{-}$ separates these two intermediate paths. However, since $P_{Ia}$ is a saddle for $s > 1/2$, all trajectories ultimately converge to the phantom attractor, which is the only late-time stable solution. The blue-shaded region indicates the physically allowed domain with $U>0$, where cosmic acceleration can occur.
• Figure 3: Evolution of the viscous DM EoS parameter, $w_{\Pi}$, and the overall effective EoS, $w_{\rm eff}$, as a function of $N$. We adopt $\xi_{0} = 10^{-3}$, for the case $s<1/2$, while varying the viscous exponent $s$ to illustrate its impact on the global dynamics: $s = 0.4$ (top-left panel), $s = 0$ (top-right panel), and $s = -1/2$ (bottom-left panel). For these three cases, we have chosen the corresponding set of initial conditions: $U(t_{\rm init}) = 0.5$, $V(t_{\rm init}) = 0.99$, $W(t_{\rm init})= 0.2$ for red curves; $U(t_{\rm init}) = 0.7$, $V(t_{\rm init}) = 0.99$, $W(t_{\rm init}) = 0.01$ for black curves; $U(t_{\rm init}) = 0.3$, $V(t_{\rm init}) = 0.01$, $W(t_{\rm init}) = 0.002$ for blue curves. Decreasing $s$ disrupts the smooth late-time evolution, hastening the onset of accelerated expansion. In the bottom-right panel, we consider $s > 1/2$ ($s=0.6$) with $\xi_{0} = 0.1$ to illustrate the phantom-like behavior, as highlighted in the inset (zoom-in) panel. For this plot, we have chosen, instead, the following set of initial conditions: $U(t_{\rm init}) = 10^{-6}$, $V(t_{\rm init}) = \sqrt{2} v$, $W(t_{\rm init})= 0.001$ for red curves; $U(t_{\rm init}) = 0$, $V(t_{\rm init}) = 0.98$, $W(t_{\rm init}) = 0.98$ for black curves; $U(t_{\rm init}) = 0$, $V(t_{\rm init}) = -\sqrt{2}v$, $W(t_{\rm init}) = -0.001$ for magenta curves. In all plots we have adopted $v = 0.71$, and $k = 0.3$, as references values.
• Figure 4: Evolution of the ratio $|\Pi|/H^{2}$, which serves as an indicator of the relative contribution of viscous effects during cosmic expansion. We consider the case $s=0.4$ with $\xi_{0}=10^{-3}$ and explore three distinct initial conditions: IC1: $U(t_{\rm init}) = 0.5$, $V(t_{\rm init}) = 0.99$, $W(t_{\rm init}) = 0.2$; IC2: $U(t_{\rm init}) = 0.7$, $V(t_{\rm init}) = 0.99$, $W(t_{\rm init}) = 0.01$; IC3: $U(t_{\rm init}) = 0.3$, $V(t_{\rm init}) = 0.01$, $W(t_{\rm init}) = 0.002$. The solid red curve (IC1) illustrates a scenario in which the effective dark matter pressure behaves as a stiff fluid due to the dominant contribution of the viscous pressure at early times. The two dashed curves (IC2 and IC3) correspond to cases in which viscosity becomes dynamically relevant only at late times, driving the accelerated expansion phase where $|\Pi|/p_{T} = 1$.
• Figure 5: Phase space portrait of the dynamical system for $s = 1/2$. Fixed points appear distributed nearly symmetrically with respect to $Y=0$, with similar dynamical character. For $Y = \pm 1$, three fixed points exist on each branch, labeled $P_{I}^{\pm}$. Among them, $P_{Ia}^{+}$, $P_{Ib}^{+}$, and $P_{Ib}^{-}$ all act as attractors within the physical region (light red). Two saddle points arises on each branch: $P_{Ic}^{+}$ and $P_{Ic}^{-}$. At $Y = 0$, two radiation-like points with $V_{II} = \pm \sqrt{2}v$ are present ($P_{II}^{\pm}$). Numerical trajectories illustrate the global behavior of the system. In the left panel, we adopt $v = 0.71$, $k = 0.3$, and $\xi_{0} = 0.9$, while in the right panel, we use $v = 0.95$, $k = 0.3$, and $\xi_{0} = 2.49$.