# Shan–Chen interacting vacuum cosmology

Natalie B. Hogg1,2 and Marco Bruni2,3
1Instituto de Física Teórica UAM/CSIC, C/Nicolás Cabrera 13-15, Universidad
Autónoma de Madrid, Cantoblanco, Madrid 28049, Spain
2Institute of Cosmology and Gravitation, University of Portsmouth, Burnaby
Road, Portsmouth PO1 3FX, United Kingdom
3INFN Sezione di Trieste, Via Valerio 2, 34127 Trieste, Italy E-mail:
<EMAIL_ADDRESS><EMAIL_ADDRESS>
(Accepted 2022 February 2. Received 2022 February 2; in original form 2021
September 23)

###### Abstract

In this paper, we introduce a novel class of interacting vacuum models, based
on recasting the equation of state originally developed in the context of
lattice kinetic theory by Shan & Chen as the coupling between the vacuum and
cold dark matter (CDM). This coupling allows the vacuum to evolve and is
nonlinear around a characteristic energy scale $\rho_{*}$, changing into a
linear coupling with a typical power law evolution at scales much lower and
much higher than $\rho_{*}$. Focusing on the simplest sub-class of models
where the interaction consists only of an energy exchange and the CDM remains
geodesic, we first illustrate the various possible models that can arise from
the Shan–Chen coupling, with several different behaviours at both early and
late times depending on the values of the model parameters selected. We then
place the first observational constraints on this Shan–Chen interacting vacuum
scenario, performing an MCMC analysis to find those values of the model and
cosmological parameters which are favoured by observational data. We focus on
models where the nonlinearity of the coupling is relevant at late times,
choosing for the reference energy scale $\rho_{*}$, the critical energy
density in $\Lambda$CDM. We show that the observational data we use are
compatible with a wide range of models which result in different cosmologies.
However, we also show that $\Lambda$CDM is preferred over all of the Shan–Chen
interacting vacuum models that we study, and comment on the inability of these
models to relax the $H_{0}$ and $\sigma_{8}$ tensions.

###### keywords:

Cosmology – dark energy – dark matter

††pubyear: 2022††pagerange: Shan–Chen interacting vacuum cosmology–11

##  1 Introduction

Reports of the death of $\Lambda$CDM are greatly exaggerated. Since the
discovery of the accelerating expansion of spacetime more than twenty years
ago (Riess et al., 1998; Perlmutter et al., 1999), numerous observations of
our Universe have yielded a standard model of cosmology that has stood up to
every challenge. This model, $\Lambda$CDM, posits that the Universe is
dominated at late times by a cosmological constant, $\Lambda$, which is
responsible for the observed accelerating expansion. Myriad observations of
structure formation and galaxy evolution (but not experimental particle
physics – yet) have provided robust evidence for the existence of the other
major component of our Universe, cold dark matter (CDM) (Zwicky, 1933; Eke et
al., 1998; Sofue & Rubin, 2001; Frenk & White, 2012; Kunz et al., 2016).

Nevertheless, problems with the cosmological constant and its use as the dark
energy abound (Martin, 2012), leading to many attempts to explain the late
time acceleration with an alternative model. Examples of alternative dark
energy models range from those in which the vacuum interacts with cold dark
matter (Wands et al., 2012; Salvatelli et al., 2014; Martinelli et al., 2019;
Hogg et al., 2020), to quintessence (Amendola, 2000; Linder, 2007; Park et
al., 2021), k-essence (Armendariz-Picon et al., 2001; de Putter & Linder,
2007) and yet more esoteric models (see e.g. Copeland et al. (2006) or Di
Valentino et al. (2021a) for a review of dark energy models).

Some alternative models of dark energy to the cosmological constant are also
motivated by their potential ability to reconcile the more than $4\sigma$
discrepancy between the value of the Hubble parameter at redshift zero
($H_{0}$) determined from cosmic microwave background (CMB) temperature and
polarisation measurements (see e.g. Aghanim et al. (2020); Di Valentino et al.
(2021b)) and that measured using the low redshift standardisable candles Type
Ia supernovae (SNIa) calibrated using Cepheid variable stars by the S$H_{0}$ES
collaboration (see e.g. Riess et al. (2019)). We note that when comparing
these two measurements alone, it is preferable to refer to this as a
“supernova absolute magnitude tension” rather than the more commonly used
“$H_{0}$ tension” – see Camarena & Marra (2021); Efstathiou (2021) for more
details on this point. For a different perspective on the discrepancy, see
Haridasu et al. (2021).

However, the S$H_{0}$ES SNIa are not the only local probe of $H_{0}$. Taking
into account the numerous complementary measurements of $H_{0}$ now available,
such as those from water masers (Pesce et al., 2020), SNIa calibrated using
stars at the tip of the red giant branch of the Hertzsprung–Russell diagram
(Sakai, 1999; Freedman et al., 2019) and strong lensing time delays (Wong et
al., 2020), the familiar tension between local determinations of $H_{0}$ and
measurements from the CMB which rely on a cosmological model remains. The
“$H_{0}$ tension” is thus rightly used by many authors as a motivation to
study beyond $\Lambda$CDM cosmologies, most commonly exploring different dark
energy models.

In this context, it is worth noting that in single scalar field quintessence
models, where the dark energy equation of state is $w\geq-1$, $H_{0}$ is
always smaller relative to the value obtained in a $\Lambda$CDM universe
(Banerjee et al., 2021), implying that the $H_{0}$ tension cannot be resolved
in these models.

A consideration sometimes made when searching for alternative dark energy
models is that the model should preferably have arisen from pre-existing
physics, and hence avoid having been overly designed to solve a specific
problem. More generally, this aesthetic standpoint can be seen as a facet of
the Bayesian school of statistical thought: “naturalness”, or our distaste for
finely tuned models which produce a desired phenomenology, is simply another
prior imposed on the model.

One promising approach along these lines is to consider dark energy and dark
matter as a single cosmological fluid, such as a unified dark matter, or UDM,
which behaves according to some equation of state. This simple idea has a long
history, dating back to the introduction of the Chaplygin gas model by
Kamenshchik et al. (2001), in which the transition from cold dark matter to
dark energy domination is achieved using a perfect fluid with a non-standard
equation of state.

However, serious objections to the Chaplygin gas model were raised when it was
found that unified models of this type result in oscillations or an
exponential blow-up in the matter power spectrum, thus ruling out the vast
majority of the viable Chaplygin gas model space (Sandvik et al., 2004). To
counter this, further models in which dark matter transitions to dark energy
via a condensation mechanism have been proposed and studied, and interest in
unified models of this type remains high – see, for example, Balbi et al.
(2007); Pietrobon et al. (2008); Piattella et al. (2010); Piattella (2010);
Gao et al. (2010); Bertacca et al. (2011); De Felice et al. (2012); Wang et
al. (2013); Li et al. (2018); Li et al. (2019). A review of fluid dark energy
models can be found in Bamba et al. (2012).

The Shan–Chen equation of state was first proposed by Shan and Chen in the
context of lattice kinetic theory in Shan & Chen (1993). A fluid with this
equation of state behaves as an ideal gas in the low and high density regimes,
but has a liquid–gas coexistence curve i.e. a region of temperature and
pressure in which the fluid can be in both the liquid and gas state. Such an
equation of state means that a phase transition can easily arise.

This model was successfully applied to the cosmological context by Bini et al.
(2013, 2016), with the finding that a dark energy fluid with a Shan–Chen
equation of state naturally evolves towards a Universe with a late time
accelerating expansion without the presence of a cosmological constant. In
those works, modifications made to the background expansion by the Shan–Chen
model were studied, and quantities such as the distance modulus in the
Shan–Chen model were compared to Type Ia supernova data.

In this paper, we cast the Shan–Chen equation of state into the form of a
coupling between the vacuum and cold dark matter. After studying this
Shan–Chen coupling in detail, we implement the interacting model in the public
codes CAMB (Lewis et al., 2000; Howlett et al., 2012) and CosmoMC (Lewis &
Bridle, 2002; Lewis, 2013) in order to obtain constraints on the cosmological
parameters and the model parameters in the Shan–Chen interacting vacuum model.

The structure of this paper is as follows: in Section 2 we outline the theory
of the Shan–Chen interacting vacuum model, in Section 3 we describe our
methodology, in Section 4 we present the results of our analysis, and in
Section 5 we discuss our results and conclude the paper.

Throughout this work, we assume General Relativity, use units in which the
speed of light $c=1$, and assume a spatially flat
Friedmann–Lemaître–Robertson–Walker (FLRW) background spacetime.

##  2 Interacting vacuum Shan–Chen model

###  2.1 Preliminary: the Shan–Chen equation of state

The nonlinear equation of state proposed in Shan & Chen (1993) is

$\displaystyle P$
$\displaystyle=w\rho_{*}\left[\frac{\rho}{\rho_{*}}+\frac{g}{2}\psi^{2}\right],$
(1) where $\rho$ is the energy density of the fluid, and $\displaystyle\psi$
$\displaystyle=1-e^{-\alpha\frac{\rho}{\rho_{*}}},$ (2)

where $\rho_{*}$ is a characteristic energy scale and $w$, $g$ and $\alpha$
are free (dimensionless) parameters of the model.

In the Shan–Chen model of dark energy introduced by Bini et al. (2013, 2016),
the matter–energy content of the Universe is assumed to be a perfect fluid
which obeys the equation of state (1), and $\rho_{*}$ is identified with the
critical density today:
$\rho_{*}=\rho_{\mathrm{crit},0}=\frac{3H_{0}^{2}}{8\pi G}$, where $H_{0}$ is
the value of the Hubble parameter at redshift zero.

We note that another choice for this energy scale could be, for example, the
value of the matter density at matter–radiation equality, as this could result
in an early dark energy-type behaviour. An example of a successful early dark
energy model can be found in Poulin et al. (2019), although that model’s
competitiveness with $\Lambda$CDM and its ability to resolve the $H_{0}$
tension was called into question in Hill et al. (2020), where the analysis of
Poulin et al. (2019) was updated using Planck 2018 data as well as additional
information from large scale structure. There have since been a number of
other works exploring different possibilities under the early dark energy
umbrella, including New Early Dark Energy (Niedermann & Sloth, 2020), Chain
Early Dark Energy (Freese & Winkler, 2021) and Chameleon Early Dark Energy
(Karwal et al., 2021). We leave the exploration of the Shan–Chen dark energy
scenario, including early dark energy models, to a future work.

While we do not need to go into the details of the Shan–Chen dark energy model
here, it is useful to point out few features for the sake of comparison with
the Shan–Chen interacting vacuum model we propose in the next section.

To this end, let us consider the energy conservation equation for a Shan–Chen
dark energy in an FLRW universe,

$\displaystyle\dot{\rho}=-3H(\rho+P(\rho)),$ (3)

where $P(\rho)$ is given by (1) and as usual $H=\frac{\dot{a}}{a}$ is the
Hubble expansion scalar, with $a$ being the cosmic scale factor.

We can then define the effective equation of state parameter
$w_{\mathrm{eff}}$,

$\displaystyle
w_{\mathrm{eff}}=\frac{P(\rho)}{\rho}=w+\frac{wg\rho_{*}}{2\rho}(1-e^{-\alpha\frac{\rho}{\rho_{*}}})^{2}\;.$
(4)

Clearly, this only depends on $x=\rho/\rho_{*}$. Although the energy scale
$\rho_{*}$ is physically meaningful, in that it determines at which energies
the nonlinear term in (4) becomes relevant, it has no influence on the
qualitative behaviour of the Shan–Chen effective equation of state.

Depending on the values of the model parameters $\alpha$, $g$ and $w$, there
may be stationary points $\rho=\rho_{\mathrm{S}}$ such that
$P(\rho_{\mathrm{S}})=-\rho_{\mathrm{S}}$, i.e. $w_{\mathrm{eff}}=-1$ and
$\dot{\rho}=0$, from (3). If they exist, in general they appear in
pairs111There can be none or two, or one when the two coincide.: they
represent cosmological constants, one stable and one unstable, and the dark
energy density can either evolve asymptotically toward the stable one, away
from the unstable one, or in between the two.

Around $\rho_{*}$ the fluid is not ideal, while in the limits
$\rho\gg\rho_{*}$ and $\rho\ll\rho_{*}$ we get $P=w\rho$ and the usual linear
equation of state is recovered. It is easy to show that $w_{\mathrm{eff}}$ has
a minimum when $wg<0$ and a maximum when $wg>0$; in addition, the position of
the minimum or maximum, given by a specific value of $x=\rho/\rho_{*}$, does
not depend on $w$ and $g$, but is shifted by $\alpha$.

In the top left panel of Figure 1, we show four examples of the effect of
different parameter values on the behaviour of $w_{\mathrm{eff}}$ as a
function of $x=\rho/\rho_{*}$: $w=\pm 1/3$ and $g=\pm 8$, keeping a fixed
value of $\alpha=2.7$, the best fit value of this parameter found in Bini et
al. (2013) in an analysis of SNIa distance modulus data. Note that there are
six general cases for $wg<0$ and six cases for $wg>0$, plus limiting cases. We
have chosen to show four examples only to ensure clarity in the plot. We
discuss them in more depth in Section 2.3 for easier comparison with the
interacting Shan–Chen scenario, but we will give a brief overview here.

Since $w_{\mathrm{eff}}$ is asymptotic to $w$, changing the value of $w$
whilst keeping the other parameters fixed simply shifts the curves up and
down, whilst changing the value of $g$ makes the “bump” in the curves more or
less pronounced. In general, the curves may or may not cross
$w_{\mathrm{eff}}=-1/3$ and $w_{\mathrm{eff}}=-1$. We can infer from this plot
that whether or not acceleration occurs in the Shan–Chen dark energy model,
and when, depends on the values of $w$, $g$ and $\alpha$. Note that for the
purposes of this plot, and for simplicity’s sake, we are assuming a dark
energy dominated Universe here, hence acceleration is achieved when
$w_{\mathrm{eff}}<-1/3$.

In general, several different dark energy models can be obtained from (4). For
instance, in some cases $w_{\mathrm{eff}}<-1$, and so the Shan–Chen dark
energy can have phantom behaviour, evolving between the two cosmological
constants, i.e. the two values $\rho_{\mathrm{S}}/\rho_{*}$ for which
$w_{\mathrm{eff}}=-1$. When $w>-1$ and $wg<0$ there are models that would
produce acceleration only below a certain energy density threshold, with the
energy density decreasing, eventually tending to a cosmological constant when
$w_{\mathrm{eff}}=-1$ is asymptotically approached. We leave a complete
analysis of the Shan–Chen dark energy scenario to future work, and now move to
a discussion of the Shan–Chen interacting vacuum model.

###  2.2 The interacting vacuum scenario

The classical vacuum is identified as a perfect fluid with equation of state
$P=-V$, where $V$ is the vacuum energy density. If non-interacting, $V$ is
equivalent to a cosmological constant. If interacting, however, $V$ can vary,
as we are now going to summarise. More details can be found in Wands et al.
(2012).

We consider the specific case of the vacuum interacting with cold dark matter
(CDM), as in previous works (Salvatelli et al., 2014; Martinelli et al., 2019;
Hogg et al., 2020). The vacuum energy–momentum tensor is

$\displaystyle\check{T}^{\mu}_{\nu}=-Vg^{\mu}_{\nu},$ (5)

where $V$ is the vacuum energy density and $g^{\mu}_{\nu}$ is the metric
tensor. By comparing (5) with the energy–momentum tensor of a perfect fluid,

$\displaystyle T^{\mu}_{\nu}=Pg^{\mu}_{\nu}+(\rho+P)u^{\mu}u_{\nu},$ (6)

where $P$ is the pressure, $\rho$ the energy density and $u^{\mu}$ the
4-velocity of the fluid, we can identify $V=-P=\rho$.

It is clear from (5) that any timelike 4-vector is an eigenvector of
$\check{T}^{\mu}_{\nu}$, with $V$ its eigenvalue: thus, $V$ is the energy
density of vacuum in any frame, i.e. for any observer.

CDM is represented by a pressureless perfect fluid, hence its energy–momentum
tensor is given by

$\displaystyle T^{\mu}_{\nu}=\rho_{c}u^{\mu}u_{\nu},$ (7)

where $\rho_{c}$ is the rest-frame energy density of CDM and $u^{\mu}$ its
4-velocity. We can introduce an energy exchange between the two components via
an energy–momentum flow 4-vector $Q^{\nu}$,

$\displaystyle\nabla_{\mu}T^{\mu}_{\nu}$ $\displaystyle=-Q_{\nu},$ (8)
$\displaystyle\nabla_{\mu}\check{T}^{\mu}_{\nu}$
$\displaystyle=-\nabla_{\nu}V=Q_{\nu},$ (9)

so that the total energy–momentum tensor is conserved, as it should be in the
framework of General Relativity as a consequence of the Bianchi identities.
This interaction 4-vector can be projected in two parts, one parallel and one
orthogonal to the cold dark matter 4-velocity $u^{\mu}$,

$\displaystyle Q^{\mu}=Qu^{\mu}+f^{\mu},$ (10)

where $Q=-Q_{\mu}u^{\mu}$ represents the energy exchange and $f^{\mu}$ the
momentum exchange between cold dark matter and the vacuum, in the frame
comoving with CDM, and $f_{\mu}u^{\mu}=0$.

As in previous works (Salvatelli et al., 2014; Martinelli et al., 2019; Hogg
et al., 2020), we now impose the geodesic condition on the interaction,
meaning that we neglect the momentum exchange by setting $f^{\mu}=0$. Since
$f^{\mu}=\rho_{c}a^{\mu}$, where $a^{\mu}$ is the 4-acceleration, this means
that there is no additional acceleration on the cold dark matter particles due
to the interaction, and they hence remain geodesic. With this choice, in the
synchronous comoving gauge used in CAMB, the interaction is unperturbed and
fully encoded in the background $Q$. However, the interaction does enter the
evolution equation for the CDM density contrast. A detailed discussion of the
linear perturbations in the interacting vacuum can be found in Wands et al.
(2012); Wang et al. (2013); Martinelli et al. (2019).

In an FLRW background, (8) and (9) reduce to

$\displaystyle\dot{\rho_{c}}+3H\rho_{c}$ $\displaystyle=-Q,$ (11)
$\displaystyle\dot{V}$ $\displaystyle=Q,$ (12)

and the Friedmann–Raychaudhuri equation for CDM interacting with vacuum is
formally the same as in $\Lambda$CDM:

$\displaystyle\frac{\ddot{a}}{a}=-\frac{4\pi G}{3}(\rho_{c}-2V),$ (13)

where in the cosmological constant case, $V=\Lambda/8\pi G$. Accordingly, an
interacting (and thus evolving) $V>0$ always contributes positively to
$\ddot{a}$, eventually accelerating the expansion if and when it becomes the
dominant component. This is unlike the case of a dark energy component with an
equation of state $P=P(\rho)$ such as the Shan–Chen dark energy model given in
(1) which, depending on its parameter values, may or may not positively
contribute to $\ddot{a}$.

###  2.3 The Shan–Chen interacting vacuum model

We can now recast the Shan–Chen model as a parameterisation of the coupling
$Q$ between the vacuum and cold dark matter, thus introducing the Shan–Chen
interacting vacuum model.

Starting from (3), we substitute $P$ from (1), replace $\rho$ with the vacuum
energy density $V$ and, to avoid confusion with the dark energy case, we
rename the parameter $w$ as $\beta$. We also introduce the dimensionless
parameter $q$ which controls the overall strength of the interaction. This
yields the final form of the coupling between the vacuum and cold dark matter,

$\displaystyle Q=\dot{V}=-3Hq\left[(1+\beta)V+\frac{\beta
g}{2}\rho_{*}\left(1-e^{-\alpha\frac{V}{\rho_{*}}}\right)^{2}\right].$ (14)

When $q=0$, there is no interaction, and we therefore return to $\Lambda$CDM.

We have a number of additional free parameters in the Shan–Chen interacting
vacuum model with respect to $\Lambda$CDM: $q$, $\alpha$, $g$, $\beta$ and
$\rho_{*}$. It is important to note that identical models can be obtained
through a remapping of some of these parameters. Specifically, for any fixed
pair $\alpha$ and $\rho_{*}$, and a triad $q$, $\beta$ and $g$, we can choose
any $\beta^{\prime}$ and then choose

$\displaystyle q^{\prime}$
$\displaystyle=\frac{q(1+\beta)}{(1+\beta^{\prime})},$ (15) $\displaystyle
g^{\prime}$ $\displaystyle=\frac{q\beta g}{q^{\prime}\beta^{\prime}},$ (16)

thus leaving (14) invariant. Therefore, in order to simplify our analysis, we
will fix the values of some of these parameters. In the following sections, we
will fix the characteristic energy scale to the value of the critical density
today in $\Lambda$CDM, i.e. $\rho_{*}=\bar{\rho}_{\mathrm{crit},0}$, as we are
interested in the effect that the interaction may have at late times. As
previously mentioned, another choice could be the value of the matter density
at matter–radiation equality. We will discuss the values chosen for the other
parameters below, after showing the effect that changing each one has on the
overall behaviour of this interacting vacuum.

Before proceeding to the main analysis, therefore, let us formally define an
effective equation of state for the Shan–Chen interacting vacuum, which we
call $w_{\mathrm{int}}$, so that we may better understand the behaviour of the
vacuum energy density $V$ with this interaction. In analogy with (3) and (4)
we can rewrite (14) as

$\displaystyle Q=\dot{V}=-3HV(1+w_{\mathrm{int}}),$ (17)

where we define222This specific definition is also motivated by its practical
use in CAMB.

$w_{\mathrm{int}}=q(1+\beta)+\frac{q\beta g}{2x}\left(1-e^{-\alpha
x}\right)^{2}-1,$ (18)

and we have now defined $x=V/\rho_{*}$. Clearly, again in analogy with the
dark energy case, there can be values $V=V_{\mathrm{S}}$ such that
$w_{\mathrm{int}}=-1$, where $\dot{V}=0$. If they exist, these stationary
points appear in pairs and they represent cosmological constants, one stable
and one unstable, separating models where $V$ decays from models where $V$
grows. Thus there is again a variety of possible models; for instance, models
where $V$ evolves between these two constants, as well as models where $V$
asymptotically decays or grows toward the stable cosmological constant
(stationary point).

In Figure 1, we show the effects of different parameter values on the
behaviour of the function $w_{\mathrm{int}}$, comparing it with the non-
interacting equation of state $w_{\mathrm{eff}}$ (4) of the original Shan–Chen
fluid dark energy models in the upper panels.

In the upper left panel we show four examples for $w_{\mathrm{eff}}$, given by
$w=\pm 1/3$ and $g=\pm 8$, with $\alpha$ fixed to the best fit value of
$\alpha=2.7$ found by Bini et al. (2013). We can see that when $wg>0$, i.e.
the light blue and yellow curves, the effective equation of state is always
outside the accelerating regime, i.e. always $w_{\rm eff}>-1/3$. However, when
$wg<0$, i.e. the red and dark blue curves, the effective equation of state can
allow for an accelerating phase. Note that in all the plots in this figure,
values of $x\gg 1$ correspond to high energies, i.e. energies which are large
with respect to the reference energy scale $\rho_{*}$.

For the case where $w$ is positive and $g$ is negative (red curve, $w=1/3$),
the Shan–Chen fluid behaves as radiation at high energies in the region
$w_{\mathrm{eff}}>-1$, decaying and transitioning through a dark energy
accelerating phase and toward a cosmological constant (the $x$ value
corresponding to the right-most point on the line $w_{\mathrm{eff}}=-1$
crossed by the red line).

On the other hand, in the region $w_{\mathrm{eff}}<-1$, the Shan–Chen fluid is
trapped in a phantom regime, growing between the two cosmological constants
corresponding to the $x$ values where the line $w_{\mathrm{eff}}=-1$ is
crossed by the red line. When $w\leq-1/3$ and $g$ is positive (dark blue
curve), the Shan–Chen fluid behaves as a dark energy even at high energies,
again evolving toward a cosmological constant; the phantom regime is
qualitatively the same as for the red curve.

In the upper right panel of Figure 1, we show the same cases (the same
parameter values corresponding to the same colours) for the Shan–Chen
interacting vacuum, i.e. the effective equation of state $w_{\mathrm{int}}$
given in (18), fixing $q=-0.1$. Note that the y-axis has a different scale to
the y-axis in the upper left panel and, because $q<0$, the behaviour of curves
with the same colour is opposite. However, as discussed after (13), in the
interacting vacuum case $V$ always gives a positive contribution to
acceleration; therefore for the Shan–Chen interacting vacuum this is
independent from the sign combination of $\beta$ and $g$ (recalling that we
renamed $w$ as $\beta$ in the interacting case). For $w_{\mathrm{int}}<-1$
there can be cases allowing for $V$ growing between two cosmological constants
(not shown), corresponding to the dark energy phantom cases in the top left
panel, as well as cases where $V$ grows indefinitely from zero (yellow and
light blue lines) or from a cosmological constant (blue and red lines). Above
$w_{\mathrm{int}}=-1$ the red and blue lines give models where $V$ decays
between the two cosmological constants.

In the rest of Figure 1, we show the effects on $w_{\rm int}$ of different
values of $q$ (middle left panel), different values of $\beta$ (middle right
panel), different values of $g$ (bottom left panel) and different values of
$\alpha$ (bottom right). If the parameter in question is not shown with
different values it is kept fixed to the best fit values of Bini et al.
(2013): $\beta=1/3$, $g=-8.0$ and $\alpha=2.7$, with $q$ fixed to $-0.1$ as an
example.

Figure 1: Upper panels: four examples of the behaviour of $w_{\rm eff}$ (upper
left) and $w_{\rm int}$ for $w=\pm 1/3$ (a.k.a. $\beta$) and $g=\pm 8$ as a
function of $x$, the energy density normalised to the energy scale $\rho_{*}$.
In the top-left panel the dashed black line marks the transition to the
accelerating regime for $w_{\rm eff}<-1/3$ and the solid black line separates
phantom models, with $w_{\rm eff}<-1$, from standard dark energy ones.
Similarly, in all other panels the black solid line separates models where the
vacuum energy density $V$ decays, with $w_{\rm int}>-1$, from models where $V$
grows. In the top-right panel we have chosen a negative $q=-0.1$, so that
curves here have the opposite behaviour of those with the same colour (same
parameter values) in the top-left panel. Middle and lower panels: the
behaviour of $w_{\rm int}$ as a function of $x=V/\rho_{*}$ for different
values of $q$ (middle left), $\beta$ (middle right), $g$ (lower left) and
$\alpha$ (lower right). When not changed, the parameters are kept fixed to
$q=-0.1$, $\beta=1/3$, $g=-8$ and $\alpha=2.7$.

Starting with the middle left panel of Figure 1, we can see the effect of
changing both the magnitude and sign of the interaction strength $q$; in
particular, changing the sign of $q$ simply mirrors curves with same parameter
values about the $w_{\rm int}=-1$ line. In all cases shown there are two
cosmological constants where the curves cross the $w_{\rm int}=-1$ line. In
general, at high energies (i.e. $x\gg 1$) $w_{\rm int}$ behaves as a constant,
then the nonlinear term in (18) becomes dominant around $x=1$. For models
above the $w_{\rm int}=-1$ line (red and yellow lines) $V$ decays from high
energies to the cosmological constant on the right, or to zero from the
cosmological constant on the left. For models below the $w_{\rm int}=-1$ line
(light blue and dark blue lines) $V$ grows from the cosmological constant on
the right to high energies, or from zero to the cosmological constant on the
left. Models in between the two cosmological constants decay or grow between
the two for $w_{\rm int}>-1$ and $w_{\rm int}<-1$, respectively.

Next, in the middle right panel of Figure 1, we see the effect that changing
$\beta$ has. Generally, increasing $\beta$ moves the curves up ($w_{\rm
int}=q(1+\beta)-1$ at $x=0$) and increases the relevance of the nonlinear term
in (18) at low energies while shifting the point at which this term becomes
dominant over the constant one to higher energies. In addition, changing the
sign of $\beta$ flips the curves upside down around the line $w_{\rm
int}=q-1$, as is also clear from the red and the yellow curves in the upper
right panel.

Third, in the lower left panel of Figure 1, we can see the effect of changing
$g$. The result is similar to that of $\beta$, with no vertical translation.
In addition, changing the sign of $g$ flips the curves upside down around the
line $w_{\rm int}=q(1+\beta)-1$.

Finally, in the lower right panel of Figure 1, we can see the effect of
changing $\alpha$. This parameter is called the saturation scale by Bini et
al. (2013), as in the fluid Shan–Chen dark energy model it represents the
typical density above which $\psi$ undergoes a saturation effect, i.e.
$\psi\approx 1$. For this reason, $\alpha$ should always be positive. We can
see that changing $\alpha$ greatly affects the relevance and steepness of the
exponential term in $w_{\rm int}$, while at the same time shifting the $x$
position of the extrema points.

In Figure 2 we summarise all the qualitatively different models that can arise
in the Shan–Chen interacting vacuum scenario, depending on the behaviour of
$w_{\rm int}(x)$. First of all, we only illustrate cases where $w_{\rm int}<0$
always; otherwise vacuum may never become dominant over CDM, or not dominant
enough, or for the right period of time. Secondly, vacuum decreases or
increases for models above or below the $w_{\rm int}=-1$ line respectively,
representing $\Lambda$CDM. Lastly, the asymptotic value of $w_{\rm int}$ for
large $x$ is the same as $w_{\rm int}(0)$, i.e. $w_{\rm int}=q(1+\beta)-1$,
and $w_{\rm int}$ has a single minimum or maximum in between. Thus the curve
$w_{\rm int}(x)$ can either be entirely above or below $w_{\rm int}=-1$, or
$w_{\rm int}=-1$ can be crossed at two points. Hence, there are six possible
different curves. For each of the two cases crossing $w_{\rm int}=-1$ there
are three models. Otherwise, each curve represents a model; there are thus ten
different types of model in total.

The two cases where $w_{\rm int}>-1$ for any $x$, i.e. the red and cyan lines,
represent two models where vacuum decays to zero in the far future, coming
from a power law decay $V\sim a^{-3q(1+\beta)}$ ($q>0,\beta>-1$) in the far
past, at high energies, for $x\gg 1$. Conversely, the two cases in which
$w_{\rm int}<-1$ for any $x$, i.e. the dark blue and yellow lines, represent
two models where vacuum grows from zero in the past, eventually growing as the
power law $V\sim a^{-3q(1+\beta)}$ ($q<0,\beta>-1$) in the far future, for
$x\gg 1$.

Finally, each line that is the same type as the pale blue and orange lines
represents three models. This is because the two points where these lines
cross $w_{\rm int}=-1$ represent stationary points of (17), i.e. cosmological
constants, that cannot be crossed by the $V$ evolution; they are asymptotic
values, either in the past or in the future. Thus, the portion of the orange
line at large $x$ represents a model where $V$ decays, initially as $V\sim
a^{-3q(1+\beta)}$ and eventually tending to a constant. The portion of the
pale blue line at large $x$ represents a model where in the past $V$ is
asymptotic to a constant, then grows in the future as $V\sim
a^{-3q(1+\beta)}$. The portion of the orange line below $w_{\rm int}=-1$
represents a model where $V$ grows between the two cosmological constants,
while the portion of the pale blue line above $w_{\rm int}=-1$ represents a
model where $V$ decays, from a constant value in the past to a smaller
constant value in the future.

The portion of the orange and pale blue lines between $x=0$ and the point
where $w_{\rm int}=-1$ is crossed represent models where $V$ either decays
from zero to a cosmological constant at $x<1$ or vice versa. These type of
models are not relevant here, where we assume $\rho_{*}=\rho_{{\rm crit},0}$,
but could be relevant when studying the case of a higher energy $\rho_{*}$.
For example, models represented by this portion of the orange line, with $V$
decaying from a past cosmological constant as a power law $V\sim
a^{-3q(1+\beta)}$, could be relevant when using the Shan–Chen interacting
vacuum to mimic the type of early dark energy model presented in Poulin et al.
(2019), where the dark energy follows precisely this type of evolution.
Finally, models represented by this portion of the pale blue line could be
used to study the case of an emerging vacuum, which initially grows from zero
as $V\sim a^{-3q(1+\beta)}$ and asymptotically tends to a cosmological
constant.

Figure 2: The six possible cases of $w_{\rm int}$, for a total of ten possible
evolution models for the vacuum energy density $V$, decreasing for $w_{\rm
int}>-1$ and increasing for $w_{\rm int}<-1$. The cyan, red, blue and yellow
lines each represent a model where $V$ varies between zero and infinity,
either growing or decaying in time. Each of the pale blue and orange lines
represent three different models, separated by the two points at the crossing
of the $w_{\rm int}=-1$ line. These points represent cosmological constants,
one stable and one unstable, i.e. asymptotic states in the evolution of each
of the models represented by these lines.

Having demonstrated the effect of all the model parameters on $w_{\rm int}$,
we now choose to reduce our parameter space somewhat by fixing two of them.
For the remainder of our analysis, we fix $g=-8.0$ and $\alpha=2.7$, both of
which are the best fit values found by Bini et al. (2013). We fix these two
parameters in particular, as changing their values is rather similar to the
effect of $\beta$, and would result in a degeneracy with $\beta$ when
sampling.

To further demonstrate that the general effect is still that of an interacting
vacuum (as seen in previous works, e.g. Martinelli et al. (2019)), in Figures
3 and 4 we plot the CMB temperature–temperature power spectrum and matter
power spectrum at $z=0$ for different values of $q$ with $\beta$ fixed, and
for different values of $\beta$ with $q$ fixed. For the purposes of these
plots, we keep the cosmological parameters fixed to the Planck 2018 best fits
(Aghanim et al., 2020). As expected, the presence of a coupling between the
vacuum and cold dark matter acts to boost or suppress the peaks of the CMB
power spectrum, with an opposite effect on the matter power spectrum. These
plots also serve as an order of magnitude guide for the priors we will set on
these parameters in our MCMC analysis which we describe in the next section.

Figure 3: The CMB temperature–temperature power spectra for different values
of $q$ (upper panel) and different values of $\beta$ (lower panel). For the
upper panel, $\beta$ is kept fixed to $-1$ and for the lower panel, $q$ is
kept fixed to $-0.1$, hence the light blue curves in both plots represent the
same case.

Figure 4: The matter power spectrum at $z=0$ for different values of $q$
(upper panel) and for different values of $\beta$ (lower panel). For the upper
panel, $\beta$ is kept fixed to $-1$ and for the lower panel, $q$ is kept
fixed to $-0.1$, hence the light blue curves in both sets of plots represent
the same case.

##  3 Method and data

In order to constrain the free parameters of the Shan–Chen interacting vacuum
model ($q$ and $\beta$) along with the cosmological parameters, we modify the
publicly available Boltzmann code CAMB and its associated MCMC sampler
CosmoMC. We implement equation (14) in CAMB, fixing $\rho_{*}$ to the value of
$\bar{\rho}_{{\rm crit},0}$, i.e. the critical energy density today in
$\Lambda$CDM, $\bar{\rho}_{\rm crit,0}=3H_{0}^{2}/8\pi G=8.5\times 10^{-27}$kg
m-3 when $H_{0}=67.4$ kms-1Mpc-1 (Aghanim et al., 2020).

Using CosmoMC, we sample the posterior distributions of the baryon and cold
dark matter densities $\Omega_{b}h^{2}$ and $\Omega_{c}h^{2}$, the amplitude
of the primordial power spectrum and the spectral index $A_{s}$ and $n_{s}$,
the value of the Hubble parameter today, $H_{0}$, as well as the Shan–Chen
parameters $q$ and $\beta$. We keep $\alpha$ fixed to 2.7 and $g$ fixed to
$-8.0$. We list the flat priors we impose on the sampled parameters in Table
1.

We use the Planck 2018 measurements of the CMB temperature and polarisation
(Aghanim et al., 2020) together with the BAO measurements from the 6dF Galaxy
Survey (Beutler et al., 2011), the SDSS Main Galaxy Sample (Ross et al., 2015)
and the SDSS DR12 consensus catalogue (Alam et al., 2017), and the Pantheon
catalogue of Type Ia supernovae (Scolnic et al., 2018).

Note that, in contrast to some previous works, we do not include the
S$H_{0}$ES measurement of $H_{0}$ in our chosen set of data. This is because
it is nonsensical to combine datasets which are in tension in a statistical
analysis of this kind. This statement is hardly original (see e.g. Hill et al.
(2020); Ivanov et al. (2020)), yet bears repeating due to the proliferation of
analyses which claim to resolve the $H_{0}$ tension while using S$H_{0}$ES
(see some of the results listed in Table B2 of Di Valentino et al. (2021a)).
If the S$H_{0}$ES result is to be used when constraining an alternative
cosmological model, it must first be checked that the value of $H_{0}$ in that
model that has been found by a combination of cosmological-model-dependent
datasets, such as CMB+BAO+SNIa, is consistent with the S$H_{0}$ES value
itself. Only then is it safe to run an analysis using e.g.
CMB+BAO+SNIa+S$H_{0}$ES. Furthermore, as recently noted, it may be more
appropriate to consider a prior on the supernova absolute magnitude rather
than $H_{0}$ when using the S$H_{0}$ES data (Benevento et al., 2020; Camarena
& Marra, 2021; Efstathiou, 2021)

Parameter | Prior
---|---
$\Omega_{b}h^{2}$ | $[0.005,0.1]$
$\Omega_{c}h^{2}$ | $[0.001,0.99]$
$\log{10^{10}A_{s}}$ | $[2.0,4.0]$
$n_{s}$ | $[0.8,1.2]$
$H_{0}$ | $[50,100]$
$q$ | $[-0.25,0.25]$
$\beta$ | $[-1.0,1.0]$
Table 1: Prior ranges of the parameters sampled in our analysis.

##  4 Results and discussion

We now proceed to the results of our MCMC analysis, beginning with case I, in
which we keep $\beta$ fixed and sample $q$, followed by case II, in which we
keep $q$ fixed and sample $\beta$, and finally case III, where we sample both.
Note that the following is not an exhaustive exploration of every parameter
combination, but serves to represent the main models under the Shan–Chen
interacting vacuum umbrella. For the sake of completeness, we include the
study of two additional sub-cases in Appendix A.

###  4.1 Case I: different values of $q$

In Figure 5, we show the constraints obtained using the full combination of
data: CMB plus BAO plus supernovae. We study a case in which we fix $q$ to
zero (i.e. the $\Lambda$CDM limit), which we call case Ia; a case in which we
fix $q=0.1$, which we call case Ib; a case in which we fix $q=-0.1$, which we
call case Ic; and a case in which we sample over $q$ along with the
cosmological parameters, which we call case Id. In all of the sub-cases shown
here, $\beta$ is fixed to $1/3$.

We see that for all these cases, even when we sample over $q$ (yellow
contours) and therefore have an additional degree of freedom, the resulting
cosmologies are constrained by the data to be rather close to $\Lambda$CDM. We
quantify this in Table 2, where we report the mean posterior values and
$1\sigma$ limits obtained using GetDist (Lewis, 2019) for this and all
subsequent cases studied. The value of $q$ found in this case (Id) is
$q=-0.007\pm 0.1$, which is completely consistent with the $\Lambda$CDM limit
of $q=0$.

Since all the present-day values of the cosmological parameters in these
models are constrained to be very close to $\Lambda$CDM, we also see no
resolution of the tensions present in $\Lambda$CDM. To resolve the $H_{0}$
tension with this combination of data we would need to see $H_{0}$ reaching
higher values, e.g. $71$ to $73$ kms-1 Mpc-1. Besides the $H_{0}$ tension
there also exists a lesser tension between the value of $\sigma_{8}$ obtained
from Planck (Aghanim et al., 2020) and from large scale structure surveys such
as the Dark Energy Survey (DES) (Abbott et al., 2020) and the Kilo Degree
Survey (KiDS) (Heymans et al., 2021). A resolution of this tension would
require a smaller value of $\sigma_{8}$ with this combination of data – again,
we do not find this.

###  4.2 Case II: different values of $\beta$

Next, in Figure 6, we show the result of fixing $q=-0.1$ and fixing $\beta$ to
two different values, $1/3$ and $-1/3$, and then sampling it in the range
$[-1.0,1.0]$. In this plot, we see that the posterior distribution of $\beta$
in the case which it is sampled over (purple contours, case IIb) describes the
allowed range of the posterior distributions that result from $\beta$ being
fixed within the range of the prior used when it is sampled,
$\beta\in[-0.1,0.1]$. The other two cases shown are $\beta=1/3$ (blue
contours; note that this case is identical to case Ic) and $\beta=-1/3$
(orange contours, case IIa). The mean posterior values for all these cases are
again shown in Table 2.

Interestingly, we see that in case IIa, a larger value of $H_{0}$ is found,
demonstrating how interacting vacuum models can in general be invoked to
resolve the tension in this parameter. However, the value of $H_{0}$ in this
case does not provide a full resolution of the tension, but merely a slight
relaxation, since its marginalised posterior value is $69.2\pm 0.5$ kms-1
Mpc-1. We can see from the very large value of $\sigma_{8}$ that the tension
in this parameter is significantly worsened in this model.

###  4.3 Case III: sampling $q$ and $\beta$

We finally consider the situation in which we sample over both $q$ and
$\beta$. The constraints on the cosmological and model parameters for this
case are shown in Figure 7. We can see in this case that $q$ is relatively
unconstrained but that $\beta$ is constrained to be $0.39\pm 0.2$. However, we
do not necessarily expect $\beta$ to be constrained to be close to zero, as
zero is not a $\Lambda$CDM limit of the model; in fact, $\beta=0$ corresponds
to the Cfix case studied in Martinelli et al. (2019).

A reason for $q$ to be less constrained than $\beta$ can be found in the best
fit value for $\beta$ that we find. Indeed, for $\beta\approx 0.35$ (and for
the given fixed values of $g$ and $\alpha$ that we use) the today’s value of
the function $w_{\mathrm{i}nt}$ in (18) is very close to $-1$, resulting is a
slow evolution of $V$ that is then relatively unaffected by the value of $q$
in the given prior.

Figure 5: Case I: constraints on $\Omega_{\mathrm{m}}$, $H_{0}$ and $\sigma_{8}$ for $\Lambda$CDM (black), and the Shan–Chen model with $q$ fixed to $0.1$ (red), $-0.1$ (blue) and sampled over the range $[-0.25,0.25]$ (yellow), using the Planck 2018 plus BAO plus Pantheon supernova catalogue. In all of these cases, $\beta=1/3$. Figure 6: Case II: constraints on $\Omega_{\mathrm{m}}$, $H_{0}$ and $\sigma_{8}$ for the Shan–Chen model with $q=-0.1$ and $\beta=1/3$ (blue; note that this case is identical to the blue contour shown in Figure 5), $\beta=-1/3$ (orange) and $\beta$ sampled over the range $[-1.0,1.0]$ (purple).  Figure 7: Case III: constraints on $\Omega_{\mathrm{m}}$, $H_{0}$, $q$, $\beta$ and $\sigma_{8}$ for the Shan–Chen model with both $q$ and $\beta$ sampled over the ranges $[-0.25,0.25]$ and $[-1.0,1.0]$ respectively. Parameter | Case | Value | Parameter | Case | Value
---|---|---|---|---|---
| Ia | $0.02244\pm 0.00013$ |  | Ia | $67.8\pm 0.4$
| Ib | $0.02243\pm 0.00013$ |  | Ib | $67.8\pm 0.4$
| Ic | $0.02244\pm 0.00013$ |  | Ic | $67.7\pm 0.4$
$\Omega_{b}h^{2}$ | Id | $0.02244\pm 0.00013$ | $H_{0}$ | Id | $67.8\pm 0.4$
| IIa | $0.02235\pm 0.00013$ |  | IIa | $69.2\pm 0.5$
| IIb | $0.02241\pm 0.00014$ |  | IIb | $68.2\pm 0.8$
| III | $0.02243\pm 0.00013$ |  | III | $67.9\pm 0.5$
| Ia | $0.119\pm 0.0010$ |  | Ia | $0.8110\pm 0.0044$
| Ib | $0.119\pm 0.0010$ |  | Ib | $0.8132\pm 0.0048$
| Ic | $0.120\pm 0.0010$ |  | Ic | $0.8086\pm 0.0039$
$\Omega_{c}h^{2}$ | Id | $0.119\pm 0.0011$ | $\sigma_{8}$ | Id | $0.8107\pm 0.0055$
| IIa | $0.098\pm 0.0013$ |  | IIa | $0.9619\pm 0.0030$
| IIb | $0.123\pm 0.0102$ |  | IIb | $0.8594\pm 0.071$
| III | $0.117\pm 0.0041$ |  | III | $0.8245\pm 0.028$
| Ia | $3.050\pm 0.005$ |  | Ia | $0.0$
| Ib | $3.050\pm 0.006$ |  | Ib | $0.1$
| Ic | $3.050\pm 0.006$ |  | Ic | $-0.1$
$\log{10^{10}A_{s}}$ | Id | $3.050\pm 0.006$ | $q$ | Id | $-0.008\pm 0.1$
| IIa | $3.053\pm 0.006$ |  | IIa | $-0.1$
| IIb | $3.051\pm 0.006$ |  | IIb | $-0.1$
| III | $3.050\pm 0.006$ |  | III | $-0.028\pm 0.1$
| Ia | $0.9667\pm 0.0035$ |  | Ia | $1/3$
| Ib | $0.9674\pm 0.0036$ |  | Ib | $1/3$
| Ic | $0.9675\pm 0.0036$ |  | Ic | $1/3$
$n_{s}$ | Id | $0.9675\pm 0.0036$ | $\beta$ | Id | $1/3$
| IIa | $0.9645\pm 0.0036$ |  | IIa | $-1/3$
| IIb | $0.9665\pm 0.0040$ |  | IIb | $0.12\pm 0.32$
| III | $0.9673\pm 0.0039$ |  | III | $0.31\pm 0.16$
Table 2: Marginalised parameter values with $1\sigma$ errors for each case.
Note that all parameters quoted are dimensionless, besides $H_{0}$ which has
dimensions of kms-1Mpc-1.

###  4.4 Statistical and physical comparison of models

In Table 3, we show the $\chi^{2}$ and $\Delta\chi^{2}=\chi^{2}_{\rm
SC}-\chi^{2}_{\Lambda{\rm CDM}}$ for each case studied. Defined in this way, a
negative $\Delta\chi^{2}$ thus implies that the model under consideration is a
better fit to the data than $\Lambda$CDM, apart from in cases Id, IIb and III,
where we must compute the significance of the change in $\chi^{2}$ using a
$\chi^{2}$ difference table (e.g. Dodge (2008)). This is because these cases
have additional degrees of freedom with respect to $\Lambda$CDM. Case Id and
IIb each have one additional degree of freedom, because $q$ is allowed to vary
in the former and $\beta$ in the latter. Case III have two additional degrees
of freedom because both $q$ and $\beta$ are allowed to vary.

From Table 3, we can see that out of the cases with the same number of degrees
of freedom as $\Lambda$CDM, we have no cases with a negative $\Delta\chi^{2}$.
This implies that none of these models are a better fit to the data than
$\Lambda$CDM. For the comparison of $\Lambda$CDM with cases with additional
degrees of freedom, we choose a significance level of 95%. This means that for
one additional degree of freedom, a $|\Delta\chi^{2}|>3.841$ would be
considered significant. The $\Delta\chi^{2}$ values for case Id and IIb are
both smaller than this value, so $\Lambda$CDM is a better fit than these
models. For two additional degrees of freedom, a $|\Delta\chi^{2}|>5.991$ is
significant. The $\Delta\chi^{2}$ for case III is $0.82$, so $\Lambda$CDM is
again a better description of the data than this model.

Furthermore, it is possible to produce some additional models which have the
same number of degrees of freedom as $\Lambda$CDM, from the best fit values of
the parameters found in cases Id, IIb and III. Specifically, we run three
additional MCMC chains with $q$ and $\beta$ fixed to their best fit values
found (see Table 2) in those cases. We allow the usual cosmological parameters
to vary. The result of these chains is therefore three more Shan–Chen models
with the same number of free parameters as $\Lambda$CDM, which we plot in
Figure 8. We can therefore compute the $\Delta\chi^{2}$ for these additional
best fit cases with no need to penalise for the additional degrees of freedom.

Figure 8: The results of sampling the cosmological parameters while fixing $q$
and $\beta$ to their best fit values found in case Id (yellow), IIb (purple)
and case III (green, open). The $\Lambda$CDM result is shown with dashed black
lines.

The $\Delta\chi^{2}$ in these best fit cases are: $\Delta\chi^{2}=0.59$ for
the best fit case Id; $\Delta\chi^{2}=0.27$ for the best fit case IIb; and
$\Delta\chi^{2}=0.61$ for the best fit case III. The best fit Id and III are
both very close to $\Lambda$CDM and indeed their original $\Delta\chi^{2}$
values, representing no improvement in fit. However, the best fit case IIb
sees a relatively good improvement in $\chi^{2}$, and the resulting
cosmological parameters are fairly distinct from those of $\Lambda$CDM –
especially $\Omega_{m}$ and $\sigma_{8}$ – as can be seen in Figure 8. This is
because in the best fit case IIb, $q$ is fixed to $-0.1$, whereas in the other
two cases $q$ is close to the $\Lambda$CDM limit of zero (recalling that these
values of $q$ came from allowing $q$ to vary), thus producing models with
cosmological parameters which are more similar to $\Lambda$CDM.

When simply looking at the $\Delta\chi^{2}$ values for the cases we have
studied, a casual reader may conclude that there is very little to choose
between the various models. However, this comparison hides a deeper truth
about these cases: they describe completely different evolution histories (and
futures) for the Universe. To demonstrate this, we can once again plot the
behaviour of $w_{\rm int}$, as we did in section 2, but this time for the best
fit models we have found through our MCMC analysis.

Case | $\chi^{2}$ | $\Delta\chi^{2}$
---|---|---
Ia ($\Lambda$CDM) | $3831.38$ | $0.0$
Ib | $3831.87$ | $0.49$
Ic | $3831.90$ | $0.52$
Id | $3831.98$ | $0.60$
IIa | $3833.95$ | $2.57$
IIb | $3832.76$ | $1.38$
III | $3832.20$ | $0.82$
Best fit Id | $3831.97$ | $0.59$
Best fit IIb | $3831.65$ | $0.27$
Best fit III | $3831.99$ | $0.61$
Table 3: The $\chi^{2}$ and $\Delta\chi^{2}$ for each case studied.

Figure 9: Top panel: the behaviour of $w_{\rm int}$ for the best fit models.
Middle panel: $V(z)/\rho_{{\rm crit},0}$ for the best fit models. Bottom
panel: the evolution of $\Omega_{V}(z)=8\pi GV(z)/3H^{2}$ (solid lines) and
$\Omega_{m}(z)=8\pi G\rho_{m}(z)/3H^{2}$ (dashed lines) for the best fit
models. Note that at $z=0$, the values of $V(z)/\rho_{{\rm crit},0}$ shown in
the middle panel are equivalent to the values of $\Omega_{V}(z)$ in the bottom
panel.

Figure 10: Left panel: $w_{\rm int}$ for case Ib, Ic, Id and III plotted with
solid lines, with current values $V_{0}=V(0)$ for the same cases,
$x_{0}=V_{0}/\bar{\rho}_{\rm crit,0}$, represented by vertical dashed lines.
For cases Ib, Ic and Id, the current value is between the two cosmological
constants where $w_{\rm int}=-1$ (the rightmost of the two is shown in the
plot); thus in these models the vacuum $V$ evolves between these two
constants. For case III the current value is to the right of $w_{\rm int}=-1$
(where the green and the black line cross), hence $V$ in this model is growing
away from this constant. Right panel: $V(z)/\rho_{{\rm crit},0}$ for case Ib,
Ic, Id and III. As expected from their current values shown in the left panel,
$V$ is growing in cases Ib and III and decaying for cases Ic and Id.

In the top panel of Figure 9 (to be compared with Figure 2), we show the
behaviour of the function $w_{\rm int}(x)$ for all the cases studied in this
work, with $q$ and $\beta$ fixed to their best fit values (shown in Table 2)
for those cases in which these parameters were sampled. From this plot, we can
see that cases Ic, Id, IIb and III have a maximum in $w_{\rm int}$, while
cases Ib and IIa have a minimum. Furthermore, cases IIa and IIb represent
models in which $w_{\rm int}<-1$ always, i.e. models where the vacuum energy
density $V$ is indefinitely growing from zero in the past. For the cases in
which $w_{\rm int}$ crosses the $-1$ line, i.e. Ib, Ic, Id and III, we can
examine the evolution of $V(z)$ in order to determine at which point in
$w_{\rm int}$ the particular model sits, i.e. whether $V$ is growing or
decaying.

In the middle panel of Figure 9, we show the evolution of the vacuum energy as
a function of redshift, $V(z)$, normalised by the critical density today for
each of these models, $\rho_{{\rm crit},0}$, so that the value of
$V(z)/\rho_{{\rm crit},0}$ at $z=0$ is equivalent to $\Omega_{V}$. We use CAMB
to calculate this evolution, using not only the best fit values of $q$ and
$\beta$ in each case, but also the best fit values of the cosmological
parameters found in our MCMC analysis. This means that $\rho_{{\rm crit},0}$
is also different for each model.

From this plot, we can see that cases IIa and IIb both give a very strong
growth in $V(z)$.On the other hand, for all of the above mentioned cases where
$w_{\rm int}$ crosses the $-1$ line, $V(z)$ is very slowly evolving, i.e.
$\rm{d}V/\rm{d}z\approx 0$. Comparing with the bottom panel of Figure 9, we
conclude that for those redshifts where $V(z)$ becomes important, eventually
dominating over matter, it is effectively close to a cosmological constant for
these four models.

Finally, in the bottom panel of Figure 9, we show the evolution of the matter
and vacuum densities in all the cases studied, with $\Omega_{V}(z)$ for each
case shown with a solid line and $\Omega_{m}(z)$ being shown with dashed
lines. For example, from this plot it is clear that case IIa is the most
different from $\Lambda$CDM (again shown in black), which is reflected in the
uncompetitive $\chi^{2}$ value for this case.

To enable us to clearly distinguish the four cases which cross $w_{\rm
int}=-1$, i.e. cases Ib, Ic, Id and III, we plot them separately in Figure 10,
where we zoom-in $w_{\rm int}$ in the left panel and we zoom-out $V(z)$ in the
right panel, i.e. we also increase the redshift range of this plot to
emphasise the evolution of $V(z)$ in each case.

In the left panel the vertical lines represent, for each model, the current
value $V_{0}=V(0)$ of the vacuum $V(z)$, normalised to $\bar{\rho}_{\rm
crit,0}$ (the reference critical energy density in $\Lambda$CDM (Aghanim et
al., 2020)), $x_{0}=V_{0}/\bar{\rho}_{\rm crit,0}$. These are the initial
values, from $z=0$, in the evolution of $V$, thus their position tells us what
the evolution is for each of the four models. We see that for cases Ib, Ic and
Id these initial values are between the two cosmological constants for these
models, where $w_{\rm int}=-1$ in the top panel of Figure 9, with the
rightmost one on the right of the vertical lines in Figure 10. Since for Ic
and Id the initial values are where $w_{\rm int}>-1$, for these models $V$ is
currently decreasing from the higher to the lower cosmological constant. For
case Ib, the initial value is where $w_{\rm int}<-1$, hence $V$ is currently
increasing from the lower to the higher cosmological constant. Finally, for
case III, the initial value is to the right of the rightmost cosmological
constant in the top panel of Figure 9, which here is the point where the green
and the black lines cross; therefore, $V$ is indefinitely growing from this
past asymptotic value.

In the right panel of in Figure 10, we plot $V(z)$ for the four slowly
evolving models. The evolution confirms what is is expected from the left
panel: case Ib is growing while case Ic is decaying; case Id is also decaying,
but at a much slower rate; case III is growing, but again at a slower rate
than Ib.

Let us now finally comment on models in which $V$ is always growing, such as
cases IIa and IIb: because $w_{\rm int}(x)<-1$ for any $x$, these should be
taken with a pinch of salt, and need to be studied more in detail in order to
establish if they are physically viable. Since for these models the coupling
(14) is positive and not proportional to the energy density of matter
$\rho_{c}$, it follows from (11) that $\rho_{c}$ can become negative; in
addition, it may well be the case that these models develop some sort of
future singularity, see Ananda & Bruni (2006a, b); Bamba et al. (2012), and
references therein. In this sense, perhaps they cannot be considered good
physical models for the entire history of the Universe, but rather a
parameterisation of the growth of the vacuum $V$ (akin to a a phantom-like
evolution).

##  5 Conclusions

The primary aim of this paper was to present a novel alternative to the
cosmological constant that could be responsible for the accelerating expansion
of the Universe, within the general framework of the interacting vacuum
scenario. In particular, we were motivated by previous models such as the
Chaplygin gas (Kamenshchik et al., 2001) or the Van der Waals equation of
state (Capozziello et al., 2005) to introduce a dark energy which has a
physical foundation, and is not purely an ad hoc phenomenological description
such as the Chevallier–Polarski–Linder (CPL) parameterisation of the dark
energy equation of state, $w=w_{0}+(1-a)w_{a}$ (Chevallier & Polarski, 2001;
Linder, 2003), in essence a Taylor expansion at low redshift (which is
actually less able to distinguish $w_{a}\neq 0$ than other so-called $w_{0}$,
$w_{a}$ parameterisations (Colgáin et al., 2021)).

Every dark energy or unified dark matter energy–momentum tensor can be
remapped into vacuum interacting with CDM (Wands et al., 2012). However,
especially in the case of unified dark matter, these models are often subject
to severe observational constraints because of the speed of sound of the dark
component (Sandvik et al., 2004), see also Balbi et al. (2007); Pietrobon et
al. (2008); Piattella et al. (2010); Piattella (2010); Gao et al. (2010);
Bertacca et al. (2011); De Felice et al. (2012); Wang et al. (2013); Li et al.
(2018); Li et al. (2019); an advantage of the interacting vacuum formulation
is that vacuum does not cluster and, in the specific sub-case of a pure energy
exchange we consider here (no momentum exchange in the rest frame of CDM) CDM
remains geodesic, as in $\Lambda$CDM. The generalised Chaplygin gas dark
energy model was reconsidered in this fashion by Bento et al. (2004) and Wands
et al. (2012), and Wang et al. (2013) considered observational constraints on
the same.

Following this approach, we cast the fluid dynamical equation of state
originally introduced by Shan & Chen (1993) in the context of lattice kinetic
theory into a CDM–vacuum interaction. This equation of state was considered as
a dark energy model in the restricted case of a purely FLRW cosmology by Bini
et al. (2013, 2016). As in previous works (Wands et al., 2012; Salvatelli et
al., 2014; Martinelli et al., 2019; Hogg et al., 2020), we started from a
general covariant formulation of the interacting vacuum scenario, including
perturbations, thereby introducing the Shan–Chen interacting vacuum class of
models.

The Shan–Chen equation of state is nonlinear and characterised by a reference
energy scale $\rho_{*}$, and so is the corresponding interaction, even if the
coupling becomes linear (Quercellini et al., 2008) at energies well above and
below $\rho_{*}$. Thanks to this nonlinearity a large variety of models can
arise, giving different dynamics, cf. Ananda & Bruni (2006a, b): in
particular, the vacuum energy density can either grow (from zero or from a
past constant value) or decay indefinitely (to zero or to a constant value),
or evolve between two cosmological constants. We studied a number of models
belonging to this class, and placed observational constraints on them.

Our findings show that the observational data we used (CMB temperature and
polarisation, BAO and SNIa measurements) are compatible with a wide range of
models which result in very different cosmologies. Models where $V$ grows
indefinitely require further theoretical analysis in order to establish if
they lead to a future singularity (Ananda & Bruni, 2006a, b; Bamba et al.,
2012) but, given that in these models the energy density of CDM becomes
negative in the future, more than anything else they can be considered a
parameterisation of the vacuum evolution.

The cases studied in this paper are only a small subset of all possible models
that exist under the umbrella of Shan–Chen dark energy. We have not given any
consideration to the original dark energy fluid as introduced in Bini et al.
(2013), focusing instead on the Shan–Chen interacting vacuum scenario which we
introduced in this work. Even within the Shan–Chen interacting vacuum model,
we have limited ourselves to fixing $\rho_{*}=\bar{\rho}_{\mathrm{crit},0}$,
using the critical energy density of the best fit Planck 2018 $\Lambda$CDM
(Aghanim et al., 2020)) as a reference energy scale. Other possibilities
exist, most obviously $\rho_{*}=\rho_{\mathrm{m}}(z_{\mathrm{eq}})$. In
addition, the Shan–Chen coupling considered here could be generalised in order
to prevent the vacuum $V$ from growing indefinitely and the CDM density
$\rho_{c}$ becoming negative. All of these remain open to exploration in
future works.

As we have previously mentioned, besides the desire to find a satisfactory
alternative to $\Lambda$, studies of alternative dark energy models such as
the Shan–Chen interacting vacuum model presented here can be motivated by
their potential to relax the $H_{0}$ tension. As we found in subsection 4.4,
cases Ib and Ic have the most competitive $\Delta\chi^{2}$, apart from the
best fit case IIb. We can see from Figure 5 that the values for $H_{0}$
obtained in Ib and Ic are virtually identical to those in $\Lambda$CDM. As
discussed earlier, there also exists a tension between the Planck 2018 value
of $\sigma_{8}$ and that obtained from large scale structure surveys. A
resolution of this tension would require a smaller value of $\sigma_{8}$. The
similarity of the cosmological parameters in these particular Shan–Chen
interacting vacuum models to those in $\Lambda$CDM means that there is also no
relaxation of the $\sigma_{8}$ tension.

The repeated failures of the interacting vacuum models to cure the $H_{0}$
tension point to the fact that perhaps this is the wrong line to continue down
to try and achieve this particular goal. While it is beyond the scope of this
work to explore these ideas, another type of dark energy such as the early
dark energy we discussed in section 2 could be a more fruitful avenue to
explore. Another possibility for future work could be the study of a model
which combines both early and late time effects – if the sole motivation is a
resolution of the $H_{0}$ tension.

Our repurposing of the Shan–Chen equation of state fluid dark energy model as
an interacting vacuum was prompted in part by its basis in pre-existing
physics, rather than being a purely phenomenological model of a vacuum – cold
dark matter interaction of the type previously considered in the literature.
This opens up an interesting philosophical question: should this kind of
argument be used more often when constructing alternative dark energy models?

On the one hand, simple models with a strong physical motivation are
preferable due to their elegance, and could yield interesting results, as we
have seen here. On the other hand, it is important not to become dogmatic when
using simplicity as a motivator. For example, a more “natural” (Weinberg,
1989) value of the cosmological constant may be preferable from an aesthetic
standpoint – but such a Universe would be inhospitable to life as we know it
(Dijkstra, 2019).

One may ask if the Shan–Chen models presented here are truly physically
motivated, or if they too fall under the label of phenomenological. Let us
note that while the concept of dark matter is nearly a century old, we still
do not have direct evidence for its existence besides the gravitational
effects which we attribute to it. Dark energy, either in the form of a
cosmological constant or something more exotic like the interacting vacuum
considered in this work, is even more puzzling. This leaves us with a huge
amount of freedom to construct models to describe these components. In other
words, while the fundamental physics behind dark matter and dark energy remain
unknown, the question of whether the Shan–Chen physical model remains valid
when transposed from its original context to that of cosmology also remains
unanswerable, even more so if the equation of state is turned into an
interaction, as we have done here. Nonetheless, the Shan-Chen modeling we have
introduced is physical at least in the sense that it is covariantly formulated
and as such is a description of matter (an interacting vacuum component in our
specific case) that at least in principle could be applied in different
contexts, e.g. neutron stars, unlike truly phenomenological parameterisations
such as CPL (Chevallier & Polarski, 2001; Linder, 2003) and other similar ad
hoc descriptions of evolving dark energy (Colgáin et al., 2021).

To conclude, the study presented here is a preliminary and non-exhaustive
examination of the Shan–Chen interacting vacuum, and many unexplored
combinations of parameters remain. In this analysis, we found that no
particular Shan–Chen interacting vacuum model of the cases studied here is
very competitive with $\Lambda$CDM when performing a model comparison.
Overall, our analysis shows it is generally difficult for these particular
interacting dark energy models to be both a better fit to the data than
$\Lambda$CDM and simultaneously resolve the tensions which persist in
cosmology.

## Acknowledgements

The Python libraries Matplotlib (Hunter, 2007) and NumPy (Harris et al., 2020)
were used for some of the visualisation and analysis presented in this paper.
NBH was supported by UK STFC studentship ST/N504245/1 for part of this work,
and partly by a postdoctoral position funded through two “la Caixa” Foundation
fellowships (ID 100010434), with fellowship codes LCF/BQ/PI19/11690015 and
LCF/BQ/PI19/11690018. NBH also received financial support from a G-Research
grant. MB is supported by UK STFC Grant No. ST/S000550/1. Numerical
computations were done on the Sciama High Performance Compute (HPC) cluster
which is supported by the ICG, SEPNet and the University of Portsmouth.

## Data Availability

The data underlying this article will be shared on reasonable request to the
corresponding author.

## CRediT authorship contribution statement

Natalie B. Hogg: Software, formal analysis, investigation, writing – original
draft, visualisation. Marco Bruni: Conceptualisation, formal analysis, writing
– review and editing, supervision.

## References

  * Abbott et al. (2020) Abbott T., et al., 2020, Physical Review D, 102, 023509
  * Aghanim et al. (2020) Aghanim N., et al., 2020, Astronomy & Astrophysics, 641, A6
  * Alam et al. (2017) Alam S., et al., 2017, Monthly Notices of the Royal Astronomical Society, 470, 2617
  * Amendola (2000) Amendola L., 2000, Physical Review D, 62, 043511
  * Ananda & Bruni (2006a) Ananda K. N., Bruni M., 2006a, Phys. Rev. D, 74, 023523
  * Ananda & Bruni (2006b) Ananda K. N., Bruni M., 2006b, Phys. Rev. D, 74, 023524
  * Armendariz-Picon et al. (2001) Armendariz-Picon C., Mukhanov V., Steinhardt P. J., 2001, Phys. Rev. D, 63, 103510
  * Balbi et al. (2007) Balbi A., Bruni M., Quercellini C., 2007, Phys. Rev. D, 76, 103519
  * Bamba et al. (2012) Bamba K., Capozziello S., Nojiri S., Odintsov S. D., 2012, Astrophys. Space Sci., 342, 155
  * Banerjee et al. (2021) Banerjee A., Cai H., Heisenberg L., Colgáin E. O., Sheikh-Jabbari M. M., Yang T., 2021, Phys. Rev. D, 103, L081305
  * Benevento et al. (2020) Benevento G., Hu W., Raveri M., 2020, Phys. Rev. D, 101, 103517
  * Bento et al. (2004) Bento M. C., Bertolami O., Sen A. A., 2004, Phys. Rev. D, 70, 083519
  * Bertacca et al. (2011) Bertacca D., Bruni M., Piattella O. F., Pietrobon D., 2011, Journal of Cosmology and Astroparticle Physics, 1102, 018
  * Beutler et al. (2011) Beutler F., et al., 2011, Monthly Notices of the Royal Astronomical Society, 416, 3017
  * Bini et al. (2013) Bini D., Geralico A., Gregoris D., Succi S., 2013, Phys. Rev. D, 88, 063007
  * Bini et al. (2016) Bini D., Esposito G., Geralico A., 2016, Phys. Rev. D, 93, 023511
  * Camarena & Marra (2021) Camarena D., Marra V., 2021, Mon. Not. Roy. Astron. Soc., 504, 5164
  * Capozziello et al. (2005) Capozziello S., Cardone V. F., Carloni S., De Martino S., Falanga M., Troisi A., Bruni M., 2005, JCAP, 04, 005
  * Chevallier & Polarski (2001) Chevallier M., Polarski D., 2001, Int. J. Mod. Phys. D, 10, 213
  * Colgáin et al. (2021) Colgáin E. O., Sheikh-Jabbari M. M., Yin L., 2021, Phys. Rev. D, 104, 023510
  * Copeland et al. (2006) Copeland E. J., Sami M., Tsujikawa S., 2006, International Journal of Modern Physics D, 15, 1753
  * De Felice et al. (2012) De Felice A., Nesseris S., Tsujikawa S., 2012, Journal of Cosmology and Astroparticle Physics, 05, 029
  * Di Valentino et al. (2021a) Di Valentino E., et al., 2021a, Class. Quant. Grav., 38, 153001
  * Di Valentino et al. (2021b) Di Valentino E., et al., 2021b, Astropart. Phys., 131, 102605
  * Dijkstra (2019) Dijkstra C. D., 2019, Master’s thesis, University of Groningen, arxiv.org/abs/1906.03036
  * Dodge (2008) Dodge Y., 2008, Chi-Square Table. Springer New York, New York, NY, pp 76–77, doi:10.1007/978-0-387-32833-1_56, https://doi.org/10.1007/978-0-387-32833-1_56
  * Efstathiou (2021) Efstathiou G., 2021, Mon. Not. Roy. Astron. Soc., 505, 3866
  * Eke et al. (1998) Eke V. R., Navarro J. F., Frenk C. S., 1998, The Astrophysical Journal, 503, 569
  * Freedman et al. (2019) Freedman W. L., et al., 2019, The Astrophysical Journal, 882, 34
  * Freese & Winkler (2021) Freese K., Winkler M. W., 2021, Chain Early Dark Energy: Solving the Hubble Tension and Explaining Today’s Dark Energy (arXiv:2102.13655)
  * Frenk & White (2012) Frenk C., White S. D., 2012, Annalen der Physik, 524, 507
  * Gao et al. (2010) Gao C., Kunz M., Liddle A. R., Parkinson D., 2010, Physical Review D, 81, 043520
  * Haridasu et al. (2021) Haridasu B. S., Viel M., Vittorio N., 2021, Phys. Rev. D, 103, 063539
  * Harris et al. (2020) Harris C. R., et al., 2020, Nature, 585, 357
  * Heymans et al. (2021) Heymans C., et al., 2021, Astron. Astrophys., 646, A140
  * Hill et al. (2020) Hill J. C., McDonough E., Toomey M. W., Alexander S., 2020, Phys. Rev. D, 102, 043507
  * Hogg et al. (2020) Hogg N. B., Bruni M., Crittenden R., Martinelli M., Peirone S., 2020, Phys. Dark Univ., 29, 100583
  * Howlett et al. (2012) Howlett C., Lewis A., Hall A., Challinor A., 2012, Journal of Cosmology and Astroparticle Physics, 1204, 027
  * Hunter (2007) Hunter J. D., 2007, Computing in Science & Engineering, 9, 90
  * Ivanov et al. (2020) Ivanov M. M., McDonough E., Hill J. C., Simonović M., Toomey M. W., Alexander S., Zaldarriaga M., 2020, Phys. Rev. D, 102, 103502
  * Kamenshchik et al. (2001) Kamenshchik A. Y., Moschella U., Pasquier V., 2001, Physics Letters B, 511, 265
  * Karwal et al. (2021) Karwal T., Raveri M., Jain B., Khoury J., Trodden M., 2021, Chameleon Early Dark Energy and the Hubble Tension (arXiv:2106.13290)
  * Kunz et al. (2016) Kunz M., Nesseris S., Sawicki I., 2016, Physical Review D, 94, 023510
  * Lewis (2013) Lewis A., 2013, Physical Review D, 87, 103529
  * Lewis (2019) Lewis A., 2019, GetDist: a Python package for analysing Monte Carlo samples (arXiv:1910.13970)
  * Lewis & Bridle (2002) Lewis A., Bridle S., 2002, Physical Review D, 66, 103511
  * Lewis et al. (2000) Lewis A., Challinor A., Lasenby A., 2000, The Astrophysical Journal, 538, 473
  * Li et al. (2018) Li H., Yang W., Wu Y., 2018, Physics of the Dark Universe, 22, 60
  * Li et al. (2019) Li H., Yang W., Gai L., 2019, Astronomy & Astrophysics, 623, A28
  * Linder (2003) Linder E. V., 2003, Phys. Rev. Lett., 90, 091301
  * Linder (2007) Linder E. V., 2007, General Relativity and Gravitation, 40, 329
  * Martin (2012) Martin J., 2012, Comptes Rendus Physique, 13, 566
  * Martinelli et al. (2019) Martinelli M., Hogg N. B., Peirone S., Bruni M., Wands D., 2019, Monthly Notices of the Royal Astronomical Society, 488, 3423
  * Niedermann & Sloth (2020) Niedermann F., Sloth M. S., 2020, Phys. Rev. D, 102, 063527
  * Park et al. (2021) Park M., Raveri M., Jain B., 2021, Phys. Rev. D, 103, 103530
  * Perlmutter et al. (1999) Perlmutter S., et al., 1999, The Astrophysical Journal, 517, 565
  * Pesce et al. (2020) Pesce D., et al., 2020, The Astrophysical Journal Letters, 891, L1
  * Piattella (2010) Piattella O. F., 2010, Journal of Cosmology and Astroparticle Physics, 03, 012
  * Piattella et al. (2010) Piattella O. F., Bertacca D., Bruni M., Pietrobon D., 2010, Journal of Cosmology and Astroparticle Physics, 1001, 014
  * Pietrobon et al. (2008) Pietrobon D., Balbi A., Bruni M., Quercellini C., 2008, Phys. Rev. D, 78, 083510
  * Poulin et al. (2019) Poulin V., Smith T. L., Karwal T., Kamionkowski M., 2019, Physical Review Letters, 122, 221301
  * Quercellini et al. (2008) Quercellini C., Bruni M., Balbi A., Pietrobon D., 2008, Phys. Rev. D, 78, 063527
  * Riess et al. (1998) Riess A. G., et al., 1998, The Astronomical Journal, 116, 1009
  * Riess et al. (2019) Riess A. G., Casertano S., Yuan W., Macri L. M., Scolnic D., 2019, The Astrophysical Journal, 876, 85
  * Ross et al. (2015) Ross A. J., Samushia L., Howlett C., Percival W. J., Burden A., Manera M., 2015, Mon. Not. Roy. Astron. Soc., 449, 835
  * Sakai (1999) Sakai S., 1999, The Tip of the Red Giant Branch as a Population II Distance Indicator. p. 48
  * Salvatelli et al. (2014) Salvatelli V., Said N., Bruni M., Melchiorri A., Wands D., 2014, Physical Review Letters, 113, 181301
  * Sandvik et al. (2004) Sandvik H., Tegmark M., Zaldarriaga M., Waga I., 2004, Physical Review D, 69, 123524
  * Scolnic et al. (2018) Scolnic D. M., et al., 2018, Astrophys. J., 859, 101
  * Shan & Chen (1993) Shan X., Chen H., 1993, Phys. Rev. E, 47, 1815
  * Sofue & Rubin (2001) Sofue Y., Rubin V., 2001, Annual Review of Astronomy and Astrophysics, 39, 137
  * Wands et al. (2012) Wands D., De-Santiago J., Wang Y., 2012, Classical and Quantum Gravity, 29, 145017
  * Wang et al. (2013) Wang Y., Wands D., Xu L., De-Santiago J., Hojjati A., 2013, Physical Review D, 87, 083503
  * Weinberg (1989) Weinberg S., 1989, Reviews of Modern Physics, 61, 1
  * Wong et al. (2020) Wong K. C., et al., 2020, Monthly Notices of the Royal Astronomical Society, 498
  * Zwicky (1933) Zwicky F., 1933, Helvetica Physica Acta, 6, 110
  * de Putter & Linder (2007) de Putter R., Linder E. V., 2007, Astroparticle Physics, 28, 263

##  Appendix A Additional results

In this appendix, we present two additional sets of results. We emphasise that
these are presented for the sake of completeness and are in line with the rest
of our findings and conclusions about the Shan–Chen interacting vacuum.

In Figure 11, we plot the results from sampling the cosmological parameters
when $q=0.05$ and $\beta=-0.05$ (light blue), and when $q=0.1$ and $\beta=0.2$
(pink). We can see that in both of these cases, the error bars on the
parameters (inferred from the 1 and 2 $\sigma$ regions of the contours) are
equivalent to those in $\Lambda$CDM (black), which is the same as what we find
in the main cases studied in this work when we keep both $q$ and $\beta$ fixed
(i.e. cases Ib, Ic and IIa). We can also see that in both of these cases, the
$H_{0}$ tension is largely unaffected, while the $\sigma_{8}$ tension is
somewhat relaxed. The cases shown here have a $\Delta\chi^{2}$ of 1.49 and
0.97 respectively, showing that they are both worse fits to the data than
$\Lambda$CDM.

Figure 11: Results of sampling the cosmological parameters while keeping
$q=0.05$ and $\beta=-0.05$ (light blue), and when $q=0.1$ and $\beta=0.2$
(pink). The $\Lambda$CDM result is shown in black.