Bose–Einstein Condensate Dark Matter That Involves Composites

A.M. Gavrilik111e-mail<EMAIL_ADDRESS>and A.V. Nazarenko

Bogolyubov Institute for Theoretical Physics of NAS of Ukraine
14b, Metrolohichna Str., Kyiv 03143, Ukraine

###### Abstract

By improving the Bose–Einstein condensate model of dark matter through the
repulsive three-particle interaction to better reproduce observables such as
rotation curves, both different thermodynamic phases and few-particle
correlations are revealed. Using the numerically found solutions of the
Gross–Pitaevskii equation for averaging the products of local densities and
for calculating thermodynamic functions at zero temperature, it is shown that
the few-particle correlations imply a first-order phase transition and are
reduced to the product of single-particle averages with a simultaneous
increase in pressure, density, and quantum fluctuations. Under given
conditions, dark matter exhibits rather the properties of an ideal gas with an
effective temperature determined by quantum fluctuations. Characteristics of
oscillations between bound and unbound states of three particles are estimated
within a simple random walk approach to qualitatively models the instability
of particle complexes. On the other hand, the density-dependent conditions for
the formation of composites are analyzed using chemical kinetics without
specifying the bonds formed. The obtain results can be extended to the models
of multicomponent dark matter consisting of composites formed by particles
with a large scattering length.

Keywords: dark matter; halo core; Bose–Einstein condensate; two- and three-
particle self-interactions; two-phase structure; first-order phase transition;
composites; random walk; chemical kinetics

##  1 Introduction

The current observable restrictions on the properties of dark matter still
admit the construction of theoretical models with different constituents [1,
2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22,
23, 24, 25, 26]. Adopting the concept of dark matter consisting of light
bosons in a condensate state [15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26],
we might also borrow the known properties of atomic Bose–Einstein condensates
in gaseous or liquid form [27, 28, 29], reproduced in the laboratory.

In the case of weakly interacting structureless bosons with or without
specification of their sort, the description of the macroscopic properties of
the dark matter halo of galaxies is usually based on the Gross–Pitaevskii
equation under the imposed conditions of diluteness and a large healing
length. Indeed, such models for the halo of dwarf galaxies give a good
qualitative description with a very short (pairwise) scattering length of
ultralight bosons [5, 20, 24]. It turns out that attempts to achieve better
agreement between theoretical estimates and astronomical observables by taking
into account additional terms of few-particle interactions in the
Gross–Pitaevskii equation result in the existence of various thermodynamic
phases of dark matter [5, 23, 24], as well as nontrivial correlations
associated with few-particle complexes (composites).

On the other hand, there are experimental data on the formation of Bose
condensate droplets [30, 31, 32] from a gas of alkali atoms with a long
scattering length, which is accompanied by the processes of particle
recombination, three-particle losses, and the formation of molecules that are
technically suppressed in the laboratory. This regime, associated with the
unitary limit of infinite scattering length, implies a strong interaction
between particles leading to the structural transformations. Much attention is
paid to its study in the literature [33, 34, 35] in connection with
experimental exploration of the Bose condensate. Theoretically, similar
processes like noted above cannot in principal be ruled out in dark matter,
which could eventually be a (multicomponent) mixture of different composites
with a resulting short scattering length, as required by current observations.
In addition, they could take place under the extreme conditions in the early
Universe [36]. Thus, the assumptions proposed give rise to many implications,
and we would like to dwell on a few. It is worth to note that dark matter,
made of several constituents, is already being actively studied [37, 38].

Summing up, we outline the tasks of this paper, which is devoted to the
analysis of density-dependent conditions for the formation of composites using
rather simple models.

In Section 2, considering a model based on the Gross–Pitaevskii equation that
takes into account two- and three-particle interactions, we try to obtain the
admissibility conditions for the occurrence of particle complexes by
estimating few-particle correlators over a relatively wide range of density,
pressure, and quantum fluctuations at zero temperature, i.e., quantities that
affect the transfer momentum in an inhomogeneous one-component system of
bosons. In Section 3, we describe the time-oscillatory instability of two- and
three-particle states that arises in a dense medium using a simple random walk
model [39]. Various scenarios for the formation of composites from three
particles are considered in Section 4 within the framework of chemical
kinetics [40] without clarifying the nature of the bonds. We operate only with
the scattering length and the density of the components as the most universal
characteristics, in order to extract the necessary conditions for the
existence of those. We end with the Discussion where we present our
conclusions on the obtained results.

##  2 Self-Gravitating Bose–Einstein Condensate: The Role of Few-Particle
Interactions

We start with an analysis of the equilibrium properties of halo dark matter,
described by the Gross–Pitaevskii equation, where both two- and three-particle
repulsive interactions are taken into account. As established in [24], such a
model reveals two-phase structure, represented by dilute and denser
(liquidlike) states, and the first-order phase transition. It is intended to
describe the resulting state of dark matter after a long-term evolution under
certain conditions, although only some of the physical options provided by
this model can be justified. After many attempts of description, with their
special features, made in the literature [20, 5, 23], here we can expect and
generally conclude that the resulting constituents of dark matter possess a
vanishingly small (pair) $s$-scattering length $a$. This does not contradict
the development of models describing particles with a longer scattering
length, because their composites (with a short scattering length $a$) are
allowed to exist in this “final” stage.

So, we first refer to our exploration of the model elaborated in [24], where
the Bose–Einstein condensate (BEC) is described by a real function $\psi(r)$
of the radial variable $r=|\mathbf{r}|$ in a ball
$B=\\{\mathbf{r}\in\mathbb{R}^{3}|\,|\mathbf{r}|\leq R\\}$. The main equation
of our model is generated by the functional $\Gamma$, depended on a chemical
potential $\mu$, and is supplemented by the Poisson equation for gravitational
field $4\pi GmV(r)$:

$\displaystyle\Gamma$ $\displaystyle=$ $\displaystyle
4\pi\int_{0}^{R}\left[\frac{\hbar^{2}}{2m}(\partial_{r}\psi(r))^{2}-\mu\psi^{2}(r)+4\pi
Gm^{2}\psi^{2}(r)V(r)\right.$ (1)
$\displaystyle\left.+\frac{U_{2}}{2}\psi^{4}(r)+\frac{U_{3}}{3}\psi^{6}(r)\right]\,r^{2}\,\mathrm{d}r,$
$\displaystyle\Delta_{r}V(r)$ $\displaystyle=$ $\displaystyle\psi^{2}(r).$ (2)

Let two- and three-particle interactions be characterized by the quantities:

$U_{2}=\frac{4\pi\hbar^{2}a}{m},\quad
U_{3}=gU_{2}\lambda^{3},\quad\lambda=\frac{\hbar}{mc},$ (3)

where $a>0$ is the $s$-scattering length of pairwise interaction; $\lambda$ is
the Compton wave length of the particle with mass $m$; $g$ is a some
dimensionless parameter introduced for convenience; $\Delta_{r}$ and
$\Delta^{-1}_{r}$ are the radial part of the Laplace operator and its inverse.

Note that the relativistic nature of three-particle interaction, because of
dependence of $\lambda$ and $U_{3}$ on speed of light $c$, is pointed out in
[5].

Applicability of the local form of interaction requires that

$n_{0}a^{3}<1,\quad n_{0}\lambda^{3}<1,$ (4)

where $n_{0}$ is the maximal particle density in the system.

Vanishing of variation $\delta\Gamma$ with respect to $\delta\psi$ gives the
stationary Gross–Pitaevskii equation with the gravitational and three-particle
interactions:

$\left[-\frac{\hbar^{2}}{2m}\Delta_{r}-\mu+4\pi
Gm^{2}V(r)\right]\psi+U_{2}\psi^{3}+U_{3}\psi^{5}=0,\quad\psi(R)=0,$ (5)

where $R$ is defined as the (first) zero of $\psi(r)$ and determines a core
boundary of BEC.

The further description of the system is conveniently carried out in terms of
dimensionless wave function $\chi$ and coordinate $\xi$:

$\psi(r)=\sqrt{n_{0}}\,\chi(\xi),\quad r=r_{0}\,\xi,\quad R=r_{0}\,\xi_{B},$
(6)

where $n_{0}$ and $r_{0}$ set the scale of the particle density and radius.

A normalization of the dimensionless wave function $\chi(\xi)$ defines the
mean particle density $n$ in the total volume $4\pi R^{3}/3$:

$\sigma\equiv\frac{3}{\xi^{3}_{B}}\int_{0}^{\xi_{B}}\chi^{2}(\xi)\,\xi^{2}\,\mathrm{d}\xi=\frac{n}{n_{0}}.$
(7)

Empirically, the value of $\sigma$ is less than unity [23, 24].

The wave function $\chi(\xi)$ determines the local particle density
$\eta(\xi)=\chi^{2}(\xi)$ and is found from the field equation

$\frac{1}{2}(\Delta_{\xi}+\nu)\chi-A\chi\Phi-Q\chi^{3}-B\chi^{5}=0,$ (8)

where dimensionless parameters of the model are

$Q=4\pi r^{2}_{0}an_{0},\quad A=Q\,\frac{Gm^{3}r^{2}_{0}}{\hbar^{2}a},\quad
B=g\,Q\,\lambda^{3}n_{0},\quad\varepsilon_{0}=\frac{\hbar^{2}}{mr^{2}_{0}}\,n_{0}.$
(9)

The effective chemical potential $\nu$ plays the role of free parameter
absorbing the gravitational potential in the center (at $\xi=0$):

$\upsilon=-A\int_{0}^{\xi_{B}}\eta(\xi)\,\xi\,\mathrm{d}\xi,$ (10)

and it leads to the following expression for $\Phi(\xi)$:

$\Phi(\xi)=-\frac{1}{\xi}\int_{0}^{\xi}\eta(s)\,s^{2}\,\mathrm{d}s+\int_{0}^{\xi}\eta(s)\,s\,\mathrm{d}s,\hskip
25.0pt\Delta_{\xi}\Phi(\xi)=\eta(\xi).$ (11)

As usual [27], we associate $r_{0}$ with the healing length:

$r_{0}=\frac{1}{\sqrt{4\pi an_{0}}},$ (12)

what leads immediately to $Q=1$ in (9). Under the diluteness conditions (4)
which guarantee applicability of the given formalism, the following inequality
holds:

$r^{3}_{0}n_{0}\gg 1.$ (13)

Let us relate the characteristics of dark particles with the parameters of our
model, using the typical mass density value $\rho_{0}=10^{-20}\
\text{kg}/\text{m}^{3}$, that gives

$n_{0}=\frac{\rho_{0}}{m}.$ (14)

Then, substituting it and (12) into (9), we obtain the relations

$\displaystyle m=\left(\frac{g}{B}\right)^{1/4}m_{0},\quad
m_{0}\equiv\left(\frac{\hbar^{3}}{c^{3}}\,\rho_{0}\right)^{1/4}=0.081\
\text{eV}/c^{2},$ (15) $\displaystyle a=\left(\frac{g}{4\pi
AB}\right)^{1/2}\,\ell_{\text{P}},\quad\ell_{\text{P}}\equiv\sqrt{\frac{G\hbar}{c^{3}}}=1.62\times
10^{-35}\ \text{m},$ (16) $\displaystyle
r_{0}=\left(\frac{B}{g}\right)^{1/8}\,\left(\frac{A}{4\pi}\right)^{1/4}\,L,\quad
L\equiv\sqrt{\frac{m_{0}}{\rho_{0}\,\ell_{\text{P}}}}=9.44\times 10^{8}\
\text{m}.$ (17)

It is required that $mc^{2}\leq 1\ \text{eV}$ accordingly to [19]. For
definiteness, let us fix $g=4.73\times 10^{-83}$ and

$Q=1,\quad A=\pi^{2}\simeq 9.87,\quad B=20,$ (18)

to obtain that $a=2.23\times 10^{-78}\ \text{m}$, $m=10^{-22}\
\text{eV}/c^{2}$, $r_{0}=0.818\ \text{kpc}$.

Of course, we admit a more suitable parametrization for concrete situations.
Besides, the core size $R=r_{0}\xi_{B}$ and the mass density
$\rho=\rho_{0}\sigma$ should be found by solving the field equations for given
restrictions. However, here we focus on analyzing the physical effects without
applying the model to specific objects.

Equations (8), (11) are numerically integrated under the following initial
conditions: $\chi^{\prime}(0)=0$ and $\chi^{\prime\prime}(0)<0$. As it was
argued in [24], the finite initial value $\chi(0)=z$ lies in the range
$z_{1}<z<z_{2}$, where

$z_{1}=\left[\sqrt{\left(\frac{3Q}{10B}\right)^{2}+\frac{\nu}{10B}}-\frac{3Q}{10B}\right]^{1/2},\quad
z_{2}=\left[\sqrt{\left(\frac{Q}{2B}\right)^{2}+\frac{\nu}{2B}}-\frac{Q}{2B}\right]^{1/2},$
(19)

and it should be minimal positive solution of equation:

$\frac{A}{2}z^{2}-\left(5Bz^{4}+3Qz^{2}-\frac{\nu}{2}\right)\left(\frac{\nu}{2}-Qz^{2}-Bz^{4}\right)=0.$
(20)

The absence of such a solution $z$ at given $(A,Q,B,\nu)$ means that
$\chi(\xi)=0$ everywhere. It happens at $\nu<\nu_{\text{min}}$ where
$\nu_{\text{min}}(A,Q,B)$ is found numerically.

Studying the statistical properties of matter, one uses the part of
oscillating wave function $\chi(\xi)$ in the range $\xi\in[0;\xi_{B}]$,
limited by its first zero $\chi(\xi_{B})=0$. The value of $\xi_{B}$ determines
the total radius $R=r_{0}\xi_{B}$ of the system and depends on the set of
parameters $(A,Q,B,\nu)$. Let us write down the internal energy density ${\cal
E}$ and the local pressure ${\cal P}$ in dimensionless units:

${\cal E}(\xi)=\frac{Q}{2}\eta^{2}(\xi)+\frac{B}{3}\eta^{3}(\xi),\quad{\cal
P}(\xi)=\frac{Q}{2}\eta^{2}(\xi)+\frac{2}{3}B\eta^{3}(\xi).$ (21)

The statistical mean of any local characteristic is defined by the integral:

$\langle(\ldots)\rangle=\frac{3}{\xi^{3}_{B}}\int_{0}^{\xi_{B}}(\ldots)\xi^{2}\mathrm{d}\xi.$
(22)

In general, there are two limiting regimes of the model (8), (11): i) the
Thomas–Fermi (TF) approximation, when the terms of quantum fluctuations like
$(\partial_{\xi}\chi)^{2}$ and $\Delta_{\xi}\chi$ are neglected; ii) the
regime of the standing spherical wave, when the chemical potential $\nu$ is
large enough and/or the interaction terms are negligibly small in (8).

Since our model generalizes the well-studied model with pair interaction
($B=0$) in the Thomas–Fermi approximation [20], it is useful to indicate its
basic equations and the values of the parameters. The set of this model
equations contains Equation (11) and

$\frac{\nu}{2}-A\Phi-Q\eta=0.$ (23)

Solution to these equations is given (see [20, 22]) by

$\eta(\xi)=\eta_{c}\,j_{0}(\kappa\xi),\quad\kappa=\sqrt{\frac{A}{Q}},\quad\eta_{c}=\frac{\nu}{2Q},$
(24)

where $j_{0}(z)=\sin{z}/z$ is the spherical Bessel function. In fact, the
central density $\eta_{c}$ is regarded as a free parameter of the model
instead of $\nu$.

Limiting the system radius by the value $r=r_{0}$ (or $\xi_{B}=1$), one
derives that $\kappa=\pi$. Hence, for arbitrary $\nu>0$, the parameters of
such a model can be fixed as $Q=1$ and $A=\pi^{2}$. This confirms our choice
of the value of $A$ in (18).

Note that close parameter values, $Q\sim 1.36$ and $A=10$, were used in [24].
However, taking into account the quantum fluctuations $\Delta_{\xi}\chi$ and
the three-particle interaction regulated by parameter $B$, the system radius
$\xi_{B}>1$ and the central density $\eta_{c}=\chi^{2}(0)$ reveal nonlinear
dependence on $\nu$ (see Figure 1 below and Figure 3 in [23]).

Let us note the relation $r_{0}(a,m)$:

$r_{0}=\pi\sqrt{\frac{\hbar^{2}a}{Gm^{3}}},$ (25)

which is discussed in detail in [20] to outline physically reasonable values
of $a$ and $m$ by describing the dark matter of dwarf galaxies.

In contrast with the Thomas–Fermi approximation, equation (8) admits the
limiting regime of the standing wave, implied by the equation:

$(\Delta_{\xi}+k^{2})\chi(\xi)=0,\hskip 25.0ptk^{2}=\nu-w>0,$ (26)

where $w$ is a correction due to small interaction, which is not derived here.

Figure 1: Dependence of the dimensionless radius $\xi_{B}$ (the first zero of
$\chi(\xi)$) on effective chemical potential $\nu$. Here $A=\pi^{2}$, $Q=1$,
$B=20$, $\nu_{\text{min}}\simeq 2.09$. The approximation (27) is built at
$w=0.6$. The chemical potential $\mu=0$ at $\nu\simeq 2.263$, defining the
lower limit of applicability of the approximation.

This equation has a simple analytic solution:

$\chi(\xi)=\chi_{0}\,j_{0}\left(\pi\frac{\xi}{\xi_{B}}\right),\hskip
25.0pt\xi_{B}=\frac{\pi}{\sqrt{\nu-w}}.$ (27)

Obvious merging of the approximate function $\xi_{B}(\nu)$ and the numerically
found one at large $\nu$ is shown in Figure 1 for the fixed values of $Q$ and
$A$.

Besides, the form of wave function (27) guarantees that
$\langle\eta^{n}(\xi)\rangle=\eta^{2}_{c}\,c_{n}$, where
$\eta_{c}\equiv\chi^{2}_{0}$, and the numeric coefficients $c_{n}$ are
independent of $\nu\gg\nu_{\text{min}}$. Although the exact situation is only
qualitatively the same, as we shall see in Figure 3 below, this effect
excludes a possibility of the composites formation because of destruction of
many-particle states at large both $\nu$ and particle momentum
$k=\sqrt{\langle(\partial_{\xi}\chi)^{2}\rangle/\sigma}$. As it is easily
seen, the averages of $\eta^{n}$ with similar properties are obtained for the
solution (24) in the Thomas–Fermi approximation.

Referring to [24], we especially note that the successful attempts of
reproducing rotational curves of dwarf galaxies require the use of
$\nu\simeq\nu_{\text{min}}$, what can be associated with a ground state of the
model (8), (11).

Moreover, the two regimes of “small” and ”large” $\nu$ are naturally separated
by condition $\mu=0$ for the physical chemical potential or its dimensionless
value $u\equiv\nu/2-A\upsilon=0$. In other words, the regime of large $\nu$
corresponds to $u>0$, while the bound states are described at $u<0$. As we
have seen in [24, 23], the point $u=0$ determines parameters of the first-
order phase transition in the system with three-particle interaction.

After matching the parameters of our model and the reference one, we would
also like to compare some dependencies. We start by examining the mean
particle density $\sigma$ as a function of the central density $\eta_{c}$, at
large values of which the transformation of particles can occur.

As already shown in [20] for model (23)-(24) with $\xi_{B}=1$, the relation
between $\sigma$ and $\eta_{c}$ is given by

$\sigma=\frac{3}{\pi^{2}}\,\eta_{c}.$ (28)

Figure 2: Mean particle density $\sigma$ (a) and system radius $\xi_{B}$ (b)
as functions of central density $\eta_{c}$ in two models. In both panels, the
orange line corresponds to the model with pair interaction in the Thomas–Fermi
approximation. The red and blue branches of the single curve correspond to the
model that takes into account both quantum fluctuations and three-particle
interaction.

A completely different relation is expected for our model (8). The
dependencies of interest are shown in Figure 2a.

The first difference between the model with pair interaction in the
Thomas–Fermi approximation (23)-(24) and our model which takes into account
both quantum fluctuations and three-particle interaction is the presence of
the minimal admissible central density in the latter, which for the given
model parameters is equal to $\eta^{*}_{c}\simeq 0.1294$.

At $\eta_{c}\geq\eta^{*}_{c}$, there are two branches of the function
$\sigma(\eta_{c})$ for the model (8), which are colored in blue and red in
Figure 2a. The behavior of the blue branch is similar to the orange line
described by (28), while the red curve looks abnormal due to three-particle
repulsion. To explain the discrepancy between the mean density curves colored
in blue and orange as well as the red branch, let us turn to the dependence of
the size $\xi_{B}$ on $\eta_{c}$.

We see in Figure 2b that an increase in the system size $\xi_{B}$ along red
branch is accompanied by a decrease in the corresponding mean density in
Figure 2a. That is, the red sectors in Figures 2a and 2b describe quite dilute
matter.

Since the blue curve in Figure 2b describes a system with a larger radius
$\xi_{B}$ than the orange line for the reference model, this is consistent
with the lower density in blue compared to the mean density in orange in
Figure 2a. As $\eta_{c}$ grows, the blue curve $\xi_{B}(\eta_{c})$ tends to
$\xi_{B}=1$. Therefore, the blue curves do show a denser state of matter.

On the other hand, the existence of composites in this dense matter becomes
problematic due to the mutual influence leading to their dissociation into
particles. Assuming that the interactions still guarantee the formation of
composites for a short time, which is similar to oscillations, we will
consider such a process in the next Section using the model of quantum random
walk.

To test our assumptions, we introduce the ratios of the dimensionless means:

$C_{1}=\frac{\langle\eta^{2}\rangle}{\langle\eta\rangle^{2}},\quad
C_{2}=\frac{\langle\eta^{3}\rangle}{\langle\eta\rangle^{3}},\quad
C_{3}=\frac{\langle\eta^{3}\rangle}{\langle\eta\rangle^{2}\langle\eta\rangle},$
(29)

where the particle density $\eta(\xi)=\chi^{2}(\xi)$ is formed from the
solution $\chi(\xi)$ of (8).

Figure 3: Behavior of the ratios $C_{1,2,3}$ of mean values
$\langle\eta^{3}\rangle$, $\langle\eta^{2}\rangle$ and $\langle\eta\rangle$
[see (29)] with varying internal pressure $P$.

To examine the behavior of these quantities, found with the same parameters
$(A,Q,B)$ as before, but varying $\nu$, we use the mean internal pressure $P$
(see [24]):

$P=\frac{Q}{2}\langle\eta^{2}\rangle+\frac{2}{3}B\langle\eta^{3}\rangle.$ (30)

The quotients $C_{1,2,3}(P)$, presented in Figure 3, reveal the splitting and
equivalence:

$\langle\eta^{2}\rangle\sim\langle\eta\rangle^{2},\quad\langle\eta^{3}\rangle\sim\langle\eta\rangle^{3},\quad\langle\eta^{3}\rangle\sim\langle\eta^{2}\rangle\langle\eta\rangle,$
(31)

described by almost horizontal plateau at relatively high pressure
$P>1.8\times 10^{-3}$. (To obtain the value of pressure and energy density in
physical units, the dimensionless quantity should be multiplied by
$\varepsilon_{0}=\hbar^{2}n_{0}/(mr^{2}_{0})$, see (9).) This means that two
and three-particle configurations behave as two and three isolated particles,
respectively. The reasons of such a behavior are explained above by
considering the standing wave solution (27). The transition to this regime is,
generally speaking, a first-order phase transition described in [23, 24].

Figure 4: Dependence of the measure of quantum fluctuations $\tau$ (a) and
mean density (b) on internal pressure $P$.

It is natural to associate the instability of many-particle configurations
(composites) with quantum fluctuations in the system. Indeed, the quantity
$\tau=\langle(\partial_{\xi}\chi)^{2}\rangle$ undergoes a jump in Figure 4a
when passing to the pressure interval of interest. Moreover, at $P>1.8\times
10^{-3}$, the almost linear proportionality between $P$ and $\tau$ resembles
the relation between pressure and temperature of an ideal gas.

In summary, we can conclude that the high pressure regime does not provide the
formation of composites in BEC, although the particle density increases with
the pressure, as shown in Figure 4b. Besides, the inequality $nr^{3}_{0}\gg 1$
is fulfilled in BEC due to large healing length $r_{0}$, which determines the
BEC core radius [27]. Despite the strong correlations, characterized by
$r_{0}$, the diluteness condition $na^{3}<1$ persists.

##  3 Quantum Random Walk Based Treatment

Based on the predictions of the previous Section, we would like to consider
the instability of two- and three-particle states in dense matter within the
framework of a simple probabilistic model. In the case of an elementary system
of three particles (monomers, denoted M), which are confined within the radius
of interactions leading to the formation of a dimer (D) and a trimer (T),
there is a possibility of oscillations, namely, a continuous transition
between three states:

$\displaystyle\psi_{1}=|3\rangle_{\mathrm{M}}\bigotimes|0\rangle_{\mathrm{D}}\bigotimes|0\rangle_{\mathrm{T}},\hskip
25.0pt\psi_{2}=|1\rangle_{\mathrm{M}}\bigotimes|1\rangle_{\mathrm{D}}\bigotimes|0\rangle_{\mathrm{T}},$
(32)
$\displaystyle\psi_{3}=|0\rangle_{\mathrm{M}}\bigotimes|0\rangle_{\mathrm{D}}\bigotimes|1\rangle_{\mathrm{T}},$

representing the unbound state and two bound ones. In the absence of external
factors, the reversibility of microscopic processes is supported by the
equality of the probabilities of forming and decaying of the composites. Of
course, the known lifetimes of each (metastable) state could point out the
dominant configuration of the particles. However, having permanently
interacting particles, we will neglect the probability of their adhesion.
Nevertheless, it seems possible to estimate the probabilities of states over
the period of oscillations.

Thus, we consider a toy model with three identical monomers, which is
described by the wave function $\Psi(t)$:

$\Psi(t)=\sum\limits_{i=1}^{3}a_{i}(t)\psi_{i},\hskip
25.0pt\sum\limits_{i=1}^{3}|a_{i}(t)|^{2}=1,$ (33)

that evolves in time $t$ according to the Schrödinger equation with initial
state (at $t=0$) assumed to be $\Psi_{0}=\psi_{1}$.

Then we define the fractions $\Omega_{i}(t)=|a_{i}(t)|^{2}$ and analyze their
means $\overline{\Omega}_{i}$ over a period of time.

Evolution of the initial state,

$|\Psi_{0}\rangle=\begin{pmatrix}1\\\ 0\\\ 0\end{pmatrix}\xrightarrow{\quad
t\quad}|\Psi(t)\rangle=\begin{pmatrix}a_{1}(t)\\\ a_{2}(t)\\\
a_{3}(t)\end{pmatrix},$ (34)

is given by the Schrödinger equation with Hamiltonian
$\hat{H}=\hat{H}_{0}+\hat{U}$, where
$\hat{H}_{0}=\mathrm{diag}(\varepsilon_{1},\varepsilon_{2},\varepsilon_{3})$:

$\mathrm{i}\hbar\frac{\mathrm{d}}{\mathrm{d}t}|\Psi\rangle=\hat{H}|\Psi\rangle.$
(35)

Actually, the matrix $\hat{H}$ with constant elements determines quantum
random walk (a Markov process). Analytically, solution of (35) for the wave
function can be immediately represented (and found) in the form:

$|\Psi(t)\rangle=\hat{S}(t)|\Psi_{0}\rangle,\hskip
25.0pt\hat{S}(t)=\exp{\left(-\frac{\mathrm{i}}{\hbar}\hat{H}t\right)}.$ (36)

Neglecting spatial evolution of localized particles and time delay in each
state, we omit $\hat{H}_{0}$ and remain with the matrix of admissible
transitions:

$\hat{U}=\begin{pmatrix}0&U_{-1}&U_{-3}\\\ U_{1}&0&U_{-2}\\\
U_{3}&U_{2}&0\end{pmatrix},$ (37)

where the quantities $U_{\pm i}$ play a role of transition amplitudes, and we
put $U_{-i}=U_{i}\in\mathbb{R}$ to make $\hat{U}$ self-adjoint with real
elements.

The conserved energy of the system is

$E\equiv\langle\Psi(t)|{\hat{H}}|\Psi(t)\rangle=U_{1}a^{*}_{2}(t)a_{1}(t)+U_{2}a^{*}_{3}(t)a_{2}(t)+U_{3}a^{*}_{1}(t)a_{3}(t)+\text{c.c.}$
(38)

Introducing auxiliary quantities:

$U^{(3)}=U_{1}U_{2}U_{3},\hskip
25.0ptU^{(2)}=\frac{1}{3}(U^{2}_{1}+U^{2}_{2}+U^{2}_{3}),\\\ $ (39)

three eigenvalues of $\hat{U}$, which are roots of equation

$0=\det\left(\hat{U}-\lambda\hat{I}\right)=-\lambda^{3}+3\lambda
U^{(2)}+2U^{(3)},$ (40)

are found to be

$\displaystyle\lambda_{1}=A+B,\hskip
25.0pt\lambda_{2}=\mathrm{e}^{-2\pi\mathrm{i}/3}A+\mathrm{e}^{2\pi\mathrm{i}/3}B,$
$\displaystyle\lambda_{3}=\mathrm{e}^{2\pi\mathrm{i}/3}A+\mathrm{e}^{-2\pi\mathrm{i}/3}B,$
(41)

where

$A=\left[U^{(3)}+\sqrt{\left(U^{(3)}\right)^{2}-\left(U^{(2)}\right)^{3}}\right]^{1/3},\
B=\frac{U^{(2)}}{A}.$ (42)

General solution to equation (35) is presented in the form:

$|\Psi(t)\rangle=\sum\limits_{i=1}^{3}C_{i}\,{\bf
u}_{i}\,\exp{\left(-\mathrm{i}\frac{\lambda_{i}t}{\hbar}\right)},\hskip
25.0pt{\bf u}_{i}=\begin{pmatrix}\lambda^{2}_{i}-U^{2}_{2}\\\
\lambda_{i}U_{1}+U_{2}U_{3}\\\ \lambda_{i}U_{3}+U_{1}U_{2}\end{pmatrix},$ (43)

where ${\bf u}_{i}$ are the non-normalized eigenvectors of matrix $\hat{U}$;
$C_{i}$ are the constants, determined by initial condition
$|\Psi(0)\rangle=|\Psi_{0}\rangle$ from (34):

$C_{1}=\mathrm{i}\frac{\sqrt{3}}{9}\frac{\lambda_{2}-\lambda_{3}}{A^{3}-B^{3}},\quad
C_{2}=\mathrm{i}\frac{\sqrt{3}}{18}\frac{\lambda_{3}-\lambda_{2}-3\lambda_{1}}{A^{3}-B^{3}},\quad
C_{3}=\mathrm{i}\frac{\sqrt{3}}{18}\frac{\lambda_{3}-\lambda_{2}+3\lambda_{1}}{A^{3}-B^{3}}.$
(44)

We omit the study of time dependence in three-parameter case (43) and restrict
ourselves to a particular case, setting $U_{3}=0$ in (42), when

$A=3^{-1/2}\,\mathrm{e}^{\mathrm{i}\pi/6}\sqrt{U^{2}_{1}+U^{2}_{2}}.$

We easily obtain the spectrum:

$\lambda_{1}=\sqrt{U^{2}_{1}+U^{2}_{2}},\ \ \lambda_{2}=0,\ \
\lambda_{3}=-\sqrt{U^{2}_{1}+U^{2}_{2}}.$ (45)

Physically, these quantities determine the energy levels in the system and
indicate the presence of the unbound state of 3M (three monomers) with
$\lambda_{1}$, the bound one of T with $\lambda_{3}$, and a transient one with
$\lambda_{2}$, associated with M+D.

It is useful to parametrize $U_{1}$ and $U_{2}$ by $s$ and $T$ as

$U_{1}=J\cos{s},\quad U_{2}=J\sin{s},\quad J=\hbar\omega=\frac{2\pi\hbar}{T}.$
(46)

Parameter $T$ defines the time scale of the problem, while the quantities
$\cos{s}/T$ and $\sin{s}/T$ represent the rates of the processes.

Under the imposed conditions, the components of the matrix $\hat{S}(t)$ are
easily calculated and yield the instant fractions:

$\displaystyle\Omega_{\mathrm{M}}(t,s)=\left[\cos^{2}{s}\,\cos{(\omega
t)}+\sin^{2}{s}\right]^{2},\hskip
25.0pt\Omega_{\mathrm{D}}(t,s)=\left[\cos{s}\,\sin{(\omega t)}\right]^{2},$
(47)
$\displaystyle\Omega_{\mathrm{T}}(t,s)=\left[\cos{s}\,\sin{s}\,\left(\cos{(\omega
t)}-1\right)\right]^{2}.$

To obtain this, we can also use in (43) that

$C_{1}{\bf u}_{1}=\frac{1}{2}\begin{pmatrix}\cos^{2}{s}\\\ \cos{s}\\\
\cos{s}\,\sin{s}\end{pmatrix},\ \ C_{2}{\bf
u}_{2}=\begin{pmatrix}\sin^{2}{s}\\\ 0\\\ -\cos{s}\,\sin{s}\end{pmatrix},\ \
C_{3}{\bf u}_{3}=\frac{1}{2}\begin{pmatrix}\cos^{2}{s}\\\ -\cos{s}\\\
\cos{s}\,\sin{s}\end{pmatrix}.$ (48)

Note that configuration with $\Omega_{\mathrm{M}}=0$ and
$\Omega_{\mathrm{T}}>\Omega_{\mathrm{D}}$ is realized at the time moment
(during the first period)

$t_{1}=\frac{T}{2}-\frac{T}{2\pi}\arccos{\left(\tan^{2}{s}\right)},$ (49)

when $s$ tends to $\pi/4$ from below.

Figure 5: Time-averaged fractions $\overline{\Omega}_{\alpha}$ as functions of
parameter $s=\mathrm{arctan}\left(U_{2}/U_{1}\right)$ at $U_{3}=0$.

Averaging $\Omega_{\alpha}(t)$ over the period $T$, we arrive at

$\displaystyle\overline{\Omega}_{\mathrm{M}}(s)=\frac{3}{2}\cos^{4}{s}-2\cos^{2}{s}+1,$
$\displaystyle\overline{\Omega}_{\mathrm{D}}(s)=\frac{1}{2}\cos^{2}{s},\quad\overline{\Omega}_{\mathrm{T}}(s)=\frac{3}{2}\cos^{2}{s}\,\sin^{2}{s}.$
(50)

Comparison of these fractions for different $s$ is made in Figure 5. We see a
slight dominance of $\overline{\Omega}_{\mathrm{T}}$ in the range
$s_{c}<s<\pi/4$, while $\overline{\Omega}_{\mathrm{D}}$ prevails for
$s<s_{c}$, where $s_{c}=\mathrm{arcsin}\left(\sqrt{3}/3\right)\simeq
3.13\pi/16$, and
$\overline{\Omega}_{\mathrm{M}}(s_{c})=\overline{\Omega}_{\mathrm{D}}(s_{c})=\overline{\Omega}_{\mathrm{T}}(s_{c})=1/3$.

Additionally, let us indicate the changes in the energy spectrum when
$\hat{H}_{0}\not=0$. Giving $\hat{H}_{0}=\mathrm{diag}(\mu,\mu,\mu)$, we get a
shift in (45): $\lambda_{i}\to\lambda_{i}+\mu$. Including
$\hat{H}_{0}=\mathrm{diag}(0,-2J\sinh(a),0)$ leads to rescaling:
$\lambda_{1}=J\mathrm{e}^{-a}$, $\lambda_{2}=0$,
$\lambda_{3}=-J\mathrm{e}^{a}$ [see (46)]. Dynamics in these cases will be
somewhat more complicated and we omit it.

In summary, we note that in a three-particle system with permanent
oscillations between quantum states there exists the possibility of the
composites dominance in a long-time picture. This effect is due to the tuning
of interaction parameters. In general, the probability of detecting monomers
(of dark matter) increases at higher density and pressure in the model with
few-particle interactions.

##  4 Chemical Kinetics Consideration

For the sake of simplicity, let us assume that the elementary constituents of
dark matter are identical spinless particles (monomers, or M). Due to an
interaction between them, the nature of which is not specified herein, we
assume the ability of formation of dimers (D) and trimers (T) in the ground
state. Therefore, we restrict ourselves to processes involving no more than
three monomers, taking into account also those contained in dimers and
trimers. Attempting to describe the formation/decay of dimers and trimers, as
well as their fractions, we clarify some physical conditions.

Account of the gravitational interaction, resulting in an inhomogeneous
distribution of particles, implies that the statistical characteristics of the
matter should, in principle, depend on space. We assume that the particle
density is highest in a core and tends to zero at halo periphery. However, in
regions of space that are small in astronomical scale, but large enough to
apply the statistical approach, the gravitational potential (together with the
chemical one) can be considered as constant one. This allows us to explore the
properties of quasi-homogeneous subsystems. However, simplifying the problem
by introducing regions of homogeneity with different particle densities, the
dominance of certain processes does depend on the ratio of the interparticle
distance to the radius of interaction (or scattering lengths).

The particle density in the core can be estimated as

$n_{0}\simeq 5.61\cdot 10^{15}\,\mathrm{m}^{-3}\,\left[\frac{mc^{2}}{1\
\mathrm{eV}}\right]^{-1}\left[\frac{\rho_{0}}{10^{-20}\
\mathrm{kg}\,\mathrm{m}^{-3}}\right],$ (51)

which is determined by the individual bosonic monomer mass $m\lesssim 1$ eV
and the typical mass density $\rho_{0}$.

Let us consider possible reactions occurring simultaneously in each subsystem:

$2\mathrm{M}\xrightleftharpoons[k_{-1}]{k_{1}}\mathrm{D},\hskip
25.0pt\mathrm{M}+\mathrm{D}\xrightleftharpoons[k_{-2}]{k_{2}}\mathrm{T},\hskip
25.0pt3\mathrm{M}\xrightleftharpoons[k_{-3}]{k_{3}}\mathrm{T},$ (52)

where the positive constants $k_{i}$ and $k_{-i}$ ($i=1,2,3$) are related to
the rates of the corresponding reactions and have the dimension of inverse
time.

The time evolution of concentrations $C_{\alpha}$ (the numbers of particles of
sort $\alpha=\mathrm{M},\mathrm{D},\mathrm{T}$ in a fixed volume) can be
described accordingly to the rules of chemical kinetics [40]:

$\displaystyle\frac{\mathrm{d}C_{\mathrm{M}}}{\mathrm{d}t}$ $\displaystyle=$
$\displaystyle-2k_{1}V_{\mathrm{2M}}C^{2}_{\mathrm{M}}+2k_{-1}C_{\mathrm{D}}-k_{2}V_{\mathrm{MD}}C_{\mathrm{M}}C_{\mathrm{D}}$
$\displaystyle+k_{-2}C_{\mathrm{T}}-3k_{3}V^{2}_{\mathrm{3M}}C^{3}_{\mathrm{M}}+3k_{-3}C_{\mathrm{T}},$
$\displaystyle\frac{\mathrm{d}C_{\mathrm{D}}}{\mathrm{d}t}$ $\displaystyle=$
$\displaystyle
k_{1}V_{\mathrm{2M}}C^{2}_{\mathrm{M}}-k_{-1}C_{\mathrm{D}}-k_{2}V_{\mathrm{MD}}C_{\mathrm{M}}C_{\mathrm{D}}+k_{-2}C_{\mathrm{T}},$
(53) $\displaystyle\frac{\mathrm{d}C_{\mathrm{T}}}{\mathrm{d}t}$
$\displaystyle=$ $\displaystyle
k_{2}V_{\mathrm{MD}}C_{\mathrm{M}}C_{\mathrm{D}}-k_{-2}C_{\mathrm{T}}+k_{3}V^{2}_{\mathrm{3M}}C^{3}_{\mathrm{M}}-k_{-3}C_{\mathrm{T}},$

with the initial conditions at $t=0$: $C^{(0)}_{\mathrm{M}}>0$,
$C^{(0)}_{\mathrm{D,T}}=0$; for $V_{\mathrm{2M}}$, $V_{\mathrm{3M}}$, and
$V_{\mathrm{MD}}$, see next proposition.

The above set of equations for $C_{\alpha}$ conserves in time the total number
of monomers in the subsystem of a fixed volume,

$C_{\mathrm{M}}(t)+2C_{\mathrm{D}}(t)+3C_{\mathrm{T}}(t)=\mathrm{const},$ (54)

and implies that the elementary act depends on constants $k_{\pm i}$ and
interaction volumes $V_{s}\simeq a^{3}_{s}$, where $a_{\mathrm{2M}}$,
$a_{\mathrm{3M}}$ and $a_{\mathrm{MD}}$ are the scattering lengths (for
$s$-wave channel, in instance). It is natural to require conditions for the
successful performing of three forward reactions:

$\displaystyle V_{\mathrm{2M}}C_{\mathrm{M}}\geq 1,\quad
V_{\mathrm{3M}}C_{\mathrm{M}}\geq 1,$ $\displaystyle
V_{\mathrm{MD}}C_{\mathrm{M}}\geq 1\quad\mathrm{and}\quad
V_{\mathrm{MD}}C_{\mathrm{D}}\geq 1.$ (55)

They indicate the presence of at least one particle in each interaction volume
$V_{s}$. Failure to meet one of these conditions may lead to the exclusion of
the corresponding process from consideration.

Associating $C_{\alpha}$ with the local particle density far from the core (of
galactic halo) center, in order to satisfy the imposed conditions (55), we
require that $n_{0}a^{3}_{s}\gg 1$ for the scattering lengths $a_{s}$ in the
model.

Strictly speaking, the validity of condition $na^{3}_{s}>1$ depends on both
$n$ and $a_{s}$. However, as shown above, an increase in monomer density $n$
(due to squeezing, for example) correlates with the destruction of composites.
Therefore, there remains the only possibility of considering interactions with
longer scattering lengths $a_{s}$. Next, we evaluate the contribution of
dimers and trimers to dark matter.

To simplify further notations, we will also use

$K_{1}=k_{1}V_{\mathrm{2M}},\ K_{2}=k_{2}V_{\mathrm{MD}},\
K_{3}=k_{3}V^{2}_{\mathrm{3M}},\ K_{-i}=k_{-i}.$

It is convenient to describe the subsystem evolution in terms of the
fractions:

$\Omega_{\alpha}=\frac{c_{\alpha}C_{\alpha}}{C_{\mathrm{M}}+c_{\mathrm{D}}C_{\mathrm{D}}+c_{\mathrm{T}}C_{\mathrm{T}}},\hskip
25.0pt\Omega_{\mathrm{M}}+\Omega_{\mathrm{D}}+\Omega_{\mathrm{T}}=1,$ (56)

where $c_{\mathrm{M}}=1$, $c_{\mathrm{D}}=2$, $c_{\mathrm{T}}=3$.

Let us analyze the equilibrium fractions $\Omega^{\mathrm{eq}}_{\alpha}$ when
$\mathrm{d}C_{\alpha}/\mathrm{d}t=0$. For equilibrium concentrations (dropping
the label ‘eq’), we find $C_{\mathrm{D}}$ and $C_{\mathrm{T}}$ as functions of
$C_{\mathrm{M}}$:

$\displaystyle
C_{\mathrm{D}}=C^{2}_{\mathrm{M}}\,f_{\mathrm{D}}(C_{\mathrm{M}}),\hskip
25.0pt\quad
C_{\mathrm{T}}=C^{3}_{\mathrm{M}}\,f_{\mathrm{T}}(C_{\mathrm{M}}),$ (57)
$\displaystyle
f_{\alpha}(x)=\frac{a^{(11)}_{\alpha}x+a^{(12)}_{\alpha}}{a^{(21)}_{\alpha}x+a^{(22)}_{\alpha}},\quad
A_{\alpha}=\left(\begin{array}[]{cc}a^{(11)}_{\alpha}&a^{(12)}_{\alpha}\\\
a^{(21)}_{\alpha}&a^{(22)}_{\alpha}\end{array}\right),$ (60)

where the elements of matrix $A_{\alpha}$ are formed from $K_{\pm i}$:

$\displaystyle
A_{\mathrm{D}}=\left(\begin{array}[]{cc}K_{3}K_{-2}&K_{1}(K_{-2}+K_{-3})\\\
K_{2}K_{-3}&K_{-1}(K_{-2}+K_{-3})\end{array}\right),$ (63) $\displaystyle
A_{\mathrm{T}}=\left(\begin{array}[]{cc}K_{2}K_{3}&K_{1}K_{2}+K_{3}K_{-1}\\\
K_{2}K_{-3}&K_{-1}(K_{-2}+K_{-3})\end{array}\right).$ (66)

For a subsystem in a fixed (unit) volume, the initial concentration
$C^{(0)}_{\mathrm{M}}$ is related to the equilibrium one $C_{\mathrm{M}}$:

$C_{\mathrm{M}}+2C^{2}_{\mathrm{M}}\,f_{\mathrm{D}}(C_{\mathrm{M}})+3C^{3}_{\mathrm{M}}\,f_{\mathrm{T}}(C_{\mathrm{M}})=C^{(0)}_{\mathrm{M}}.$
(67)

That means, the initial state $C^{(0)}_{\mathrm{M}}$ can be determined by
using the $C_{\mathrm{M}}$.

Note the possibility to generalize our considerations by specifying the
concentration $C^{(q)}_{\alpha}$ of composites of sort $\alpha$ in certain
quantum state $q$ at a given temperature $T$, that is, the occupation of
energy level $E^{(q)}_{\alpha}$ with the degeneracy $g_{q}$. For the processes
that do not violate the equilibrium distribution of particles over states, one
writes [40]:

$C^{(q)}_{\alpha}=C_{\alpha}p_{\alpha}(q),\quad
p_{\alpha}(q)=\frac{g_{q}}{Z_{\alpha}}\exp{\left(-\frac{E^{(q)}_{\alpha}}{T}\right)},\quad\sum\limits_{q}p_{\alpha}(q)=1,$
(68)

where $Z_{\alpha}$ is the partition function over states, and $C_{\alpha}$ is
the total concentration. From equation (53), with account of Boltzmann factors
$p_{\alpha}(q)$, one can deduce the dependence of the reaction rates $k_{\pm
i}$ on non-zero temperature.

At finite $K_{\pm i}$, there are two limiting regimes: i) $C_{\mathrm{M}}\to
0$ (low concentration of monomers) admits approximation:
$C^{(0)}_{\mathrm{M}}=C_{\mathrm{M}}+O(C^{2}_{\mathrm{M}})$, and

$\displaystyle\Omega^{\mathrm{eq}}_{\mathrm{M}}=1-2\frac{K_{1}}{K_{-1}}\,C_{\mathrm{M}}+O(C^{2}_{\mathrm{M}}),$
$\displaystyle\Omega^{\mathrm{eq}}_{\mathrm{D}}=2\frac{K_{1}}{K_{-1}}\,C_{\mathrm{M}}+O(C^{2}_{\mathrm{M}}),$
(69)
$\displaystyle\Omega^{\mathrm{eq}}_{\mathrm{T}}=3\frac{K_{1}K_{2}+K_{3}K_{-1}}{K_{-1}(K_{-2}+K_{-3})}\,C^{2}_{\mathrm{M}}+O(C^{3}_{\mathrm{M}});$

ii) asymptotic expressions at large $C_{\mathrm{M}}$ yield:

$\displaystyle
C^{(0)}_{\mathrm{M}}=\frac{3K_{3}}{K_{-3}}\,C^{3}_{\mathrm{M}}+o(C^{2}_{\mathrm{M}}),$
$\displaystyle\Omega^{\mathrm{eq}}_{\mathrm{M}}=\frac{K_{-3}}{3K_{3}}\,C^{-2}_{\mathrm{M}}+O(C^{-3}_{\mathrm{M}}),\quad\Omega^{\mathrm{eq}}_{\mathrm{D}}=\frac{2K_{-2}}{3K_{2}}\,C^{-1}_{\mathrm{M}}+O(C^{-2}_{\mathrm{M}}),$
(70)
$\displaystyle\Omega^{\mathrm{eq}}_{\mathrm{T}}=1-\frac{2K_{-2}}{3K_{2}}\,C^{-1}_{\mathrm{M}}+O(C^{-2}_{\mathrm{M}}).$

The formulas (69)-(70) are obtained by expanding up to the first non-trivial
terms. Thus, the normalization
$\Omega^{\mathrm{eq}}_{\mathrm{M}}+\Omega^{\mathrm{eq}}_{\mathrm{D}}+\Omega^{\mathrm{eq}}_{\mathrm{T}}=1$
is valid for $\Omega^{\mathrm{eq}}_{\alpha}$ in the same approximation.

As seen from (69)-(70), the case i) indicates a rather low probability of
obtaining composites D and T, while the mode ii) shows the dominant role of
trimers (due to reaction $\mathrm{M}+\mathrm{D}\rightarrow\mathrm{T}$) and the
suppression of the monomer number because of converting
$3\mathrm{M}\rightarrow\mathrm{T}$. In this case we also see that
$C_{\mathrm{M}}\sim\left(C^{(0)}_{\mathrm{M}}\right)^{1/3}$.

Since the admissible value of $\Omega^{\mathrm{eq}}_{\mathrm{T}}$ varies in
the interval $(0;1)$, we can get $\Omega^{\mathrm{eq}}_{\mathrm{T}}$ higher
than the expected threshold $\Omega^{\mathrm{thr}}_{\mathrm{T}}\in(0;1)$ by
requiring $C_{\mathrm{M}}>C^{\mathrm{thr}}_{\mathrm{M}}$. For a given
$\Omega^{\mathrm{thr}}_{\mathrm{T}}$, the value of
$C^{\mathrm{thr}}_{\mathrm{M}}$ is determined by a real positive root $z$ of
the quartic equation:

$3z^{3}f_{\mathrm{T}}(z)\left(1-\Omega^{\mathrm{thr}}_{\mathrm{T}}\right)-\Omega^{\mathrm{thr}}_{\mathrm{T}}\left(z+2z^{2}f_{\mathrm{D}}(z)\right)=0.$
(71)

A similar problem arises when expressing $C_{\mathrm{M}}$ through
$C^{(0)}_{\mathrm{M}}$. Besides, one can also find the threshold value of
$C^{(0)}_{\mathrm{M}}$ when
$\Omega^{\mathrm{eq}}_{\mathrm{T}}=\Omega^{\mathrm{thr}}_{\mathrm{T}}$. To
demonstrate some analytical solutions, we restrict ourselves to simpler
special cases.

Let $K_{-1}=K_{-2}=K_{3}=0$. It means that the allowed processes are

$3\mathrm{M}\xrightarrow{k_{1}}\mathrm{M}+\mathrm{D}\xrightarrow{k_{2}}\mathrm{T}\xrightarrow{k_{-3}}3\mathrm{M}.$
(72)

The fractions in such a model is given by

$\displaystyle\Omega^{\mathrm{eq}}_{\mathrm{M}}=\frac{K_{2}K_{-3}}{K_{2}K_{-3}+2K_{1}K_{-3}+3K_{1}K_{2}C_{\mathrm{M}}},$
$\displaystyle\Omega^{\mathrm{eq}}_{\mathrm{D}}=\frac{2K_{1}K_{-3}}{K_{2}K_{-3}+2K_{1}K_{-3}+3K_{1}K_{2}C_{\mathrm{M}}},$
(73)
$\displaystyle\Omega^{\mathrm{eq}}_{\mathrm{T}}=\frac{3C_{\mathrm{M}}K_{1}K_{2}}{K_{2}K_{-3}+2K_{1}K_{-3}+3K_{1}K_{2}C_{\mathrm{M}}},\quad
C_{\mathrm{M}}={\cal C}_{+}\left(C^{(0)}_{\mathrm{M}}\right),$

where ${\cal C}_{+}$ is a positive solution among the two solutions ${\cal
C}_{\pm}$ to the quadratic equation obtained from (67) under imposed
restrictions:

${\cal
C}_{\pm}=\sqrt{\frac{K_{-3}}{3K_{1}}}\left[\pm\sqrt{C^{(0)}_{\mathrm{M}}+{\cal
C}}-\sqrt{{\cal C}}\right],\quad{\cal
C}=\frac{K_{-3}}{12K_{1}}\left(1+2\frac{K_{1}}{K_{2}}\right)^{2}.$ (74)

Therefore, $C_{\mathrm{M}}\sim\sqrt{C^{(0)}_{\mathrm{M}}}$ within this model.

To guarantee
$\Omega^{\mathrm{eq}}_{\mathrm{T}}\geq\Omega^{\mathrm{thr}}_{\mathrm{T}}$ for
$\Omega^{\mathrm{thr}}_{\mathrm{T}}=0.5$, we arrive at the condition:

$C_{\mathrm{M}}\geq\frac{K_{-3}}{3K_{1}}+\frac{2K_{-3}}{3K_{2}}\quad\mathrm{or}\quad
C^{(0)}_{\mathrm{M}}=4{\cal C}\,\gamma(\gamma+1),\quad\gamma\geq 1,$ (75)

with definite $\gamma$.

Rewriting this as

$\frac{1}{3}\frac{k_{-3}}{k_{1}}\frac{1}{C_{\mathrm{M}}a^{3}_{\mathrm{2M}}}+\frac{2}{3}\frac{k_{-3}}{k_{2}}\frac{1}{C_{\mathrm{M}}a^{3}_{\mathrm{MD}}}\leq
1,$ (76)

one has two auxiliary inequalities: $C_{\mathrm{M}}a^{3}_{\mathrm{2M}}\geq
k_{-3}/k_{1}$ and $C_{\mathrm{M}}a^{3}_{\mathrm{MD}}\geq k_{-3}/k_{2}$, which
are obviously satisfied at $k_{-3}/k_{1}\leq 1$, $k_{-3}/k_{2}\leq 1$, and
$C_{\mathrm{M}}a^{3}_{s}\geq 1$ from (55).

Now let $K_{3}=K_{-3}=0$ in (52). Denoting $\varkappa_{i}=K_{-i}/K_{i}$, we
obtain the fractions:

$\displaystyle\Omega^{\mathrm{eq}}_{\mathrm{M}}=\frac{\varkappa_{1}\varkappa_{2}}{\varkappa_{1}\varkappa_{2}+2\varkappa_{2}C_{\mathrm{M}}+3C^{2}_{\mathrm{M}}},\quad\Omega^{\mathrm{eq}}_{\mathrm{D}}=\frac{2\varkappa_{2}C_{\mathrm{M}}}{\varkappa_{1}\varkappa_{2}+2\varkappa_{2}C_{\mathrm{M}}+3C^{2}_{\mathrm{M}}},$
$\displaystyle\Omega^{\mathrm{eq}}_{\mathrm{T}}=\frac{3C^{2}_{\mathrm{M}}}{\varkappa_{1}\varkappa_{2}+2\varkappa_{2}C_{\mathrm{M}}+3C^{2}_{\mathrm{M}}},\quad
C_{\mathrm{M}}={\cal C}_{0}\left(C^{(0)}_{\mathrm{M}}\right),$ (77)

where ${\cal C}_{0}$ denotes a real solution among three roots ${\cal C}_{m}$
($m=0,\pm 1$) of a cubic equation:

$\displaystyle{\cal C}_{m}={\cal
A}\,\mathrm{e}^{\mathrm{i}\frac{2\pi}{3}m}-\frac{{\cal A}^{(2)}}{{\cal
A}}\,\mathrm{e}^{-\mathrm{i}\frac{2\pi}{3}m}-\frac{2}{9}\varkappa_{2},\quad{\cal
A}=\left[{\cal A}^{(3)}+\sqrt{\left({\cal A}^{(3)}\right)^{2}+\left({\cal
A}^{(2)}\right)^{3}}\right]^{1/3},$ (78) $\displaystyle{\cal
A}^{(2)}=\frac{\varkappa_{2}}{9}\left(\varkappa_{1}-\frac{4}{9}\varkappa_{2}\right),\quad{\cal
A}^{(3)}=\frac{\varkappa_{1}\varkappa_{2}}{6}\left[C^{(0)}_{\mathrm{M}}+\frac{2}{9}\varkappa_{2}\left(1-\frac{8}{27}\frac{\varkappa_{2}}{\varkappa_{1}}\right)\right].$
(79)

To provide $\Omega^{\mathrm{eq}}_{\mathrm{T}}\geq 0.5$, we demand
$3C^{2}_{\mathrm{M}}\geq\varkappa_{1}\varkappa_{2}+2\varkappa_{2}C_{\mathrm{M}}$
[see (71)], or, equivalently,

$C_{\mathrm{M}}\geq\frac{\varkappa_{2}}{3}\left(1+\sqrt{1+3\frac{\varkappa_{1}}{\varkappa_{2}}}\right).$
(80)

Therefore, $C_{\mathrm{M}}a^{3}_{\mathrm{MD}}>k_{-2}/k_{2}$, where
$k_{-2}/k_{2}\leq 1$ to permit the synthesis of T. Thus, condition (80) is not
stronger than (55).

Summarizing, the applied analysis enables the leading role of composites under
the derived conditions (55) and the definite values of monomer density
$C^{(0)}_{\mathrm{M}}$ and $C_{\mathrm{M}}$ [see (75) and (80)]. However, the
mechanism of the formation of composites needs detailed study.

##  5 Discussion

In this paper, we have considered three models that deal with the problem of
formation and stability of the two- and three-particle composites in dark
matter.

Starting with the Bose-condensate dark matter model with respective self-
interaction on the base of the Gross–Pitaevskii equation extended by an
additional three-particle interaction, we find that the few-particle complexes
can be formed at relatively low density and pressure, when quantum
fluctuations and the transfer momentum are small enough and cannot destroy the
condensate and the desired structures. Such a conclusion is deduced on the
base of calculating the averages of the powers of the local density $\eta$ and
thermodynamic functions at zero temperature. However, increasing the pressure
(and density), the multiparticle averages are reduced to the product of
single-particle averages which determine the dependence on thermodynamic
parameters. This enables to relate the thermodynamic functions by relations
for an ideal gas, where the mean value of quantum fluctuations is regarded as
the effective temperature of the inhomogeneous system.

Note that a similar splitting of many-particle correlators into single-
particle functions is observed in the model with pairwise interaction in the
Thomas-Fermi approximation. This emphasizes the importance of the model with
three-particle interaction along with explicit accounting of quantum
fluctuations.

We associate the revealed changes in properties with the first-order phase
transition described in our previous work [24]. Besides, as found in Ref.
[24], a description of rotational curves using this model is realized in the
vicinity of the smallest value of the chemical potential (the ground state)
and at a short pairwise scattering length $a$, the magnitude of which we
reproduced above. Then, we can compare our expression (21) for the internal
energy with the well-known one, where the term of three-particle interaction
arises due to quantum corrections in a system with pairwise interaction and
the use of renormalization group calculations [33, 41, 42, 43]:

$\displaystyle{\cal E}=\frac{{\cal Q}(n)}{2}n^{2}+\frac{{\cal B}(n)}{3}n^{3},$
(81) $\displaystyle{\cal
Q}(n)=\frac{4\pi\hbar^{2}a}{m}\left(1+\frac{128}{15\sqrt{\pi}}\sqrt{na^{3}}\right),$
(82) $\displaystyle{\cal
B}(n)=\frac{16\pi\hbar^{2}a^{4}}{m}\left(4\pi-3\sqrt{3}\right)\left[\ln{\left(na^{3}\right)}+2d\right].$
(83)

Here $n$ is the particle density of homogeneous system of bosons, $a$ is their
pairwise $s$-scattering length, $d$ is a numeric parameter.

Note that the internal energy ${\cal E}$ is proportional to inverse of the
mass of (dark matte) particle.

However, we also admit that the origin of the three-body interaction can be
explained by relativistic effects [5, 6] or by the very nature of dark matter
particles.

It is interesting to note the attempts to generalize the expression (81) in
the case of the unitary regime with an infinitely long scattering length $a$
[44, 45]. This is motivated, in particular, by the study of Bose condensates
of alkali atoms, the strong interaction between which can cause the molecules
synthesis, recombination, and three-particle losses. Hypothetically, similar
processes can result in multicomponent structure of dark matter , if the
resulting composites are stable enough and have a short scattering length.
Although, due to the greater mass, the Compton wavelength of the composite
would already become shorter than that of individual constituents.

Within the framework of an efficient model of chemical kinetics in Section 4,
we show that a long scattering length, but not a high density, is necessary
for the formation of few-particle complexes. Even without knowledge of the
nature of the bonds, it is clear that the appearance of dimers and trimers
therein cannot be compared with the synthesis of molecules from atoms that are
complicated by themselves. Moreover, there is no evidence that the processes
under consideration are similar to the production of baryons: quark-antiquark
mesons and three-quark nucleons or hyperons.

We assume that the Efimov effect could be appropriate mechanism of forming the
trimers from dark bosons. The existence of such an effect is successfully
confirmed by experiments with cold-atom systems [46] and is investigated [35].
Comparison of manifestations of this effect in the two different situations
shows a number of similar features.

To simulate unstable composites consisting of no more than three particles, we
consider a random walk model in Section 3. By adjusting the coupling
constants, it is able to obtain configurations with the dominance of
composites over long time intervals (long oscillation period). This model is
quite useful to visualize the decay of correlators in the Bose condensate with
a three-body interaction.

We expect that detailed quantum analysis may provide results applicable to the
description of dark matter, which will be published elsewhere.

## Acknowledgments

Both authors acknowledge support from the National Academy of Sciences of
Ukraine by its priority project “Fundamental properties of the matter and its
manifestation in the micro world, astrophysics and cosmology in relativistic
collisions of nuclei and in the early Universe”.

## References

  * [1] Arcadi, G. et al. The waning of the WIMP? A review of models, searches, and constraints. Eur. Phys. J. C 2018, 78, 203.
  * [2] Guth, A.H.; Hertzberg, M.P.; Prescod-Weinstein, C. Do dark matter axions form a condensate with long-range correlation? Phys. Rev. D 2015, 92, 103513.
  * [3] Valle, D.C.; Mielke, E.W. Relativistic soliton collisions of axion type dark matter. Phys. Lett. B 2016, 758, 93.
  * [4] Schwabe, B.; Niemeyer, J.C.; Engels, J.F. Simulations of solitonic core mergers in ultralight axion dark matter cosmologies. Phys. Rev. D 2016, 94, 043513.
  * [5] Chavanis, P.H. Phase transitions between dilute and dense axion stars. Phys. Rev. D 2018, 98, 023009.
  * [6] Chavanis, P.H. K-essence Lagrangians of polytropic and logotropic unified dark matter and dark energy models. arXiv: 2109.05953 2021.
  * [7] Boyarsky, A.; Drewes, M.; Lasserre, T.; Mertens, S.; Ruchayskiy, O. Sterile neutrino dark matter. Prog. Part. N. Phys. 2019, 104, 1.
  * [8] Belotsky, K.M. et al. Signatures of primordial black hole dark matter. Mod. Phys. Lett. A 2014, 29, 1440005.
  * [9] Carr, B.; Kühnel, F.; Sandstad, M. Primordial black holes as dark matter. Phys. Rev. D 2016, 94, 083504.
  * [10] Gani, V.A.; Khlopov, M.Yu.; Voskresensky, D.N. Double charged heavy constituents of dark atoms and superheavy nuclear objects. Phys. Rev. D 2019, 99, 015024.
  * [11] Hu, W.; Barkana, R.; Gruzinov, A. Cold and fuzzy dark matter. Phys. Rev. Lett. 2000, 85, 1158.
  * [12] Lee, J.W. Brief history of ultra-light scalar dark matter models. EPJ Web Conf. 2018, 168, 06005.
  * [13] Schive, H.-Yu ; Chiueh, T.; Broadhurst, T. Soliton random walk and the cluster-stripping problem in ultralight dark matter. Phys. Rev. Lett. 2020, 124, 201301.
  * [14] Matos, T.; Guzman, F.S.; Urena-Lopez, L.A. Scalar field as dark matter in the universe. Class. Quantum Grav. 2000, 17, 1707.
  * [15] Gavrilik, A.M.; Kachurik, I.I.; Khelashvili, M.V.; Nazarenko, A.V. Condensate of $\mu$-Bose gas as a model of dark matter. Physica A: Stat. Mech. Applic. 2018, 506, 835.
  * [16] Gavrilik, A.M.; Kachurik, I.I.; Khelashvili, M.V. Galaxy rotation curves in the $\mu$-deformation based approach to dark matter. Ukr. J. Phys. 2019, 64, 1042.
  * [17] Sin, S.J. Late-time phase transition and the galactic halo as a Bose liquid. Phys. Rev. D 1994, 50, 3650.
  * [18] Bohmer, C.G.; Harko, T. Can dark matter be a Bose–Einstein condensate? J. Cosmol. Astropart. Phys. 2007, 06, 025.
  * [19] Fukuyama, T.; Morikawa, M.; Tatekawa, T. Cosmic structures via Bose–Einstein condensation and its collapse. J. Cosmol. Astropart. Phys. 2008, 06, 033.
  * [20] Harko, T. Bose–Einstein condensation of dark matter solves the core/cusp problem. J. Cosmol. Astropart. Phys. 2011, 05, 022.
  * [21] Suarez, A.; Roblez, V.; Matos, T. A review on the scalar field/Bose-Einstein condensate dark matter model. Astrophys. Space Sci. Proc. 2014, 38, 107-142.
  * [22] Nazarenko, A.V. Partition function of the Bose–Einstein condensate dark matter and the modified Gross–Pitaevskii equation. Int. J. Mod. Phys. D 2020, 29, 2050018.
  * [23] Gavrilik, A.M.; Khelashvili, M.V.; Nazarenko, A.V. Bose–Einstein condensate dark matter model with three-particle interaction and two-phase structure. Phys. Rev. D 2020, 102, 083510\.
  * [24] Gavrilik, A.M.; Nazarenko, A.V. Phases of the Bose–Einstein condensate dark matter model with both two- and three-particle interactions. Universe 2021, 7, 359.
  * [25] L.A. Urena-Lopez, Scalar field dark matter with a cosh potential, revisited. J. Cosmol. Astropart. Phys. 2019, 06, 009.
  * [26] Ferreira, E.G.M. Ultra-light dark matter. Astron. Astrophys. Rev. 2021, 29, 7.
  * [27] Lifshitz, E.M.; Pitaevskii, L.P. Statistical physics: Theory of condensate state (Volume 9 of Course of Theoretical Physics); Pergamon Press: New York, USA, 1980.
  * [28] Ketterle, W.; Durfee, D.S.; Stamper-Kurn D.M. Making, probing and understanding Bose-Einstein condensates. arXiv: cond-mat/9904034 1999.
  * [29] Leggett, A.J. Bose–Einstein condensation in the alkali gases: Some fundamental aspects. Rev. Mod. Phys. 2001, 73, 307-356.
  * [30] Cabrera, C.R. et al. Quantum liquid droplets in a mixture of Bose-Einstein condensates. Science 2017, 359, 301-304.
  * [31] Roati, G. et al. ${}^{39}\mathrm{K}$ Bose-Einstein condensate with tunable interactions. Phys. Rev. Lett. 2007, 99, 010403.
  * [32] Bulgac, A. Dilute quantum droplets. Phys. Rev. Lett. 2002, 89, 050402.
  * [33] Braaten, E.; Nieto, A. Quantum corrections to the energy density of a homogeneous Bose gas. Eur. Phys. J. B 1999, 11, 143–159.
  * [34] Braaten, E.; Hammer, H.-W. Universality in few-body systems with large scattering length. Phys. Rep. 2006, 428, 259–390.
  * [35] Musolino, S.; Kurkjian, H.; Van Regemortel, M.; Wouters, M.; Kokkelmans, S.J.J.M.F.; Colussi, V.E. Bose-Einstein condensation of Efimovian triples in the unitary Bose gas. arXiv: 2106.10991 2021.
  * [36] Braaten, E.; Kang, D.; Laha, R. Production of dark-matter bound states in the early universe by three-body recombination. JHEP 2018, 11, 084.
  * [37] Poulin, A.; Godfrey, S. Multicomponent dark matter from a hidden gauged SU(3). Phys. Rev. D 2019, 99, 076008.
  * [38] Tellez-Tovar, L.O.; Matos, T.; Vazquez, J.A. Cosmological constraints on the Multi Scalar Field Dark Matter model. arXiv: 2112.09337 2021.
  * [39] Kempe, J. Quantum random walks - an introductory overview. Contemporary Physics 2003, 44, 307–327.
  * [40] Zamarev, K.I. Chemical kinetics. Part I; Novosibirsk Univ. Publ.: Novosibirsk, Russia, 2004 (in Russian).
  * [41] Lee, T.D.; Yang, C.N. Many-body problem in quantum mechanics and quantum statistical mechanics. Phys. Rev. 1957, 105, 1119.
  * [42] Wu, T.T. Ground state of a Bose system of hard spheres. Phys. Rev. 1959, 115, 1390-1404.
  * [43] Bedaque, P.F.; Bulgac, A.; Rupak, G. Quantum corrections to dilute Bose liquids. Phys. Rev. A 2003, 68, 033606.
  * [44] Cowell, S. et al. Cold Bose gases with large scattering lengths. Phys. Rev. Lett. 2002, 88, 210403.
  * [45] Castin, Y.; Werner, F. Single-particle momentum distribution of an Efimov trimer. Phys. Rev. A 2011, 83, 063614.
  * [46] Ferlaino, F.; Grimm, R. Forty years of Efimov physics: How a bizarre prediction turned into a hot topic. Physics 2010, 3, 9.