# Critical point influenced by Bose-Einstein condensation

V. A. Kuznietsov  Bogolyubov Institute for Theoretical Physics, 03680 Kyiv,
Ukraine  Physics Department, Taras Shevchenko National University of Kyiv,
03022 Kyiv, Ukraine  O. Savchuk  GSI Helmholtzzentrum für Schwerionenforschung
GmbH, Planckstr. 1, D-64291 Darmstadt, Germany  Bogolyubov Institute for
Theoretical Physics, 03680 Kyiv, Ukraine  Frankfurt Institute for Advanced
Studies, Giersch Science Center, Ruth-Moufang-Str. 1, D-60438 Frankfurt am
Main, Germany  O. S. Stashko  Frankfurt Institute for Advanced Studies,
Giersch Science Center, Ruth-Moufang-Str. 1, D-60438 Frankfurt am Main,
Germany  M. I. Gorenstein  Bogolyubov Institute for Theoretical Physics, 03680
Kyiv, Ukraine  Frankfurt Institute for Advanced Studies, Giersch Science
Center, Ruth-Moufang-Str. 1, D-60438 Frankfurt am Main, Germany

(August 27, 2024)

###### Abstract

A system of bosons studied within the mean field framework has two fascinating
phenomena: a liquid-gas first order phase transition and Bose-Einstein
condensation. Interplay between these two phenomena is being investigated.
Depending on the mean-field potential parameters one can observe two types of
critical points, called ”Boltzmann” and ”Bose”, that belong to different
universality classes with distinct sets of critical exponents. As examples of
Bose and Boltzmann CPs pion and $\alpha$ matter are considered, respectively.
In general, the phase diagram can have one of the CPs or both of them.

Phase transition, mean-field model, Bose-Einstein condensation

##  I Introduction

When both attractive and repulsive interactions are present, systems of
classical particles demonstrate a liquid-gas first order phase transition
(FOPT), see, e.g., Ref. Landau and Lifshitz (1975). This transition ends with
a critical point (CP) at temperature $T=T_{c}$ where different thermodynamic
quantities show abnormal behavior specific for a second order phase
transition. The concept of critical exponents splits a great diversity of CP
phenomena into several universality classes with the same critical exponents
(see, e.g., Ref. Greiner  _et al._ (2012)).

Bose-Einstein condensation (BEC) in an ideal gas of bosons was predicted many
years ago Bose (1924); Einstein (1925) and experimentally confirmed for cold
atomic gases in magnetic traps Anderson  _et al._ (1995); Bradley  _et al._
(1995); Davis  _et al._ (1995); Dalfovo  _et al._ (1999). The theory of the
BEC phenomenon for interacting particles has been extensively discussed
Kapusta and Gale (2006); Andersen (2004); Griffin  _et al._ (1996); Watabe and
Ohashi (2013); Watabe (2019). In particular, modifications of the line of
onset of the BEC (further referred to as the BEC line) due to the small
repulsive interactions between particles were discussed Baym  _et al._ (1999a,
2000); Holzmann and Krauth (1999); Holzmann  _et al._ (2001); Savchuk  _et
al._ (2020a); Baym  _et al._ (1999b); Huang (1999). Connections between the
BEC and percolation phenomena are pointed out in Refs. Schakel (2001); Feynman
(1953a, b).

BEC phenomena have been studied in condensed matter physics, nuclear physics,
astrophysics, and cosmology (see, e.g., Refs. Satarov  _et al._ (2017); Begun
and Gorenstein (2007, 2008); Strinati  _et al._ (2018); Nozieres and Schmitt-
Rink (1985); Funaki  _et al._ (2008); Chavanis and Harko (2012); Mishustin
_et al._ (2019); Padilla  _et al._ (2019); Harko and Madarassy (2022)). In
most cases, particle interactions need to be taken into account. If the
attractive and repulsive forces are present, then both phenomena – the BEC and
a liquid-gas FOPT with a CP – are simultaneously expected in systems of
interacting bosons. An important property of the CP is infinitely high
particle number fluctuations (see Refs. Vovchenko (2017, 2020)). Some special
features of the mixed liquid-gas phase within the van der Waals model were
discussed in Refs. Satarov  _et al._ (2021) and Poberezhnyuk  _et al._ (2021).

The phase diagram of strongly interacting matter is one of the most important
problems in physics. A liquid-gas FOPT with formation of the Bose condensate
(BC) in the liquid phase at $T<T_{c}$ was considered for interacting $\alpha$
particles in Ref. Satarov  _et al._ (2017) and for interacting pions in Refs.
Stashko  _et al._ (2021a) and Kuznietsov  _et al._ (2021). It was also pointed
out Fedotkin  _et al._ (2019, 2021); Sanzhur (2022) that the CP temperature
$T_{c}$ of classical Boltzmann particles increases due to the Bose statistics
while the Fermi statistics leads to the decrease of $T_{c}$. The CP location
changes a by few percent due to the small quantum correction in $\alpha$
matter. On the other hand, in the system of interacting pions the location of
a CP found in Refs. Stashko  _et al._ (2021a) and Kuznietsov  _et al._ (2021)
lies just on the BEC line. Thus, much stronger Bose effects are expected in
the case of pions.

In the present paper we investigate the phase diagram of interacting Bose
particles within the mean-field model. The strength of the quantum statistics
effects in the vicinity of a CP depends on the value of the system parameters.
By changing these parameters we find two qualitatively different scenarios for
CP properties, which are denoted as the Boltzmann and Bose CPs. The critical
exponents for these CPs belong to different universality classes. For a
special region of model parameters, double phase transition with simultaneous
presence of the two CPs (Boltzmann and Bose), becomes possible. We find that
relativistic effects appear to be crucially important for these double phase
transitions.

The paper is organized as follows. The mean-field model results are presented
in Sec. II. This section includes the ground state properties, Boltzmann
approximation, and effects of the Bose statistics. In Sec. III the Boltzmann
CP and Bose CP are defined. The particle number fluctuations and the critical
exponents are calculated for the Bose CP and compared to those for the
Boltzmann CP. Scenarios with double phase transitions are considered in Sec.
IV. Summary in Sec. V closes the paper.

##  II mean-field model

The statistical system in the thermodynamic limit is defined in the grand
canonical ensemble in terms of the pressure $p$ as a function of temperature
$T$ and chemical potential $\mu$. All thermodynamic functions can then be
found from the pressure function using thermodynamic identities. The particle
number density $n$, entropy density $s$, and energy density $\varepsilon$ are
equal to

$\displaystyle n\equiv\left(\frac{\partial
p}{\partial\mu}\right)_{T}\leavevmode\nobreak\ ,\leavevmode\nobreak\
\leavevmode\nobreak\ s=\left(\frac{\partial p}{\partial
T}\right)_{\mu}\leavevmode\nobreak\ ,\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \varepsilon=Ts+\mu n-p\leavevmode\nobreak\ .$ (1)

For the ideal Bose gas one finds

$\displaystyle p_{\rm id}(T,\mu)$
$\displaystyle=\frac{g}{6\pi^{2}}\int^{\infty}_{0}f_{k}(T,\mu)\frac{k^{4}dk}{\sqrt{k^{2}+m^{2}}}\leavevmode\nobreak\
,$ (2) $\displaystyle n_{\rm id}(T,\mu)$
$\displaystyle\equiv\left(\frac{\partial p_{\rm
id}}{\partial\mu}\right)_{T}=\frac{g}{2\pi^{2}}\int^{\infty}_{0}f_{k}(T,\mu)k^{2}dk\leavevmode\nobreak\
,$ (3) $\displaystyle f_{k}(T,\mu)$
$\displaystyle=\left[\exp{\left(\frac{\sqrt{k^{2}+m^{2}}-\mu}{T}\right)-1}\right]^{-1}\leavevmode\nobreak\
,$ (4)

where $m$ and $g$ are the particle mass and degeneracy factor, respectively
111We use the units with $\hbar=c=1$. .

In relativistic systems the chemical potential $\mu$ regulates the conserved
charge, i.e., the number of particles minus the number of antiparticles. The
thermodynamic functions of antiparticles correspond to a substitution of
$\mu\rightarrow-\leavevmode\nobreak\ \mu$. In what follows we are mostly
interested in the region of large chemical potentials, $\mu/T>1$. Thus, a
contribution from antiparticles to the system thermodynamics can be neglected.
The number of particles plays then the role of the conserved charge. To have
closed form expressions we also often exploit the validity of the non-
relativistic approximation, $\sqrt{k^{2}+m^{2}}\cong m+k^{2}/(2m)$, accurate
for $T/m\ll 1$.

The mean-field model of interacting bosons is given by the following set of
self-consistent equations (see, e.g., Ref. Anchishkin and Vovchenko (2015)):

$\displaystyle p(T,\mu)$ $\displaystyle=p_{\rm
id}(T,\mu^{*})+\int^{n}_{0}\frac{dU}{dn^{\prime}}n^{\prime}dn^{\prime}\leavevmode\nobreak\
,$ (5) $\displaystyle n(T,\mu)$ $\displaystyle=n_{\rm id}(T,\mu^{*})+n_{\rm
BC}\leavevmode\nobreak\ ,$ (6) $\displaystyle\mu^{*}$
$\displaystyle=\mu-U(n)\leavevmode\nobreak\ ,$ (7)

where the density dependent mean-field potential $U$ will be taken in the
simple form

$\displaystyle U(n)=-An+Bn^{2},\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ A>0,\leavevmode\nobreak\ \leavevmode\nobreak\
B>0\leavevmode\nobreak\ .$ (8)

The constants $A$ and $B$ correspond to the attractive and repulsive
interactions, respectively. The quantity $n_{\rm BC}\geq 0$ in Eq. (6) is the
density of the BC. It corresponds to particles at the zero-momentum state
$k=0$. Non-zero values $n_{\rm BC}>0$ can only appear at $\mu^{*}=m$. Note
that values of $\mu^{*}>m$ are forbidden as they lead to negative values of
particle numbers $f_{k}$ at small $k$. The relation

$\displaystyle\mu^{*}(T,n)\leavevmode\nobreak\ =\leavevmode\nobreak\ m-0$ (9)

defines the BEC line in the $(n,T)$ plane. At this line, one still has $n_{\rm
BC}=0$, and it corresponds to an onset of the BE condensation.

###  II.1 Ground state at $T=0$

At $T\rightarrow 0$ the thermal pressure (2) and particle number density (3)
go to zero, so the system can only exist in a form of the BC. The condition of
the BEC, $\mu^{*}=m$, can be then rewritten as

$\displaystyle Bn^{2}_{\rm BC}-An_{\rm BC}-(\mu-m)=0\leavevmode\nobreak\ .$
(10)

This yields the BC density $n_{\rm BC}$ as a function of the chemical
potential $\mu$

$\displaystyle n_{\rm BC}(T=0,\mu)$
$\displaystyle=\frac{A+\sqrt{A^{2}+4B(\mu-m)}}{2B}\leavevmode\nobreak\
\leavevmode\nobreak\ .$ (11)

At $T=0$ and

$\displaystyle\mu\leavevmode\nobreak\ <\leavevmode\nobreak\
\mu_{0}=m-\frac{3A^{2}}{16B}$ (12)

the stable thermodynamic solution corresponds to $n_{\rm BC}=0$ and $p=0$, and
at $\mu\geq\mu_{0}$ the pressure behaves as

$\displaystyle p=-\frac{An^{2}_{\rm BC}}{2}+\frac{2Bn^{3}_{\rm BC}}{3}\geq
0\leavevmode\nobreak\ .$ (13)

The first-order phase transition for $T=0$ takes place at $\mu=\mu_{0}$ (see,
e.g., Ref. Kuznietsov  _et al._ (2021)) with jumps of $n$ and $\varepsilon$
from their zero values to

$\displaystyle n_{\rm BC}(T=0,\mu=\mu_{0})=\frac{3A}{4B}\leavevmode\nobreak\
,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ $ (14)
$\displaystyle\varepsilon(T=0,\mu=\mu_{0})=\mu_{0}n_{\rm
BC}=\frac{3A}{4B}\left(m-\frac{3A^{2}}{16B}\right)\leavevmode\nobreak\ .$ (15)

The quantities (14) and (15) are referred as the ground state values of the
considered system.Note that at all values of $\mu$ they are in agreement with
the third law of thermodynamics.

###  II.2 Boltzmann approximation

The Boltzmann approximation corresponds to a change of the $f_{k}$ function
(4) to

$\displaystyle f_{k}(T,\mu)\leavevmode\nobreak\ =\leavevmode\nobreak\
\exp\left(-\leavevmode\nobreak\
\frac{\sqrt{k^{2}+m^{2}}-\mu}{T}\right)\leavevmode\nobreak\ .$ (16)

The system pressure in the mean-field model in then reduced to a simple
analytical form Vovchenko  _et al._ (2020); Savchuk  _et al._ (2020b), namely

$\displaystyle
p(T,n)=Tn-\frac{An^{2}}{2}+\frac{2Bn^{3}}{3}\leavevmode\nobreak\ .$ (17)

The classical equation of state (17) depends on the interaction parameters $A$
and $B$, but is insensitive to the particle mass $m$ and degeneracy factor
$g$. The statistical system exhibits the liquid-gas FOPT with a CP defined by
the following conditions

$\displaystyle\left(\frac{\partial p}{\partial
n}\right)_{T}=0\leavevmode\nobreak\ ,\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \left(\frac{\partial^{2}p}{\partial
n^{2}}\right)_{T}=0\leavevmode\nobreak\ .$ (18)

From Eqs.(17) and (18) one finds CP parameters

$\displaystyle T_{\rm c}^{0}=\frac{A^{2}}{8B}\leavevmode\nobreak\
,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ n_{\rm
c}^{0}=\frac{A}{4B}\leavevmode\nobreak\ ,\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
p_{c}^{0}=\frac{A^{3}}{96B^{2}}\leavevmode\nobreak\ .$ (19)

In the reduced variables $\tilde{T}=T/T_{c}^{0}$, $\tilde{n}=n/n_{c}^{0}$, and
$\tilde{p}=p/p_{c}^{0}$, Eq. (17) takes a universal form,

$\displaystyle\tilde{p}\leavevmode\nobreak\ =\leavevmode\nobreak\
3\leavevmode\nobreak\ \tilde{T}\leavevmode\nobreak\
\tilde{n}\leavevmode\nobreak\ -\leavevmode\nobreak\ 3\leavevmode\nobreak\
\tilde{n}^{2}\leavevmode\nobreak\ +\leavevmode\nobreak\
\tilde{n}^{3}\leavevmode\nobreak\ ,$ (20)

independent of the parameters $A$ and $B$.

Figure 1:  (a) Phase diagram of the mean-field model in the $(n,P)$ plane
within the Boltzmann approximation. Solid lines correspond to the stable parts
of the pressure isotherms while dashed and dotted lines show their metastable
and unstable parts, respectively. (b) Same phase diagram in the $(n,T)$ plane.

Figure 1 demonstrates several pressure isotherms
$\tilde{p}=\tilde{p}(\tilde{n})$ at $\tilde{T}\leq 1$. These isotherms include
the metastable and unstable parts denoted in Fig. 1 (a) by dashed and dotted
lines, respectively. A region of the unstable isotherms with
$\partial\tilde{p}(\tilde{n},\tilde{T})/\partial\tilde{n}<0$ shown with grey
color in Fig. 1 is restricted by the left and right spinodals where
$\partial\tilde{p}(\tilde{n},T)/\partial\tilde{n}=0$. In the
$(\tilde{n},\tilde{T})$-regions with
$\partial\tilde{p}(\tilde{n},\tilde{T})/\partial\tilde{n}<0$ the system can
not exist in the homogeneous state. Instead, it behaves as an inhomogeneous
mixture of the gas and liquid with particle number densities $n_{\rm g}$ and
$n_{\rm l}$, respectively. According to the Gibbs criteria these two phases
should have the same temperatures, chemical potentials, and pressures. A
region occupied by the liquid-gas mixed phase is restricted by the left
(gaseous) and right (liquid) binodals with particle number densities $n_{\rm
g}$ and $n_{\rm l}$, respectively. These binodals are shown by the blue solid
lines in Fig. 1 (b). The pressure isotherms inside the mixed phase shown by
the horizontal solid lines in Fig. 1 (a) are obtained from the Gibbs
conditions or, equivalently, with the Maxwell construction of equal areas
Landau and Lifshitz (1975).

Critical exponents define the behavior of the thermodynamic function in the
vicinity of a CP. The critical exponents $\alpha,\beta,\gamma$, and $\delta$
are defined as:

$\displaystyle\left.\left(\frac{\partial\varepsilon}{\partial
T}\right)\right|_{n=n_{\rm c}}\leavevmode\nobreak\ $
$\displaystyle\sim\leavevmode\nobreak\ |T_{\rm
c}-T|^{-\alpha}\leavevmode\nobreak\ ,$ (21) $\displaystyle n_{\rm
l}\leavevmode\nobreak\ -\leavevmode\nobreak\ n_{\rm g}\leavevmode\nobreak\ $
$\displaystyle\sim\leavevmode\nobreak\ \left(T_{\rm
c}-T\right)^{\beta}\leavevmode\nobreak\ ,$ (22)
$\displaystyle\left[\left(\frac{\partial p}{\partial
n}\right)_{T}\right]^{-1}_{n=n_{\rm c}}\leavevmode\nobreak\ $
$\displaystyle\sim\leavevmode\nobreak\ |T_{c}-T|^{-\gamma}\leavevmode\nobreak\
,$ (23) $\displaystyle\left.\left(p-p_{\rm c}\right)\right|_{T=T_{\rm
c}}\leavevmode\nobreak\ $ $\displaystyle\sim\leavevmode\nobreak\ |n-n_{\rm
c}|^{\delta}\leavevmode\nobreak\ {\rm sgn}(n-n_{\rm c})\leavevmode\nobreak\ .$
(24)

In the Boltzmann approximation (17) one can easily find the following critical
exponents in the mean-field model

$\displaystyle\alpha\leavevmode\nobreak\ =\leavevmode\nobreak\
0\leavevmode\nobreak\ ,\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \beta\leavevmode\nobreak\ =\leavevmode\nobreak\
\frac{1}{2}\leavevmode\nobreak\ ,\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \gamma\leavevmode\nobreak\ =\leavevmode\nobreak\
\leavevmode\nobreak\ 1\leavevmode\nobreak\ ,\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \delta\leavevmode\nobreak\
=\leavevmode\nobreak\ 3\leavevmode\nobreak\ .$ (25)

In the present paper we use the Skyrme-like mean-field potential $U(n)$ (8).
Its simple form makes possible straightforward analytical calculations for the
first order liquid-gas phase transition with the end CP for classical
particles, i.e., in the absence of the Bose statistics. The CP in the
Boltzmann approximation with its critical exponents (25) belongs to the
universality class of the the so called ”classical models” (or the van der
Waals type models). This approach gives however only an approximate
description of the liquid-gas CP in real systems. A more realistic description
of the CP phenomena based on the renormalization group methods can be found in
Kos  _et al._ (2016) and references there in.

The Boltzmann approximation is not valid at very low temperatures
$T\rightarrow 0$ and finite particle number density $n$. The Boltzmann gas
entropy $s=s(T,n)$ in this limit becomes negative in a contradiction with the
third law of thermodynamics. Thus the quantum statistics is needed to describe
correctly the ground state properties of the system at $T=0$.

###  II.3 Effects of Bose statistics

The thermodynamic properties within considered mean-field model depend on four
parameters: $A$, $B$, $m$, and $g$. The strength of the Bose effects in a
vicinity of a CP depends on the system parameters. We consider first the small
quantum statistics corrections (see Refs. Fedotkin  _et al._ (2019, 2021)).
The pressure (5) and particle number density (6) can be presented as

$\displaystyle p=\alpha\leavevmode\nobreak\
T\sum_{k=1}^{\infty}K_{2}\left(\frac{km}{T}\right)\frac{\exp(k\mu^{*}/T)}{k^{2}}-\frac{An^{2}}{2}+\frac{2Bn^{3}}{3},$
(26) $\displaystyle
n=\alpha\sum_{k=1}^{\infty}K_{2}\left(\frac{km}{T}\right)\frac{\exp(k\mu^{*}/T)}{k}\leavevmode\nobreak\
+\leavevmode\nobreak\ n_{\rm BC}\leavevmode\nobreak\ ,$ (27)

where $\alpha\equiv gTm^{2}/(2\pi^{2})$ and $K_{2}$ is the modified Bessel
function. The BEC line (9) obtained from (27) at $\mu^{*}\rightarrow m-0$ and
$n_{\rm BC}=0$ reads

$\displaystyle n\leavevmode\nobreak\ =\leavevmode\nobreak\
\alpha\sum_{k=1}^{\infty}K_{2}\left(\frac{km}{T}\right)\leavevmode\nobreak\
\exp\left(\frac{k\leavevmode\nobreak\ m}{T}\right).$ (28)

Note that in the mean-field model this line in the $(n,T)$ plane is identical
to that in the ideal Bose gas. At low temperatures $T\ll m$, one can use the
non-relativistic approximation, i.e.,
$K_{2}(x)\approx(\pi/(2x))^{1/2}\exp(-x)$, at $x\gg 1$. It leads to the well
known text-book expression for the BEC line Landau and Lifshitz (1975):

$\displaystyle T_{\rm BC}\leavevmode\nobreak\ =\leavevmode\nobreak\
\frac{2\pi}{m}\leavevmode\nobreak\
\left(\frac{n}{g\zeta(3/2)}\right)^{2/3}\leavevmode\nobreak\
,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \frac{T_{\rm BC}}{m}\ll 1.$ (29)

In the ultra-relativistic limit $T/m\gg 1$, one uses $K_{2}(x)\approx 2/x^{2}$
at $x\ll 1$ and finds Begun and Gorenstein (2008)

$\displaystyle T_{\rm
BC}=\left[\frac{\pi^{2}n}{g\zeta(3)}\right]^{1/3}\leavevmode\nobreak\
,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \frac{T_{\rm
BC}}{m}\gg 1\leavevmode\nobreak\ .$ (30)

In Eqs. (29) and (30), $\zeta(l)$ is the Riemann zeta function,
$\zeta(3/2)\approx 2.6$ and $\zeta(3)\approx 1.2$.

In the non relativistic limit and $z\equiv\exp[(\mu^{*}-m)/T]\ll 1$, one finds

$\displaystyle p$ $\displaystyle\approx
gT\left(\frac{mT}{2\pi}\right)^{3/2}\left[z+\frac{z^{2}}{2^{5/2}}\right]\leavevmode\nobreak\
-\frac{An^{2}}{2}+\frac{2Bn^{3}}{3}\leavevmode\nobreak\ ,$ (31) $\displaystyle
n$ $\displaystyle\approx
g\left(\frac{mT}{2\pi}\right)^{3/2}\left[z+\frac{z^{2}}{2^{3/2}}\right]\leavevmode\nobreak\
.$ (32)

The Boltzmann approximation corresponds to $z\rightarrow 0$ when only the
first terms $k=1$ contribute to power series in Eq. (26) and (27).

Taking into account the next $k=2$ terms one obtains Eqs. (31) and (32) with
corrections due to the Bose statistics that remain small at $z\ll 1$. By
inverting $n(z)$ to $z=z(n)$ in Eq. (32) and substituting it into (31) one
finds the pressure function

$\displaystyle p(T,n)\approx\leavevmode\nobreak\ Tn-\left(\frac{A}{2}+\gamma
T^{-1/2}\right)n^{2}+\frac{2B}{3}n^{3}\leavevmode\nobreak\ ,$ (33)

where $\gamma\equiv\pi^{3/2}(2gm^{3/2})^{-1}$.

The additional $n^{2}$ term in the pressure (33) in comparison to its
Boltzmann approximation (17) comes from the (small) effects of the Bose
statistics. It corresponds effectively to an increase of the attractive
interactions. Thus, one expects an increase of the CP parameters
$n_{c}^{0}\sim A$ and $T_{c}^{0}\sim A^{2}$. Indeed, Eqs. (18) for the
pressure function (33) read

$\displaystyle T-(A+2\gamma T^{-1/2})n+2Bn^{2}$
$\displaystyle=0\leavevmode\nobreak\ ,$ (34) $\displaystyle-(A+2\gamma
T^{-1/2})+4Bn$ $\displaystyle=0\leavevmode\nobreak\ .$ (35)

Solving the system of Eqs. (34) and (35) for the CP parameters, one obtains

$\displaystyle n^{1}_{\rm c}$
$\displaystyle=\frac{n_{c}^{0}}{2}\left(\leavevmode\nobreak\
1\leavevmode\nobreak\ +\leavevmode\nobreak\ \sqrt{1\leavevmode\nobreak\
+\leavevmode\nobreak\ 4\sigma}\leavevmode\nobreak\ \right)\approx n^{0}_{\rm
c}(1+\sigma),$ (36) $\displaystyle T_{\rm c}^{1}$
$\displaystyle\leavevmode\nobreak\ =\leavevmode\nobreak\
2B\leavevmode\nobreak\ (n_{\rm c}^{1})^{2}\leavevmode\nobreak\
\approx\leavevmode\nobreak\ T_{c}^{0}(1+2\sigma)\leavevmode\nobreak\ ,$ (37)

where

$\displaystyle\sigma\leavevmode\nobreak\ \equiv\leavevmode\nobreak\
(2\pi)^{3/2}\leavevmode\nobreak\ \frac{B^{1/2}}{A^{2}}\leavevmode\nobreak\
g^{-1}\leavevmode\nobreak\ m^{-3/2}\leavevmode\nobreak\
\equiv\leavevmode\nobreak\
(2\pi)^{3/2}\frac{\tilde{B}^{1/2}}{\tilde{A}^{2}}\leavevmode\nobreak\ ,$ (38)

with the dimensionless parameters

$\displaystyle\tilde{A}\equiv gm^{2}\,A\leavevmode\nobreak\
,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \tilde{B}\equiv g^{2}m^{5}\,B\leavevmode\nobreak\ .$ (39)

Note that Eqs. (36) and (37) are numerically accurate at small $\sigma\ll 1$.

Equations (36) and (37) demonstrate an increase of $n_{c}^{1}$ and $T_{c}^{1}$
in comparison with their $n_{c}^{0}$ and $T_{c}^{0}$ values obtained within
the Boltzmann approximation. The size of these Bose effects is regulated by a
single parameter $\sigma$ (38). This parameter increases with the repulsive
interactions, $\sigma\sim B^{1/2}$, and decreases with $A$, $g$, and $m$ as
$\sigma\sim A^{-2}g^{-1}m^{-3/2}$.

To illustrate the role of Bose effects, two examples will be discussed. The
first example concerns a system of interacting $\alpha$ particles
($g=1,m=m_{\alpha}=3727$ MeV). Following Ref. Satarov  _et al._ (2017) we fix
the ground state properties of the $\alpha$ matter as

$\displaystyle n_{0}=\frac{3A}{4B}=0.036\leavevmode\nobreak\ {\rm
fm}^{-3}\leavevmode\nobreak\ ,\leavevmode\nobreak\ \leavevmode\nobreak\
W_{0}=\frac{3A^{2}}{16B}=12\leavevmode\nobreak\ {\rm MeV}.$ (40)

From Eq. (40) one finds $A\cong 1.35$ GeV fm3 and $B\cong 28$ GeV fm6. This
leads to the CP parameters presented in Table 1.

The phase diagram of $\alpha$ matter is shown in Fig. 2 (a) in terms of the
variables $\tilde{n}=n/n_{c}^{0}$ and $\tilde{T}=T/T_{c}^{0}$. One observes
rather small Bose effects near the CP ($\sigma\approx 0.049\ll 1$).

Within the reduced variables $\tilde{n}$ and $\tilde{T}$ the phase diagram is
defined in terms of the single parameter $\sigma\ll 1$. As seen from Table 1,
the values of $n_{\rm c}^{1}$ and $T_{\rm c}^{1}$ obtained in Eqs. (36) and
(37) give a good approximations to their exact values $n_{\rm c}$ and $T_{\rm
c}$.

$T^{0}_{\rm c}$, MeV | $n^{0}_{\rm c}\ {\rm fm}^{-3}$ | $T^{1}_{\rm c}$ | $n^{1}_{\rm c}$ | $T_{\rm c}$ | $n_{\rm c}$ | $\sigma$
---|---|---|---|---|---|---
$8.0$ | $0.012$ | $8.8$ | $0.0126$ | $8.9$ | $0.013$ | $0.049$
Table 1: Results for the CP parameters for $\alpha$ matter: $m=3.73$ GeV, $g=1$, $A=1.35$ GeV fm3, and $B=28$ GeV fm6.  $T^{0}_{\rm c}$, MeV | $n^{0}_{\rm c}\ {\rm fm}^{-3}$ | $T_{\rm c}$ | $n_{\rm c}$ | $\sigma$
---|---|---|---|---
$0.034$ | $0.0035$ | $29.2$ | $0.007$ | $7\cdot 10^{3}$
Table 2: Results for the CP parameters for the pion matter: $m=0.14$ GeV,
$g=1$, $A=0.0196\ {\rm GeV\cdot fm}^{3}$, and $\ B=1.426\ {\rm GeV\cdot
fm}^{6}$ Figure 2: The phase diagram of $\alpha$ matter (a) and pion matter
(b) in terms of the reduced variables $n/n_{c}^{0}$ and $T_{c}/T_{c}^{0}$.
Dashed lines present the BEC lines and solid lines are the gaseous and liquid
binodals. The Boltzmann CP and Bose CP are shown by yellow and blue stars,
respectively, the triple point is shown by an open circle, and red stars
correspond to the CP in the Boltzmann approximation.

Our second example corresponds to the interacting pion matter
($g=1,m=m_{\pi}=140$ MeV). We take $A=19.6\ {\rm MeV\,fm}^{3}$ and $\ B=1426\
{\rm MeV\,fm}^{6}$ obtained in Ref. Kuznietsov  _et al._ (2021) from fitting
the lattice data for pion Bose condensate at $T=0$ as a function of electric
chemical potential. The results for $T_{\rm c}$ and $n_{\rm c}$ values are
presented in Table 2, and the $(n,T)$ phase diagram is shown in Fig. 2 (b). In
contrast to $\alpha$ matter, the effects of Bose statistics for pions are very
large. One finds $\sigma\gg 1$ and thus $T_{\rm c}/T_{\rm c}^{0}\gg 1$.
Approximations (36) and (37) lose their validity at large $\sigma\gg 1$, and
thus the $n_{\rm c}^{1}$ and $T_{\rm c}^{1}$ values have no physical meaning.
The BEC line goes through the CP in pion matter, as seen in Fig. 2(b). One
finds a limiting value of the critical particle number density,

$\displaystyle n_{c}^{\rm lim}\leavevmode\nobreak\ =\leavevmode\nobreak\
\frac{A}{2B}\leavevmode\nobreak\ =\leavevmode\nobreak\ 2n_{c}^{0}$ (41)

that is exactly two times larger than $n_{\rm c}^{0}$ obtained in the
Boltzmann approximation. Note that the ground state BC density given by Eq.
(14) equals $n_{\rm BC}^{0}=3n_{c}^{0}$. The limiting critical density
$n_{c}^{\rm lim}$ (41) corresponds to the $n$ value where $dU/dn=0$, and it is
independent of parameters $m$ and $g$.

On the other hand, the $T_{c}$ value at $n_{c}=n_{\rm c}^{\rm lim}$ is still
dependent on $m$ and $g$. One finds

$\displaystyle T_{\rm c}=T_{\rm BC}\leavevmode\nobreak\
\approx\leavevmode\nobreak\ \frac{2\pi}{m}\leavevmode\nobreak\
\left[\frac{n_{\rm c}^{\rm lim}}{g\leavevmode\nobreak\
\zeta(3/2)}\right]^{2/3}\leavevmode\nobreak\ ,\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ T/m\ll 1\leavevmode\nobreak\ .$ (42)

Interacting $\alpha$ matter and pion matter are considered within the same
mean-field model. We are interested in possible changes of the system
parameters $A$, $B$, $m$, and $g$ in the both systems. Our goal is to study an
interplay between the interaction parameters $A$ and $B$, and the Bose
statistics effects that are also sensitive to the particle mass $m$ and
degeneracy factor $g$. In Ref. Satarov  _et al._ (2017) the mean-field
potential $U(n)$ for $\alpha$ matter was considered as a function of three
parameters: $a$, $n_{0}$, and $\gamma$. In the present study we use a
particular case of that potential with $\gamma=1,\leavevmode\nobreak\
\leavevmode\nobreak\ n_{0}=3A/(4B)$, and $a=A/2$. This simplifies the model
consideration and the both $\alpha$ matter and pion matter are considered
within the same model. Note also that pure $\alpha$ matter is not a realistic
system. A minimal extension of this model would require one to add the
interacting nucleons and some light nuclei like $d$, $t$, etc. The presence of
these additional particles in the nuclear mixture can strongly affects the
$\alpha$ particles Röpke (2009). We do not insist on the existence of pure
$\alpha$ matter in nature. The pure $\alpha$ matter in our consideration is
only a useful and well studied toy model example to illustrate a limiting case
of the Boltzmann-like CP.

##  III Boltzmann and Bose CP regimes

The CPs shown in Figs. 2(a) and 2(b) correspond to two limiting cases of small
and large effects of quantum statistics. In what follows we denote these two
scenarios as the Boltzmann CP and Bose CP, respectively. The reduced variables
$n/n^{0}_{\rm c}$ and $T/T^{0}_{\rm c}$ appear to be rather useful for the
analysis. As discussed in Appendix A all possible sets of the model parameters
$A$, $B$, $m$, and $g$ can be combined to the two dimensionless parameters
$\tilde{A}$ and $\tilde{B}$ (39) which define the type of the CP.

In the non-relativistic approximation this is further simplified to a
dependence on the single dimensionless parameter $\sigma$ (38). At fixed $A$
and $B$ the value of $\sigma$ remains a function of $m$ and $g$. One can lead
the pion system to a region of small values of $\sigma$ by artificially
increasing the pion values of $m$ and/or $g$ with fixed $A$ and $B$ values for
pions. Similarly one can lead the $\alpha$ system to a region of large values
of $\sigma$ by artificially decreasing of the $m$ value of $\alpha$ with fixed
$A$ and $B$ values for $\alpha$ matter. These procedures will move the system
parameters for pions and $\alpha$ particles between the regions of small and
large values of $\sigma$.

The phase boundaries and their CPs for different sets of model parameters
obtained by changing the $m$ and $g$ values for pion and $\alpha$ matter are
presented for several values of $\sigma$ in Fig. 3(a). The $\sigma$ value work
as a single universal parameter. A common feature for these physical systems
is a validity of the non-relativistic approximation, i.e., $T_{\rm c}/m\ll 1$
in all these cases. A transition between Boltzmann CP and Bose CP regimes
takes place at $\sigma=\sigma_{\rm cr}\cong 1.62$. This is discussed in
Appendix B.

In Figs. 2(a) and 3(a) open circles denote the triple points, namely the
points of intersections of the BEC line with the liquid binodals at $T_{\rm
tr}<T_{\rm c}$. The liquid component of the mixed phase at $T<T_{\rm tr}$
includes the BC, $n_{\rm BC}>0$, which affects part of the liquid binodal. In
Ref. Stein  _et al._ (1998) a notion of “binodal anomaly” was introduced to
indicate a special form of liquid binodal. It is a result of the influence of
the Bose statistics effects. This corresponds to $T_{\rm tr}<T_{\rm c}$ in our
notation, in a case when the CP is still the Boltzmann type one. When $T_{\rm
tr}=T_{\rm c}$ the BEC affects the whole liquid binodal. In addition, this has
anomalous influence to the end point of the binodal: the Boltzmann CP is
transformed to the Bose CP.

In Fig. 3(b) the plane of the parameters $\tilde{A}$ and $\tilde{B}$ is
presented. The solid red line of $\sigma_{\rm cr}\cong 1.62$, defines the
boundary between the Boltzmann and Bose regimes. It is calculated in the non-
relativistic approximation. The original sets of the $\pi$ and $\alpha$
parameters are shown in the $(\tilde{A},\tilde{B})$ plane by the corresponding
symbols in Fig. 3 (b). The dotted lines in Fig. 3 (b) illustrate the changes
of the pion system parameters by increasing $m$ or $g$ and the changes of the
$\alpha$ system parameters by decreasing $m$. For the pion and $\alpha$ matter
parameters $A$ and $B$ the transitions between the Boltzmann and Bose CPs at
$\sigma\cong 1.62$ take place at $T_{\rm c}/m\ll 1$, thus, the non-
relativistic approximation is valid.

Figure 3: (a) Phase diagrams at different values of the parameter $\sigma$.
(b) The $(\tilde{A},\tilde{B})$ diagram calculated in the non relativistic
approximation. The red solid line shows the boundary $\sigma_{\rm cr}=1.62$
between two regimes (see Appendix B). Pion matter and alpha matter are
indicated with symbols $\pi$ and $\alpha$. Dotted lines with arrows illustrate
pion matter with increasing $m$ and $g$ values, and $\alpha$ matter with
decreasing $m$ value.

###  III.1 Fluctuations

The particle number $N$ fluctuates in the grand canonical ensemble. These
fluctuations can be characterized by the scaled variance

$\displaystyle\omega\leavevmode\nobreak\ $
$\displaystyle\equiv\leavevmode\nobreak\ \frac{\langle N^{2}\rangle-\langle
N\rangle^{2}}{\langle N\rangle}\leavevmode\nobreak\ =\leavevmode\nobreak\
T\leavevmode\nobreak\ \left[\frac{\partial p(T,n)}{\partial
n}\right]_{T}^{-1}$ $\displaystyle=\omega_{\rm
id}(T,\mu^{*})\left[1+\frac{n}{T}\leavevmode\nobreak\
\frac{dU}{dn}\leavevmode\nobreak\ \omega_{\rm
id}(T,\mu^{*})\right]^{-1}\leavevmode\nobreak\ ,$ (43)

where

$\displaystyle\omega_{\rm id}(T,\mu^{*})\leavevmode\nobreak\
=\leavevmode\nobreak\ \frac{T}{n}\leavevmode\nobreak\ \left(\frac{\partial
n}{\partial\mu^{*}}\right)_{T}\leavevmode\nobreak\ $ (44)

is the scaled variance for the ideal Bose gas.

In the Boltzmann approximation, it follows that $\omega_{\rm id}=1$, and one
obtains from (III.1) the simple analytical expression

$\displaystyle\omega_{0}(T,n)=\left[1+\frac{n}{T}U^{\prime}(n)\right]^{-1}=\left[1-\frac{2\tilde{n}-\tilde{n}^{2}}{\tilde{T}}\right]^{-1}.$
(45)

The gaseous and liquid spinodals in the Boltzmann approximation shown in Fig.
1 are defined by the equation

$\displaystyle\left[\frac{\partial\tilde{p}}{\partial{\tilde{n}}}\right]\leavevmode\nobreak\
=\leavevmode\nobreak\ \tilde{n}^{2}-2\tilde{n}+\tilde{T}\leavevmode\nobreak\
=\leavevmode\nobreak\ 0\leavevmode\nobreak\ .$ (46)

Therefore, along the spinodals $\omega_{0}=\infty$, and $\omega_{0}<0$ in the
unstable $(n,T)$ region between spinodals.

As seen from Eq. (III.1), the infinite values of $\omega$ emerge in the
Boltzmann CP where

$\displaystyle\frac{T_{\rm c}}{n_{\rm c}}\leavevmode\nobreak\
=\leavevmode\nobreak\ -\leavevmode\nobreak\ \leavevmode\nobreak\
\frac{dU}{dn}\leavevmode\nobreak\ \omega_{\rm id}\leavevmode\nobreak\ ,$ (47)

i.e., this requires $dU/dn<0$. In the Boltzmann CP at $n_{c}<2n_{c}^{0}$ the
infinite values of $\omega$ appear due to the interplay between repulsive and
attractive interactions between particles. The presence of attractive
interactions leads to $dU/dn<0$ at $n=n{\rm{}_{c}}<2n_{\rm c}^{0}$. Equation
(47) then leads to $\omega\rightarrow\infty$ at the CP $(n_{c},T_{c})$. This
mechanism does not require Bose statistics, and it takes place already within
the Boltzmann approximation, i.e., for $\omega_{\rm id}=1$.

The infinite fluctuations in the Bose CP at $n_{\rm c}=2n_{\rm c}^{0}$ have a
different origin. On the BEC line one finds Begun and Gorenstein (2008)
$\omega_{\rm id}(T,\mu^{*}\rightarrow m-0)\rightarrow\infty$ and therefore

$\displaystyle\omega=\frac{T}{n}\leavevmode\nobreak\
\left(\frac{dU}{dn}\right)^{-1}\leavevmode\nobreak\ .$ (48)

The particle number fluctuations are infinite, $\omega=\infty$, in the Bose CP
at $n_{\rm c}=2n_{\rm c}^{0}$ because of $dU/dn=0$. At $n>2n_{\rm c}^{0}$,
$dU/dn>0$ and Eq. (48) lead to finite values of $\omega$ on the BEC line.

The scaled variances $\omega$ for the $\alpha$ matter and pion matter are
shown in Figs. 4 (a) and (b), respectively, on the $(\mu,T)$ plane. They give
examples of the Boltzmann and Bose CP. In both cases, infinite values of the
scaled variance appear at the CP. However, the regions of large fluctuations
look rather different. For the Boltzmann CP shown in Fig. 4 (a) the large
values of $\omega$ are localized in a narrow region of the $(\mu,T)$ plane
that looks like a continuation of the FOPT line. The region with the BEC does
not generate large fluctuations. For the Bose CP shown in Fig. 4 (b) a region
of large fluctuations is much wider. It is localized under the BEC line and
the FOPT line in the vicinity of the Bose CP.

###  III.2 Critical exponents

The Bose CP corresponds to $n_{\rm c}=2n_{\rm c}^{0}$ and $T_{\rm c}=T_{\rm
BC}(n_{\rm c})$. The heat capacity (21) in the considered mean-field model can
be calculated in terms of the ideal Bose gas quantities

$\displaystyle\left(\frac{\partial\varepsilon}{\partial T}\right)_{n=n_{\rm
c}}\cong c_{\rm id}(T=T_{\rm c},\mu^{*}=m)+c(T_{\rm c}-T),\leavevmode\nobreak\
$ (49)

that corresponds to the critical exponent $\alpha=-1$.

At $T<T_{\rm c}$ the chemical and mechanical equilibrium between the gaseous
phase with density $n_{\rm g}$ and the liquid phase with density $n_{\rm l}$
are given by the Gibbs conditions:

$\displaystyle\mu_{\rm g}=\mu_{\rm l}\leavevmode\nobreak\
,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ p(T,n_{\rm
g})=p(T,n_{\rm l})\leavevmode\nobreak\ .$ (50)

Equation (50) can be rewritten as

$\displaystyle m-\mu^{*}=U(n_{\rm g})-U(n_{\rm l})\leavevmode\nobreak\ ,$ (51)
$\displaystyle p_{\rm id}(T,\mu^{*})\leavevmode\nobreak\ -p_{\rm
id}(T,m)\leavevmode\nobreak\ =\leavevmode\nobreak\
\int_{n_{g}}^{n_{l}}dn\,n\,dU/dn\leavevmode\nobreak\ .$ (52)

Using

$\displaystyle p_{\rm id}(T,\mu^{*})\leavevmode\nobreak\
\cong\leavevmode\nobreak\ p_{\rm id}(T,m)\leavevmode\nobreak\
+\leavevmode\nobreak\ n_{\rm id}(T,m)\leavevmode\nobreak\
\left(\mu^{*}\leavevmode\nobreak\ -\leavevmode\nobreak\ m\right)$ (53)

one finds at small values of $\delta n_{\rm g}=n_{\rm c}-n_{\rm g}$ and
$\delta n_{\rm l}=n_{\rm l}-n_{\rm c}$

$\displaystyle n_{\rm l}-n_{\rm g}\equiv\delta n_{\rm g}+\delta n_{\rm
l}\sim(T_{c}-T)$ (54)

that means $\beta=1$. The details are presented in Appendix C. In particular,
an asymmetry between the gaseous and liquid binodals with $\delta n_{\rm
g}/\delta n_{\rm l}\geq 2$ observed. This is different from the properties of
the Boltzmann CP where $\delta n_{\rm g}/\delta n_{\rm l}=1$.

Figure 4: The scaled variance $\omega$ in the $(\mu/m,T)$ plane for $\alpha$
matter (a) and pion matter (b).

Equation (23) at the Bose CP can be presented as

$\displaystyle\left[\left(\frac{\partial p}{\partial
n}\right)_{T}\right]^{-1}_{n=n_{\rm c}}=\frac{\omega}{T}=\frac{\omega_{\rm
id}(T,\mu^{*}=m))}{T}\sim\leavevmode\nobreak\ |T_{c}-T|^{-1}$ (55)

that corresponds to $\gamma=1$, the same as in the Boltzmann CP.

In the vicinity of the CP at $T=T_{\rm c}$ and $n\rightarrow n_{\rm c}-0$, one
has $\mu^{*}\rightarrow m-0$ and finds (see also Ref. Begun and Gorenstein
(2008))

$\displaystyle p(T,\mu^{*})\cong p_{\rm id}(T,\mu^{*}=m)+n_{\rm
id}(T,\mu^{*}=m)\,\left(m-\mu^{*}\right)$
$\displaystyle+\int_{0}^{n}dn^{\prime}n^{\prime}\frac{dU}{dn^{\prime}}=p_{\rm
c}-\int_{n}^{n_{\rm c}}dn^{\prime}n^{\prime}\frac{dU}{dn^{\prime}}+n_{\rm
c}\,\left(m-\mu^{*}\right)$ $\displaystyle\cong
p_{c}+\left(B-\frac{2\pi^{2}}{T_{c}^{2}g^{2}m^{3}}\right)n_{c}(n-n_{c})^{2}\leavevmode\nobreak\
.$ (56)

In the non relativistic regime when $T_{\rm c}/m\ll 1$ the value of $T_{\rm
c}/T_{\rm c}^{0}$ is a function of $\sigma$ given by Eq. (77). This leads to

$\displaystyle p(T_{\rm c},n)\cong p_{\rm
c}+\left[1-\left(\frac{\sigma}{\sigma_{\rm cr}}\right)^{2/3}\right]Bn_{\rm
c}(n-n_{\rm c})^{2}\leavevmode\nobreak\ .$ (57)

Equation (57) corresponds to the Bose CP in the non-relativistic approximation
and, thus, it requires $\sigma>\sigma_{\rm cr}\cong 1.62$. In the vicinity of
the CP at $T=T_{\rm c}$ and $n\rightarrow n_{\rm c}+0$, one has $\mu^{*}=m$
and finds

$\displaystyle p(T_{\rm c},n)\cong p_{\rm c}+Bn_{\rm c}(n-n_{\rm
c})^{2}\leavevmode\nobreak\ .$ (58)

Equations (57) and (58) reveal the critical exponent $\delta=2$ at the Bose
CP. Thus, in contrast to the Boltzmann CP, the value of
$(\partial^{2}p/\partial n^{2})_{T=T_{\rm c}}$ at $n=n_{c}$ is not equal to
zero. Instead, this second derivative of the pressure has a discontinuity
equal to $A(\sigma/\sigma_{\rm cr})^{2/3}$ at the Bose CP.

To summarize, the set of the critical exponents at the Bose CP,

$\displaystyle\alpha\leavevmode\nobreak\ =\leavevmode\nobreak\
-1\leavevmode\nobreak\ ,\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \beta\leavevmode\nobreak\ =1\leavevmode\nobreak\
,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\gamma\leavevmode\nobreak\ =\leavevmode\nobreak\ 1\leavevmode\nobreak\
,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\delta=2\leavevmode\nobreak\ $ (59)

is different from that (25) at the Boltzmann CP. The set of critical exponents
(59) corresponds to a universality class of ordinary percolation in $6^{+}$
dimensions and directed percolation in $4^{+}$ dimensions (see Refs. Stauffer
and Aharony (2018); Saberi (2015); Stauffer (1979)). Each set of the critical
exponents (25) and (59) satisfies the scaling relations:

$\displaystyle
2-\alpha=2\beta+\gamma=\beta(\delta+1)=\gamma\frac{\delta+1}{\delta-1}\leavevmode\nobreak\
.$ (60)

In the special case of $\sigma=\sigma_{\rm cr}$ both sets of the critical
exponents (25) and (59) are simultaneously realized. The Boltzmann set (25)
corresponds to $n<n_{\rm c}=2n_{\rm c}^{0}$ and the Bose set (59) to $n>n_{\rm
c}=2n_{\rm c}^{0}$. This is a new feature of the critical exponents. For this
single value $\sigma=\sigma_{\rm cr}$ the CP keeps a memory of the Boltzmann
like behavior at $n<n_{\rm c}$ and the Bose one at $n>n_{\rm c}$.

##  IV double phase transitions

Our discussion up to now in most cases assumes the non relativistic
approximation $T_{\rm c}^{0}/m\ll 1$. If the relativistic effects become
important, another interesting possibility with the two CPs emerges: the
Boltzmann CP at $n_{\rm c}<2n_{\rm c}^{0}$ and the Bose CP at $n_{\rm
c}=2n_{\rm c}^{0}$. Both conditions (47) and (48) are then simultaneously
satisfied, and at the both CPs the particle number fluctuations are singular
with $\omega\rightarrow\infty$. For very large attractive interactions such a
possibility was observed at zero chemical potential in Ref. Stashko  _et al._
(2021b). The double phase transitions at $\mu>0$ will be discussed in this
section.

Let us consider first the small Bose effects in the ultrarelativistic limit
$m/T\ll 1$. This permits us to obtain simple analytical approximations. Taking
into account only the $k=1$ and $k=2$ terms in the power series (26) and (27)
one obtains in the $m/T\ll 1$ limit

$\displaystyle p_{\rm
id}\approx\frac{gT^{4}}{\pi^{2}}\left[X+\frac{X^{2}}{4}\right]\leavevmode\nobreak\
,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
n\approx\frac{gT^{3}}{\pi^{2}}\left[X+\frac{X^{2}}{2}\right]\leavevmode\nobreak\
,$ (61)

where $X=\exp[\mu^{*}/T]$. By inverting $n(X)$ to $X=X(n)$ and substituting it
into $p_{\rm id}(X)$ one finds

$\displaystyle p_{\rm id}\leavevmode\nobreak\ \approx\leavevmode\nobreak\
Tn\leavevmode\nobreak\ -\leavevmode\nobreak\
\frac{3\pi^{2}}{4gT^{2}}\leavevmode\nobreak\ n^{2}\leavevmode\nobreak\ ,$ (62)

and then

$\displaystyle p(T,n)\approx\leavevmode\nobreak\
Tn-\left(\frac{A}{2}+\frac{\Gamma}{T^{2}}\right)n^{2}+\frac{2B}{3}n^{3}\leavevmode\nobreak\
.$ (63)

where $\Gamma=3\pi^{2}/(4g).$ This is the ultrarelativistic analog of non
relativistic approximation (33). From Eq. (63) one finds condition for the CP
parameters

$\displaystyle T-(A+2\Gamma T^{-2})n+2Bn^{2}$ $\displaystyle=0,$ (64)
$\displaystyle-(A+2\Gamma T^{-2})+4Bn$ $\displaystyle=0.$ (65)

From Eqs. (64) and (65) one obtains approximate relations

$\displaystyle n_{c}^{1}\approx n_{c}^{0}(1+\kappa),\leavevmode\nobreak\
\leavevmode\nobreak\ T_{c}^{1}\approx
T_{c}^{0}(1+2\kappa),\leavevmode\nobreak\ \leavevmode\nobreak\ $ (66)

valid at

$\displaystyle\kappa=\frac{96\pi^{2}\tilde{B}^{2}}{\tilde{A}^{5}}\ll
1\leavevmode\nobreak\ ,\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \frac{m}{T_{\rm c}^{0}}\ll
1\leavevmode\nobreak\ .$ (67)

With increasing $\kappa$ the critical density $n_{\rm c}$ reaches its limiting
value of $n_{\rm c}^{\rm lim}=2n_{\rm c}^{0}$, and the critical temperature
behaves as

$\displaystyle T_{\rm c}=T_{\rm BC}\leavevmode\nobreak\
\approx\left[\frac{\pi^{2}\leavevmode\nobreak\ n_{\rm c}^{\rm
lim}}{g\zeta(3)}\right]^{1/3}\leavevmode\nobreak\ .$ (68)

Let us denote by $n_{\rm tr}$ the particle number density at the triple point,
i.e., the point of an intersection of the BEC line with a liquid binodal. This
point is shown by an open circle in Fig. 2 (a). The values of $n_{\rm c}$ and
$n_{\rm tr}$ satisfy the following inequalities:

$\displaystyle n_{\rm c}^{0}\leavevmode\nobreak\ <\leavevmode\nobreak\ n_{\rm
c}\leavevmode\nobreak\ \leq\leavevmode\nobreak\ n^{\rm lim}=2n_{\rm
c}^{0}\leavevmode\nobreak\ \leq\leavevmode\nobreak\ n_{\rm
tr}\leavevmode\nobreak\ <\leavevmode\nobreak\ 3n_{\rm
c}^{0}\leavevmode\nobreak\ .$ (69)

Figure 5:  The relativistic $(\tilde{A},\tilde{B})$ diagram in the mean-field
model. The region of the double CPs is shown by yellow color between the Bose
and Boltzmann regions. See text for details. Figure 6:  (a) An example of the
phase transitions with two CPs in the $(n,T)$ plane. It corresponds to
$\tilde{A}=12.5$ and $\tilde{B}=60$. The Boltzmann CP corresponds to
$\tilde{A}=8$ and $\tilde{B}=\ 1.3$, and the Bose CP to $\tilde{A}=9.7$ and
$\tilde{B}=59$. (b) The $(\mu,T)$ plane with the double phase transition.
Yellow and blue stars correspond to the Boltzmann CP and Bose CP,
respectively. A horizontal dotted line illustrates two serial transitions with
increasing $\mu$ at fixed $T$.

In the non-relativistic approximation a transition to the Bose CP takes place
as a merger of the Boltzmann CP and the triple point at $n=2n_{\rm c}^{0}$
that happens at $\sigma=\sigma_{\rm cr}\cong 1.62$ (see Appendix B). In a
presence of the relativistic effects a transition of the Boltzmann CP to the
Bose CP takes place through a double phase transition. This happens in a
region of the parameters $\tilde{A}$ and $\tilde{B}$ noted as ”double CP” in
Fig. 5 and shown with yellow color. Several examples of the simultaneous
presence of the two CPs are shown in Fig. 6. A transition from the Boltzmann
CP to the Bose CP has the following general features seen in Fig. 6. First,
the triple point reaches the limiting density $2n_{\rm c}^{0}$. Thus, both the
Boltzmann CP at $n=n_{\rm c}<2n_{\rm c}^{0}$ and the Bose CP at $n_{\rm
c}=2n_{\rm c}^{0}$ simultaneously exist. Then the Boltzmann CP continues to
move to larger $n_{\rm c}$ values, and finally disappears at $n_{\rm
c}<2n^{0}_{\rm c}$.

A big variety of the mutual arrangements of the both CP temperatures are found
at different sets of the $\tilde{A}$ and $\tilde{B}$ values. As seen from Fig.
5, this $(\tilde{A},\tilde{B})$ region is not far away from the $T_{\rm
c}^{0}=m$ line, i.e., a presence of the double CPs requires moderate
relativistic effects, and it is absent in the both non-relativistic and ultra-
relativistic limits.

The double phase transitions require a ”quadruple point” where the two phase
transition lines in the $(\mu,T)$ plane intersect. This point at $T=T_{qd}$ is
shown as a full circle in Fig. 6(b) where an example of a double phase
transitions is presented in the $(\mu,T)$ phase diagram. The horizontal dotted
line in Fig. 6(b) corresponds a fixed temperature $T$ that is larger than
$T_{\rm qd}$ and smaller than both $T_{\rm c}({\rm Boltz})$ and $T_{\rm
c}({\rm Bose})$. With increasing of $\mu$ at fixed $T$ the system moves right
along the dotted line in Fig. 6(b). Two consecutive phase transitions take
place. The first one is between the gaseous and liquid phases, where both
phases do not include the BC. The second transition takes place between the
liquid phase without the BC and the liquid with the non zero BC density
$n_{\rm BC}>0$.

Note that a contribution of antiparticles to the thermodynamic functions is
small in comparison with that from particles at $\mu/T>1$. In the Boltzmann
approximation one finds for the CP values

$\displaystyle\frac{\mu_{\rm c}^{0}}{T_{\rm c}^{0}}\leavevmode\nobreak\
=\leavevmode\nobreak\ \ln\left[4\pi^{2}\tilde{A}^{-1}\leavevmode\nobreak\
K_{2}^{-1}(2\tilde{B}/\tilde{A}^{2})\right]\leavevmode\nobreak\
-\leavevmode\nobreak\ 1\leavevmode\nobreak\ .$ (70)

The line $\mu_{\rm c}^{0}/T_{\rm c}^{0}=1$ is shown in the
$(\tilde{A},\tilde{B})$ diagram in Fig. 5 by a dashed line. From Fig. 5 one
observes that both relations, $T_{\rm c}^{0}\approx m$ and $\mu_{\rm
c}^{0}\approx m$, are approximately satisfied in the region of the
$(\tilde{A},\tilde{B})$ parameters with the double CP.

##  V Summary

Thermodynamic properties of interacting Bose particles are studied with a
relativistic mean-field model. The mean-filed potential is chosen to be
$U(n)=-An+Bn^{2}$, where $A>0$ and $B>0$ terms correspond to attractive and
repulsion interactions, respectively. The whole range of the system parameters
$(A,B,m,g)$ is considered. The CP temperature $T_{\rm c}^{0}$ and particle
number density $n_{\rm c}^{0}$ in the Boltzmann approximation are given by the
functions (19) of the interaction parameters $A$ and $B$, and they are
independent of the values of $m$ and $g$. The Bose effects are sensitive to
the values of $m$ and $g$, and they lead to increase of the $T_{\rm c}$ and
$n_{\rm c}$ values in comparison to their $T_{\rm c}^{0}$ and $n_{\rm c}^{0}$
Boltzmann limits.

We found two different implementations of the CP, denoted as the Boltzmann CP
and Bose CP, in the system of interacting bosons. These two scenarios
correspond to small and large Bose effects, respectively. The diagram of the
dimensionless parameters $(\tilde{A},\tilde{B})$ gives us a possibility to
study a diverse physical systems of Bose particles. The interacting $\alpha$
particles and pions are presented as illustrative examples of the Botzmann CP
and Bose CP, respectively. The Bose CP has the limiting density $n_{\rm
c}=2n_{\rm c}^{0}$ and the BEC line touches the Bose CP. The Boltzmann CP and
Bose CP have rather different properties. They have two different sets of the
critical exponents and, thus, belong to distinct universality classes.

On the $(\tilde{A},\tilde{B})$ diagram a special region of the system
parameters was found. In this region $T_{\rm c}^{0}\approx m$ and $\mu_{\rm
c}^{0}\approx m$ corresponding to the double phase transitions with two CPs –
the Boltzmann CP at $n_{\rm c}<2n_{\rm c}^{0}$ and the Bose CP at $n_{\rm
c}=2n_{\rm c}^{0}$.

###### Acknowledgements.

The authors thank M. Gazdzicki, R. Poberezhnyuk, L. Satarov, H. Stoecker, and
V. Vovchenko for fruitful comments and discussions. This work is supported by
the National Academy of Sciences of Ukraine, Grant No. 0121U112254. O.S.S.
acknowledges support from National Research Foundation of Ukraine (Project No.
2020.02/0073). M.I.G. acknowledges support from the Alexander von Humboldt
Foundation.

##  Appendix A Dimensional analysis

Introducing dimensionless variables

$\displaystyle\hat{T}=\frac{T}{m},\leavevmode\nobreak\ \leavevmode\nobreak\
\hat{\mu}^{*}=\frac{\mu^{*}}{m}\leavevmode\nobreak\ ,\leavevmode\nobreak\
\leavevmode\nobreak\ \hat{n}=\frac{n}{gm^{3}},\leavevmode\nobreak\
\leavevmode\nobreak\ \hat{p}=\frac{p}{gm^{4}}$ (71)

one can rewrite Eqs. (26) and (27) as the following functions

$\displaystyle\hat{p}=\hat{p}(\hat{T},\hat{\mu}^{*};\leavevmode\nobreak\
\tilde{A},\tilde{B})\leavevmode\nobreak\ ,\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\
\hat{n}=\hat{n}(\hat{T},\hat{\mu}^{*};\leavevmode\nobreak\
\tilde{A},\tilde{B})\leavevmode\nobreak\ ,$ (72)

where $\tilde{A}$ and $\tilde{B}$ are given by Eq. (39). Presenting
$\hat{\mu}^{*}$ as $\hat{\mu}^{*}(\hat{T},\hat{n};\leavevmode\nobreak\
\tilde{A},\tilde{B})$ one obtains

$\displaystyle\hat{p}=\hat{p}(\hat{T},\hat{n};\leavevmode\nobreak\
\tilde{A},\tilde{B})\leavevmode\nobreak\ .$ (73)

At fixed temperature $\hat{T}>\hat{T}_{\rm c}$ the system isotherm is a
monotonic function of $\hat{n}$, i.e.,
$(\partial\hat{p}/\partial\hat{n})_{\hat{T}}>0$. When $\hat{T}$ decreases, one
comes to the critical isotherm $\hat{T}_{\rm c}$ for which
$(\partial\hat{p}/\partial\hat{n})_{\hat{T}_{\rm c}}=0$ at the critical
density $\hat{n}_{\rm c}$. From these critical quantities one can then find
the values of $T_{\rm c}/T_{\rm c}^{0}$ and $n_{\rm c}/n_{\rm c}^{0}$ as
functions of the two parameters, namely, the functions of $\tilde{A}$ and
$\tilde{B}$. Therefore, the type of the CP, either Boltzmann CP or Bose CP, is
defined by the values of the dimensionless parameters $\tilde{A}$ and
$\tilde{B}$. As shown in Appendix B, in the non-relativistic approximation
this is reduced to a dependence on the single parameter
$\sigma=\sigma(\tilde{A},\tilde{B})$.

##  Appendix B Non relativistic approximation

The spinodal is defined by the condition $(\partial p/\partial n)_{T}=0$, that
can be written as

$\displaystyle T+n\frac{dU}{dn}\omega_{\rm
id}(T,\mu^{*})=0\leavevmode\nobreak\ .$ (74)

The functions $T_{\rm s}=T_{\rm s}(n)$ given by Eq. (74) are shown in Fig. 7.
These spinodal lines intersect with the BEC lines $T_{\rm BC}(n)$ at $n=n^{\rm
lim}=2n_{\rm c}^{0}$ as shown in Fig. 7. The analytical approximations valid
near the BEC line can be used for $\omega_{\rm id}$ in Eq. (74) (see Refs.
Begun and Gorenstein (2008, 2007)):

$\displaystyle\omega_{\rm id}(T,\mu^{*})$
$\displaystyle=\frac{gm^{3/2}T^{2}}{2\sqrt{2}\pi
n}(m-\mu^{*})^{-1/2}\leavevmode\nobreak\ ,$ (75) $\displaystyle m-\mu^{*}(T)$
$\displaystyle=\frac{9\zeta^{2}(3/2)}{16\pi}\frac{(T-T_{\rm BC})^{2}}{T_{\rm
BC}}\leavevmode\nobreak\ .$ (76)

The BEC temperature $T_{\rm BC}$ in the non-relativistic approximation (29)
can be written as

$\displaystyle\tilde{T}_{\rm BC}=\frac{T_{\rm BC}}{T^{0}_{\rm
c}}=2\left[\frac{2\tilde{n}\sigma}{\zeta(3/2)}\right]^{2/3}\leavevmode\nobreak\
,$ (77)

where $\sigma$ is given by Eq. (38).

Figure 7: The spinodal lines $T_{\rm s}=T_{\rm s}(n)$ for different values of
the parameter $\sigma$. The liquid spinodal becomes the vertical line at
$n=2n^{0}_{\rm c}$.

By substituting Eqs. (75-77) into Eq. (74) and rewriting it in the reduced
variables, one finds the function $\tilde{T}_{\rm s}(\tilde{n})$ at
$\tilde{n}\rightarrow 2-0$. The Boltzmann CP corresponds to $d\tilde{T}_{\rm
s}/d\tilde{n}<0$ and the Bose CP to $d\tilde{T}_{\rm s}/d\tilde{n}>0$ at
$\tilde{n}=2$ (see Fig. 7). Thus, the condition

$\displaystyle\left.\frac{d\tilde{T}_{\rm
s}}{d\tilde{n}}\right|_{\tilde{n}=2}\leavevmode\nobreak\ =\leavevmode\nobreak\
0$ (78)

corresponds to the transition between the Boltzmann CP and the Bose CP.

The solution of Eq. (78) gives the critical value for the parameter
$\sigma=\sigma_{\rm cr}$. The straightforward analytical solution of Eq. (78)
gives

$\displaystyle\sigma_{\rm cr}\leavevmode\nobreak\ =\leavevmode\nobreak\
\frac{2\pi^{3/2}}{\zeta^{2}(3/2)}\leavevmode\nobreak\
\cong\leavevmode\nobreak\ 1.62\leavevmode\nobreak\ .$ (79)

##  Appendix C Binodals near the Bose CP

In a case of the Bose CP we introduce

$\displaystyle\frac{n_{\rm c}-n_{\rm g}}{n_{\rm c}}$
$\displaystyle\equiv\frac{\delta n_{\rm g}}{n_{\rm c}}=C_{\rm
g}\leavevmode\nobreak\ \frac{T_{\rm c}-T}{T_{\rm c}}\leavevmode\nobreak\ ,$
(80) $\displaystyle\frac{n_{\rm l}-n_{\rm c}}{n_{\rm c}}$
$\displaystyle\equiv\frac{\delta n_{\rm l}}{n_{\rm c}}=C_{\rm
l}\leavevmode\nobreak\ \frac{T_{\rm c}-T}{T_{\rm c}}\leavevmode\nobreak\ .$
(81)

At $T<T_{\rm c}-0$ one finds $\mu_{\rm g}^{*}(T)<m$ on the gaseous binodal and
$\mu_{\rm l}^{*}=m$ on the liquid one.

At $T\rightarrow T_{\rm c}$ one has $\mu_{\rm g}^{*}\rightarrow m$ and, using
Eq. (76), obtains the following system for dimensional binodal slope
parameters $C_{\rm g}$ and $C_{\rm l}$:

$\displaystyle\begin{split}&2\left(C_{g}^{3}+C_{l}^{3}\right)=3\left(C_{g}^{2}-C_{l}^{2}\right)\frac{T_{\rm
c}}{n_{\rm c}}\left.\frac{dn_{\rm id}(T,m)}{dT}\right|_{T=T_{\rm c}}\\\
&+2R(C_{g}^{2}-C_{l}^{2})^{3/2}\leavevmode\nobreak\ ,\end{split}$ (82)
$\displaystyle C_{g}=\frac{T_{\rm c}}{n_{\rm c}}\left.\frac{dn_{\rm
id}(T,m)}{dT}\right|_{T=T_{\rm
c}}+R(C_{g}^{2}-C_{l}^{2})^{1/2}\leavevmode\nobreak\ ,$ (83)

where $R=\sqrt{B}gT_{c}m^{3/2}/(4\sqrt{2}\pi)$.

In the non relativistic approximation the derivative will be

$\displaystyle\left.\frac{dn_{\rm id}(T,m)}{dT}\right|_{T=T_{\rm
c}}=\frac{3n_{\rm c}}{2T_{\rm c}}\leavevmode\nobreak\ .$ (84)

Using Eqs. (76), (82) and (83) one finds

$\displaystyle\begin{split}&(C_{g}^{3}+C_{l}^{3})=\frac{9}{4}(C_{g}^{2}-C_{l}^{2})\\\
&+\left(\frac{\sigma_{cr}}{\sigma}\right)^{1/3}(C_{g}^{2}-C_{l}^{2})^{3/2}\leavevmode\nobreak\
,\end{split}$ (85) $\displaystyle
C_{g}=\frac{3}{2}+\left(\frac{\sigma_{cr}}{\sigma}\right)^{1/3}\sqrt{C_{g}^{2}-C_{l}^{2}}\leavevmode\nobreak\
.$ (86)

From Eqs. (85) and (86) at $\sigma\to\sigma_{\rm cr}+0$

$\displaystyle C_{g}\to\infty,\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ C_{l}\to\infty,\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\frac{C_{g}}{C_{l}}\to\infty\leavevmode\nobreak\ ,$ (87)

and

$\displaystyle C_{g}=\frac{3}{2},\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\
C_{l}=\frac{3}{4},\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\
\frac{C_{g}}{C_{l}}=2\leavevmode\nobreak\ ,$ (88)

at $\sigma\gg\sigma_{\rm cr}$. Thus, there is an asymmetry between the gaseous
and liquid binodals. This asymmetry becomes very large when $\sigma$
approaches its critical value $\sigma_{\rm cr}$.

At $\sigma=\sigma_{\rm cr}$ Eq. (80) is broken and should be substituted by
$\delta n_{\rm g}/n_{\rm c}=C_{\rm g}\sqrt{(T_{\rm c}-T)/T_{\rm c}}$. It means
that for $\sigma=\sigma_{\rm cr}$ the critical exponent $\beta$ has the
Boltzmann value $\beta=1/2$ at $n<n_{\rm c}=2n_{\rm c}^{0}$ and the Bose value
$\beta=1$ at $n>n_{\rm c}=2n_{\rm c}^{0}$

With increasing $T_{c}$ value in Eqs. (82) and (86) one obtains a limiting
behavior

$\displaystyle C_{\rm g}\cong\frac{T_{\rm c}}{n_{\rm c}}\,\left.\frac{dn_{\rm
id}(T,m)}{dT}\right|_{T=T_{\rm c}}\leavevmode\nobreak\ ,\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ C_{\rm
l}\leavevmode\nobreak\ =\leavevmode\nobreak\ \frac{1}{2}\,C_{\rm
g}\leavevmode\nobreak\ $ (89)

that corresponds to the large $\sigma$ behavior (88) in the non relativistic
approximation.

In the ultra relativistic limit one has

$\displaystyle\left.\frac{dn_{\rm id}(T,m)}{dT}\right|_{T=T_{\rm
c}}=\frac{3n_{\rm c}}{T_{\rm c}}\leavevmode\nobreak\ $ (90)

instead of Eq. (84). Using Eq. (89) in the ultra-relativistic limit,

$\displaystyle C_{g}=3,\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\
C_{l}=\frac{3}{2},\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\
\frac{C_{g}}{C_{l}}=2\leavevmode\nobreak\ .$ (91)

## References

  * Landau and Lifshitz (1975) L. D. Landau and E. M. Lifshitz,  _Statistical Physics_ (Pergamon, Oxford, 1975).
  * Greiner  _et al._ (2012) W. Greiner, L. Neise, and H. Stöcker,  _Thermodynamics and statistical mechanics_ (Springer Science & Business Media, 2012).
  * Bose (1924) S. N. Bose, Z. Phys. 26, 178 (1924).
  * Einstein (1925) A. Einstein, Kgl. Preuss. Akad. Wiss 1, 137 (1925).
  * Anderson  _et al._ (1995) M. H. Anderson, J. R. Ensher, M. R. Matthews, C. E. Wieman, and E. A. Cornell, Science 269, 198 (1995).
  * Bradley  _et al._ (1995) C. C. Bradley, C. A. Sackett, J. J. Tollett, and R. G. Hulet, Phys. Rev. Lett. 75, 1687 (1995).
  * Davis  _et al._ (1995) K. B. Davis, M. O. Mewes, M. R. Andrews, N. J. van Druten, D. S. Durfee, D. M. Kurn, and W. Ketterle, Phys. Rev. Lett. 75, 3969 (1995).
  * Dalfovo  _et al._ (1999) F. Dalfovo, S. Giorgini, L. P. Pitaevskii, and S. Stringari, Rev. Mod. Phys. 71, 463 (1999).
  * Kapusta and Gale (2006) J. I. Kapusta and C. Gale,  _Finite-Temperature Field Theory: Principles and Applications_, 2nd ed., Cambridge Monographs on Mathematical Physics (Cambridge University Press, 2006).
  * Andersen (2004) J. Andersen, Reviews of Modern Physics 76, 599–639 (2004).
  * Griffin  _et al._ (1996) A. Griffin, D. W. Snoke, and S. Stringari,  _Bose-einstein condensation_ (Cambridge University Press, 1996).
  * Watabe and Ohashi (2013) S. Watabe and Y. Ohashi, Phys. Rev. A 88, 053633 (2013).
  * Watabe (2019) S. Watabe, Acta Phys. Polon. A 135, 1222–1230 (2019).
  * Baym  _et al._ (1999a) G. Baym, J.-P. Blaizot, M. Holzmann, F. Laloë, and D. Vautherin, Phys. Rev. Lett. 83, 1703 (1999a).
  * Baym  _et al._ (2000) G. Baym, J.-P. Blaizot, and J. Zinn-Justin, EPL 49, 150 (2000).
  * Holzmann and Krauth (1999) M. Holzmann and W. Krauth, Phys. Rev. Lett. 83, 2687 (1999), arXiv:cond-mat/9905198 [cond-mat.stat-mech] .
  * Holzmann  _et al._ (2001) M. Holzmann, G. Baym, J.-P. Blaizot, and F. Laloe, Phys. Rev. Lett. 87, 120403 (2001), arXiv:cond-mat/0103595 [cond-mat.stat-mech] .
  * Savchuk  _et al._ (2020a) O. Savchuk, Y. Bondar, O. Stashko, R. V. Poberezhnyuk, V. Vovchenko, M. I. Gorenstein, and H. Stoecker, Phys. Rev. C 102, 035202 (2020a).
  * Baym  _et al._ (1999b) G. Baym, J.-P. Blaizot, M. Holzmann, F. Laloe, and D. Vautherin, Phys. Rev. Lett. 83, 1703 (1999b), arXiv:cond-mat/9905430 [cond-mat] .
  * Huang (1999) K. Huang, Phys. Rev. Lett. 83, 3770 (1999).
  * Schakel (2001) A. M. J. Schakel, Phys. Rev. E 63, 026115 (2001).
  * Feynman (1953a) R. P. Feynman, Phys. Rev. 91, 1291 (1953a).
  * Feynman (1953b) R. P. Feynman, Phys. Rev. 91, 1301 (1953b).
  * Satarov  _et al._ (2017) L. M. Satarov, M. I. Gorenstein, A. Motornenko, V. Vovchenko, I. N. Mishustin, and H. Stoecker, Journal of Physics G: Nuclear and Particle Physics 44, 125102 (2017).
  * Begun and Gorenstein (2007) V. V. Begun and M. I. Gorenstein, Phys. Lett. B653, 190 (2007), arXiv:hep-ph/0611043 [hep-ph] .
  * Begun and Gorenstein (2008) V. V. Begun and M. I. Gorenstein, Phys. Rev. C77, 064903 (2008), arXiv:0802.3349 [hep-ph] .
  * Strinati  _et al._ (2018) G. C. Strinati, P. Pieri, G. Röpke, P. Schuck, and M. Urban, Physics Reports 738, 1–76 (2018).
  * Nozieres and Schmitt-Rink (1985) P. Nozieres and S. Schmitt-Rink, J. Low. Temp. Phys. 59, 195 (1985).
  * Funaki  _et al._ (2008) Y. Funaki, T. Yamada, H. Horiuchi, G. Röpke, P. Schuck, and A. Tohsaki, Phys. Rev. Lett. 101, 082502 (2008).
  * Chavanis and Harko (2012) P.-H. Chavanis and T. Harko, Phys. Rev. D86, 064011 (2012), arXiv:1108.3986 [astro-ph.SR] .
  * Mishustin  _et al._ (2019) I. Mishustin, D. Anchishkin, L. Satarov, O. Stashko, and H. Stoecker, Phys. Rev. C 100, 022201 (2019), arXiv:1905.09567 [nucl-th] .
  * Padilla  _et al._ (2019) L. E. Padilla, J. A. Vázquez, T. Matos, and G. Germán, Journal of Cosmology and Astroparticle Physics 2019, 056 (2019).
  * Harko and Madarassy (2022) T. Harko and E. J. Madarassy, Eur. Phys. J. C 82, 401 (2022), arXiv:2205.00297 [gr-qc] .
  * Vovchenko (2017) V. Vovchenko, Phys. Rev. C 96, 015206 (2017), arXiv:1701.06524 [nucl-th] .
  * Vovchenko (2020) V. Vovchenko, Int. J. Mod. Phys. E 29, 2040002 (2020), arXiv:2004.06331 [nucl-th] .
  * Satarov  _et al._ (2021) L. M. Satarov, R. V. Poberezhnyuk, I. N. Mishustin, and H. Stoecker, Phys. Rev. C 103, 024301 (2021), arXiv:2009.13487 [nucl-th] .
  * Poberezhnyuk  _et al._ (2021) R. V. Poberezhnyuk, O. Savchuk, M. I. Gorenstein, V. Vovchenko, and H. Stoecker, Phys. Rev. C 103, 024912 (2021), arXiv:2011.06420 [hep-ph] .
  * Stashko  _et al._ (2021a) O. S. Stashko, O. V. Savchuk, R. V. Poberezhnyuk, V. Vovchenko, and M. I. Gorenstein, Phys. Rev. C 103, 065201 (2021a).
  * Kuznietsov  _et al._ (2021) V. A. Kuznietsov, O. S. Stashko, O. V. Savchuk, and M. I. Gorenstein, Phys. Rev. C 104, 055202 (2021).
  * Fedotkin  _et al._ (2019) S. N. Fedotkin, A. G. Magner, and M. I. Gorenstein, Phys. Rev. C 100, 054334 (2019).
  * Fedotkin  _et al._ (2021) S. N. Fedotkin, A. G. Magner, and U. V. Grygoriev, (2021), arXiv:2112.13704 [nucl-th] .
  * Sanzhur (2022) A. I. Sanzhur, (2022), arXiv:2202.09749 [nucl-th] .
  * Anchishkin and Vovchenko (2015) D. Anchishkin and V. Vovchenko, Journal of Physics G: Nuclear and Particle Physics 42, 105102 (2015).
  * Vovchenko  _et al._ (2020) V. Vovchenko, C. Greiner, V. Koch, and H. Stoecker, Phys. Rev. D 101, 014015 (2020).
  * Savchuk  _et al._ (2020b) O. Savchuk, V. Vovchenko, R. V. Poberezhnyuk, M. I. Gorenstein, and H. Stoecker, Phys. Rev. C 101, 035205 (2020b).
  * Kos  _et al._ (2016) F. Kos, D. Poland, D. Simmons-Duffin, and A. Vichi, Journal of High Energy Physics 2016 (2016), 10.1007/jhep08(2016)036.
  * Röpke (2009) G. Röpke, Phys. Rev. C 79 (2009), 10.1103/physrevc.79.014002.
  * Stein  _et al._ (1998) H. Stein, C. Porthun, and G. Röpke, The European Physical Journal B 2, 393 (1998).
  * Stauffer and Aharony (2018) D. Stauffer and A. Aharony,  _Introduction to percolation theory_ (Taylor & Francis, 2018).
  * Saberi (2015) A. A. Saberi, Physics Reports 578, 1 (2015).
  * Stauffer (1979) D. Stauffer, Physics Reports 54, 1 (1979).
  * Stashko  _et al._ (2021b) O. S. Stashko, D. V. Anchishkin, O. V. Savchuk, and M. I. Gorenstein, Journal of Physics G: Nuclear and Particle Physics 48, 055106 (2021b).