# Nonrelativistic Treatment of Fully-Heavy Tetraquarks as Diquark-Antidiquark
States

Halil Mutuk<EMAIL_ADDRESS>Department of Physics, Faculty of Science and
Letters, Ondokuz Mayis University, 55139, Samsun, Turkey

###### Abstract

The goal of the present work is to obtain a reliable estimate of the masses of
the ground and radially excited states of fully-heavy tetraquark systems. In
order to do this, we use a nonrelativistic model of tetraquarks which are
assumed to be compact and consist of diquark-antidiquark pairs. This
nonrelativistic model is composed of Hulthen potential, a linear confining
potential and spin-spin interaction. We computed ground, first, and second
radially excited $cc\bar{c}\bar{c}$ and $bb\bar{b}\bar{b}$ tetraquark masses.
It was found that predicted masses of ground states of $cc\bar{c}\bar{c}$ and
$bb\bar{b}\bar{b}$ tetraquarks are significantly higher than the thresholds of
the fall-apart decays to the lowest allowed two-meson states. These states
should be broad and are thus difficult to observe experimentally. First
radially excited states are considerably lower than their corresponding
(2S-2S) two-meson thresholds. We hope that our study may be helpful to the
experimental search for ground and excited $cc\bar{c}\bar{c}$ and
$bb\bar{b}\bar{b}$ tetraquark states.

##  I Introduction

Quark model describes ordinary mesons as ($q\bar{q}$) systems and baryons as
($qqq$) systems in terms of quarks $q$ and antiquarks $\bar{q}$. In addition
to quark model, the existence of multiquark states such as tetraquarks
($q\bar{q}q\bar{q}$ or $qq\bar{q}\bar{q}$), pentaquarks ($q\bar{q}qqq$), and
structures with more quarks was proposed decades ago [1, 2, 3, 4]. These new
structures present quantum numbers, masses, flavours, decay channels and
widths in the experiments which cannot be fitted into conventional quark
model. Now they are called exotic hadrons. Concerning tetraquarks, the first
discovery of exotic states was made in 2003 by Belle Collaboration in the
charmonium sector [5]. This state was named as $X(3872)$ (now referred to name
as $\chi_{c1}(3872)$) with a mass of ($3872.0\pm 0.6$ MeV). Its extraordinary
properties, such as isospin violation [6] and radiative decays [7, 8] still
make its structure unsolved. Many new exotic hadron candidates have been
observed by decaying into final states of a charm and anticharm quarks. These
states are called $XYZ$ states. The $X$ states are neutral and have positive
parity quantum number and generally seen in $J/\psi+\text{pions}$ decays.
Typical example is the observation of $\chi_{c1}(3872)$. The $Y$ states are
neutral and have negative parity quantum number, and seen in $e^{+}e^{-}$
annihilation with or without initial state radiation. $Y(4260)$ is an example
of this family [9]. The $Z$ states are mostly charged or neutral, have
typically positive parity and decay into $J/\psi+\pi$, $h_{c}(1P)+\pi$,
$\chi_{c}(1P)+\pi$. The charged $Z_{c}(3900)$ state is the famous example of
$Z$ family [10, 11]. Observation of pentaquark states in 2019 made exotic
hadrons more interesting [12].

The physics of exotic hadrons is a thorough piece of research which involves
both short and long distance behaviors of QCD. At one side, increasing values
of radius $r$ are needed when considering the spectroscopy. This is the place
where nonperturbative effects are at the stage. On the other side hard
processes, such as decays, occur at short distances, i.e., short values of
radius $r$. This is the place where perturbative effects are at the stage.
There are two different regimes which make theoretical predictions difficult.
Many theoretical and phenomenological models are being studied to understand
and interpret these exotic states such as lattice QCD, dynamically generated
resonances, QCD sum rules, coupled channel effects and nonrelativistic
effective field theories (see Ref. [13] and references therein).

Generally two pictures are taken into account in the approaches mentioned
above: molecular and compact tetraquark pictures. Hadronic molecules are
loosely bounded systems together by the exchange of pions and other light
mesons. This scenario has received a lot of interest due to the masses of
several $XYZ$ hadrons are very close to the related meson-antimeson
thresholds. In the case of $\chi_{c1}(3872)$ it has been suggested that if it
has a binding energy of less than 200 keV with respect to the
$D^{0}\bar{D}^{\ast 0}$, according to the $R=1/\sqrt{2\mu E_{B}}$, where $\mu$
is the reduced mass of the two-hadron system and $E_{B}$ is the binding
energy, it would be at least as large as 10 fm [14]. From this perspective,
hadronic molecules can be seen as extended objects. Compact tetraquarks are
bound states of color nonsinglet diquark-antidiquarks, tightly bound by
gluons, very much along the same lines as colored quark-antiquark pairs are
bound into color-neutral mesons [15]. A diquark is a bound quark-quark $(qq)$
pair, whereas an antidiquark is a bound antiquark-antiquark $(\bar{q}\bar{q})$
pair. These pairs are colored, i.e., have non-zero color charges and can have
colorless combinations which turns out to be an ansatz of tetraquark paradigm.
Tetraquark configurations are not ruled out by QCD. Indeed this context opens
a new window of compact hadrons in QCD, which are even more numeorus than the
conventional quark-antiquark mesons.

A very recent study of LHCb present a $J/\psi$-pair invariant mass spectrum by
using pp collision data [16]. A narrow structure around 6.9 $\text{GeV/c}^{2}$
matching the lineshape of a resonance and a broad structure just above twice
the $J/\psi$ mass are observed. The deviation of the data from nonresonant
$J/\psi$-pair production is above five standard deviations in the mass region
between 6.2 and 7.4 $\text{GeV/c}^{2}$, covering predicted masses of states
composed of four charm quarks [16]. This energy range lie well above the
experimentally known range for charmonium which is in the range of 3 - 4.5 GeV
and makes fully-charm tetraquark state $cc\bar{c}\bar{c}$ very interesting.
The reason for this is that, the energy range of the $XYZ$ states are in the
same mass range of the conventional charmonium states and this can yield a
confusion on these structures [17]. The observation of a possible
$cc\bar{c}\bar{c}$ structure triggered many theoretical studies [18, 19, 20,
21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37].

In the side of $bb\bar{b}\bar{b}$ structures, there is no observation up to
now. An experimental study claimed the existence of a full-bottom tetraquark
states $bb\bar{b}\bar{b}$, with a global significance of 3.6 $\sigma$ and a
mass around 18.4 GeV, almost 500 MeV below the threshold of $\Upsilon\Upsilon$
[38]. However, LHCb Collaboration presented an intriguing analysis looking for
the exotic $bb\bar{b}\bar{b}$ tetraquark in the $\Upsilon(1S)\mu^{+}\mu^{-}$
final state and announced that no observation was made [39]. Prior to and
after these experimental studies, the possible existence of $bb\bar{b}\bar{b}$
state was investigated on the theoretical basis [40, 41, 42, 43, 44, 45, 46,
47, 48].

Motivated with this charm sector observation and a possible open window for
fully-bottomed tetraquark states, in the present work we will use a
nonrelativistic model to study tetraquarks as composed of diquarks and
antidiquarks, which interact much like ordinary quarkonia. Owing to the masses
of the valence degree of freedoms, fully-heavy four quark states can be
investigated by nonrelativistic approach. We calculate mass spectra of fully-
heavy tetraquark systems assuming that they are composed of doubly heavy quark
($QQ$) and antidiquark ($\bar{Q}\bar{Q}$). Choosing this configuration is not
just an assumption. Firstly, in light quark systems, the binding mechanism is
maintained by the light-meson (such as pion, $\pi$) and gluon exchange. The
binding mechanism in fully-heavy systems is probably dominated by the gluon-
exchange forces since the typical gluon mass scale is $m_{g}\sim
0.5~{}\text{GeV}$. This mass value is much lighter than the possible force
carriers of heavy-mesons that could be exchanged between the heavy diquark
($QQ$) and antidiquark ($\bar{Q}\bar{Q}$) of tetraquark structure. Secondly,
there is some evidence of diquark clustering in baryons. In 2017, the LHCb
Collaboration reported the observation of a doubly charmed and doubly charged
baryon, $\Xi_{cc}^{++}$ a $ccu$ state where the charm diquark may play a role
in the structure, has lead to further attention on heavy-quark systems as the
description of exotic hadrons [49].

This paper is organized as follows: In Section II, we describe nonrelativistic
potential model of this work. In Section III, the masses of
diquark/antidiquark and tetraquark systems are calculated. Detailed
comparisons of diquark and tetraquark masses with previous studies within
different approaches are given. Section IV is reserved for conclusion and
summary of the obtained results.

##  II Potential Model

The fundamental assumption in the quark potential model is that, if one
integrate out the gluon fields in the QCD action, he/she can hopefully obtain
an effective Hamiltonian with suitable potentials which describes the physics
of hadrons fairly good. The quark masses in the resulting effective
Hamiltonian might not be the same quark masses in the QCD Lagrangian, that is
why they are named as constituent quark masses.

The use of quark potential models to describe the energy spectra of mesonic
and baryonic systems gave reliable results. Phenomenological models with a
simple relativistic kinetic energy term plus a scalar potential term which
incorporates the so called linear confinement plus a term related to short
distance which incorporates color-Coulomb interaction stemmed from QCD, give
good results and descriptions of the observed spectra of both heavy and light
quark mesons and baryons.

In the present work, we have considered the following potential:

$V(r)=V_{H}(r)+V_{C}(r)+V_{0}.$ (1)

Here, $V_{H}$ denotes the Hulthen potential:

$V_{H}(r)=-\frac{h}{\text{exp}(br)-1},~{}b>0,~{}h\geq 0,$ (2)

where $h$ is the coupling strength and $b$ is the range. Hulthen potential is
extensively used in atomic, molecular, plasma and mathematical physics [50,
51, 52, 53, 54, 55], chemical physics [56, 57, 58], nuclear physics [59, 60,
61], and particle physics [62, 63, 64].

In the limit $b\to 0$, the Hulthen potential approaches from above the
Coulomb-like potential

$V_{H}(r)=-\frac{h}{br}\approx V_{\text{Coulomb}}(r)=-\frac{\kappa}{r}.$ (3)

In other words, the Hulthen potential behaves like the Coulomb potential as
$r\to 0$ but decreases exponentially in the asymptotic region when $r$ is
sufficiently large. $V_{C}$ is the confining part of the potential with

$V_{C}(r)=cr,$ (4)

where $c$ is a constant. $V_{0}$ is also a constant which would act as a zero-
point energy. The parameters for this potential are given in the Table 1 [65]:

Table 1: Potential model parameters. Parameter | Numerical value
---|---
$h$ | $0.20~{}\text{GeV}$
$b$ | $0.4~{}\text{GeV}$
$c$ | $0.193~{}\text{GeV}^{2}$
$V_{0}$ | $-0.223~{}\text{GeV}$
$m_{c}$ | $1.4~{}\text{GeV}$
$m_{b}$ | $4.812~{}\text{GeV}$
$\alpha_{s}(m_{c}^{2})$ | 0.37
$\alpha_{s}(m_{b}^{2})$ | 0.26

To see whether or not $V(r)$ (Eq. 1) coincides with the properties of Cornell
potential exhibit, we can plot them. At first step, it can be shown that for
large $r$, the Hulthen potential approaches to zero faster than the Coulomb-
like potential, $V(r)=-\frac{4}{3}\frac{\alpha_{s}}{r}$ with
$\alpha_{s}=0.5202$ [17]. This can be seen in Figure 1.

Figure 1:  Comparison of Coulomb potential and Hulthen potential for large
$r$. $V(r)$ is in GeV and $r$ in fm units.

The choice of Hulthen potential could be explained as follows. As can be seen
from Figure 1, the contribution of large $r$ value of Hulthen potential is
less than the Coulomb potential as $r$ increases. This is due to the nature of
short range potentials. Furthermore, the general theory of scattering is not
immediately applicable to the Coulomb potential case because it decreases too
slowly as the distance increases, which can be seen in Figure 1. The Coulomb
potential falls off so slowly that it continues to influence the particles
even as they move apart [66]. In order to apply the general theory of
scattering to the Coulomb potential case, some modifications of the Coulomb
potential have been considered. For example replacing Coulomb potential by the
Yukawa potential and subsequently making the exponential factor approach unity
is one example of this modification. Hulthen potential is a good approximation
of Yukawa potential, $V(r)=-V_{0}\frac{e^{-\alpha r}}{r}$. Furthermore,
Hulthen potential combined with a linear term, $\sim r^{n}$, can be a
theoretical playground for the different forms of linear part of full
potentials.

In addition to above discussion, the spectrum generated by the Hulthen
potential could be investigated. Since this is done in literature excessively,
it will be good to mention the well-known operator inequality. It was shown in
[67] that

$V_{\text{Coulomb}}(r)=-\frac{\kappa}{r}\leq-\frac{h}{br}\leq-\frac{h}{\text{exp}(br)-1}=V_{H}(r),$
(5)

for $\frac{h}{b}\leq\kappa$. This operator inequality between the Coulomb and
Hulthen potentials is shared by the corresponding Hamiltonians and thus
carries over to the entailed (ground state) energy eigenvalues.

On the other side, the plots of Cornell potential
$V(r)=-\frac{4}{3}\frac{\alpha_{s}}{r}+br$ with $\alpha_{s}=0.5202$ and
$b=0.1463~{}\text{GeV}^{2}$ [17], and potential model of this work are given
in Figure 2.

Figure 2:  Comparison of Cornell potential and potential model of this work.
$V(r)$ is in GeV and $r$ in fm units.

It can be seen from Figure 2 that in the range of $0<r<1~{}\text{fm}$, both
potentials have similar behaviour. In literature, all the phenomenological
potentials have almost similar behaviour in the range of confinement,
$r<1~{}\text{fm}$. This is the characteristic interval of heavy quarkonium
systems, charmonium $c\bar{c}$ and bottomonium $b\bar{b}$. At its present
status, it is not possible to obtain an exact interquark potential in the
whole range of distances from the first principle of QCD.

A nonrelativistic approach with static potentials can be reasonable under the
condition that the kinetic energy is much less than the rest masses of the
constituents, which is usually the case when considering heavy quark bound
states. Another advantage of assuming tetraquarks are composed from the heavy
diquark ($QQ$) and heavy antidiquark ($\bar{Q}\bar{Q}$) is that it
significantly simplifies calculations, since the four-body problem of a four
quark structure is very complicated. In the context of diquark-antidiquark
configuration, the problem turns out to be a two two-body problem. In a two-
body problem which involves a central potential, it is legitimate to work in
the center-of-mass frame (CM), where one can use spherical coordinates to
separate the radial and angular parts of the wave function. In terms of the
reduced mass, let $\mu\equiv m_{1}m_{2}/(m_{1}+m_{2})$ where $m_{1}$ and
$m_{2}$ are the constituent masses of quark 1 and quark 2, respectively. Since
we are dealing equal masses of quarks, $m\equiv m_{1}=m_{2}$ yields reduced
mass as $\mu=m/2$. Hence, the time-independent radial Schrödinger equation can
be written as

$\left[\frac{1}{2\mu}\left(-\frac{d^{2}}{dr^{2}}+\frac{L(L+1)}{r^{2}}\right)+V(r)\right]y(r)=Ey(r),$
(6)

with the orbital quantum number $L$ and the energy eigenvalue $E$. In order to
account full spectra and splittings in energy between states with different
quantum numbers, one can include spin-dependent terms to the potential. Due to
the two nature of energies (high distances and short distances), it is not
clear how the spin dependent forces change with distance. Assuming that spin
dependent forces stem from short distance potential, observed hadron masses
can be explained within well accuracy. For more about the nature of these spin
dependent forces see Refs. [68, 69]. The spin-dependent parts of the potential
are sensitive to the Lorentz structure of the interquark potential. Hence, one
has to split the total potential as the sum of Lorentz vector $V_{v}(r)$ and
Lorentz scalar $V_{s}(r)$. Based on the Breit–Fermi interaction for one-gluon-
exchange [70, 71, 72], the following spin dependent terms can be added to the
potential. For equal masses $m=m_{1}=m_{2}$, with $m$ being the constituent
mass of the two-body problem (charm quark, bottom quark, or diquark) spin-spin
interaction can be defined as

$V_{SS}=\frac{2}{3m^{2}}\nabla^{2}V_{v}(r)\textbf{S}_{1}\cdot\textbf{S}_{2}.$
(7)

This interaction is related to the $S-$wave splittings. This term has no
effect for $\ell\neq 0$ states. The expectation value of the operator for the
spin-spin interaction can be calculated in terms of the spin quantum numbers
by using

$\langle\textbf{S}_{1}\cdot\textbf{S}_{2}\rangle=\langle\frac{1}{2}(S^{2}-S_{1}^{2}-S_{2}^{2})\rangle$
(8)

where $S$, $S_{1}$, and $S_{2}$ is the total spin, the spin of quark 1, and
the spin of quark 2, respectively.

The contribution of spin-dependent terms can be calculated by writing total
potential Eq. (1) as sum of Lorentz vector structure and Lorentz scalar
structure:

$V(r)=V_{v}(r)+V_{s}(r).$ (9)

The Lorentz structures are well known, for example, of Cornell potential but
it is not clear in the Eq. (1). To account these Lorentz structures, it was
assumed that the terms in Eq. (1) have partly vector and partly scalar [73].
In our formalism, as an alternative we use a different manner. The spin-spin
correction to the nonrelativistic potential can be obtained as [74]:

$\displaystyle\langle V_{SS}\rangle$ $\displaystyle=$
$\displaystyle\frac{8\pi}{3}\frac{\kappa_{s}\alpha_{s}}{m^{2}}\int\Psi^{\ast}(\textbf{r})\Psi(\textbf{r})\delta(\textbf{r})\langle\textbf{S}_{1}\cdot\textbf{S}_{2}\rangle
d^{3}r$ (10) $\displaystyle=$
$\displaystyle\frac{8\pi}{3}\frac{\kappa_{s}\alpha_{s}}{m^{2}}|\Psi(0)|^{2}\langle\textbf{S}_{1}\cdot\textbf{S}_{2}\rangle.$

If the spin-spin interaction was treated as a first-order perturbation without
the Gaussian smearing, as it is the case in this work, it would be
proportional to modulus of the wavefunction at the origin, $|\Psi(0)|^{2}$.
Since spin-spin interaction only occurs in $S-$wave, only $S-$wave states
(i.e., orbital angular momentum $\ell=0$) have non-zero value of the
wavefunction at the origin. Therefore for $S-$wave state we have [70]:

$|\Psi(0)|^{2}=|Y_{0}^{0}(\theta,\phi)R_{n,\ell}(0)|^{2}=\frac{|R_{n,\ell}(0)|^{2}}{4\pi}.$
(11)

$|R_{n,\ell}(0)|^{2}$ can be obtained directly from the numerical calculations
and is related to the radial potential as

$|\Psi(0)|^{2}=\frac{\mu}{2\pi}\langle\frac{d}{dr}V(r)\rangle\Rightarrow|R_{n,\ell}(0)|^{2}=2\mu\langle\frac{d}{dr}V(r)\rangle,$
(12)

where $\mu$ is the reduced mass, $\mu=m_{1}m_{2}/(m_{1}+m_{2})$. Then it is
possible to replace $|\Psi(0)|^{2}$ in Eq. (10) with $|R_{n,\ell}(0)|^{2}$.

In QCD, the coupling constant $\alpha_{s}$ is not actually a “constant”. It
changes according to the energy scale of each bound state. Therefore it
changes with respect to energy scale and that is why it is called “running”
coupling constant. We have used the nonrelativistic model of Ref. [65] where
$\alpha_{s}$ is constant and has different values for fitting quarkonium
spectra. Choosing $\alpha_{s}$ constant is a common approach in many of the
nonrelativistic quark potential models.

By solving Schrödinger equation, one can obtain the energy eigenvalue $E$. The
result will depend on the number of nodes of the wave function $n$ (or
principal quantum number $N=n+1$) and the orbital angular momentum number
$\ell$. The mass calculated in this case is called spin-averaged mass. If one
include spin-dependent terms in the potential, the result will also depend on
the total spin $S$ of the constituents spins $S_{1}$ and $S_{2}$. The
Schrödinger equation has no analytical solution for the potential given in
this work. So we solve it numerically inspired from Ref. [75].

After having solved Schrödinger equation, the mass $M$ of a particular state
can be written as

$M=2m_{q/d}+E+\langle V_{ss}\rangle,$ (13)

where $m_{q}$ is the mass of corresponding quark/antiquark and $m_{d}$ is the
mass of corresponding diquark/antidiquark. We have considered the diquarks and
antidiquarks as constituents of the tetraquarks to predict the tetraquark
masses. Therefore, we first calculate the diquark masses and then tetraquark
masses.

##  III Numerical Analysis and Discussion

To study tetraquark systems in terms of diquark states, we will elaborate
four-particle system as two two-body systems, i.e., diquark-antidiquark
system. We assumed that the interaction between the diquarks and antidiquarks
is to be effectively the same for ordinary quarkonia. In this picture, a
diquark and antidiquark interact as a whole which means interactions between
quarks from a diquark system ($QQ$) with antiquarks from an antidiquark system
($\bar{Q}\bar{Q}$) are not considered. Since quark and antiquark masses are
same one does not have to calculate antidiquark masses separately. We first
calculate diquark masses (antidiquark masses at the same time) and use these
diquark (antidiquark) masses to calculate tetraquark masses in diquark-
antidiquark picture. The motivation for this factorization is the color
structure of diquarks.

Hadrons can exist only when their total color charge of constituent quarks are
zero. Technically, this means that every naturally occuring hadron is a color
singlet under the group symmetry SU(3). A diquark is composed of two quarks
($qq$) whereas conventional quark-antiquark (also called quarkonium) is
composed of ($q\bar{q}$). The difference between quark-quark systems and
quark-antiquark systems is mainly due to the color structure.

Regarding the Cornell potential, $V(r)=\kappa\frac{\alpha_{s}}{r}+br$, where
$\kappa$ is the color factor, $\alpha_{s}$ is the QCD fine structure constant
and $b$ is string tension and related to the strength of the confinement, the
color factor change in the diquark configuration. The SU(3) color symmetry of
QCD implies that, when we combine a quark and an antiquark in the fundamental
color representation, we obtain
$|q\bar{q}\rangle:\textbf{3}\bigotimes\bar{\textbf{3}}=\textbf{1}\bigoplus\textbf{8}$.
This representation gives the color factor for the color singlet as
$\kappa=-4/3$ of the quark-antiquark system. When we combine two quarks in the
fundamental color representation, it reduces to
$|qq\rangle:\textbf{3}\bigotimes\textbf{3}=\bar{\textbf{3}}\bigoplus\textbf{6}$,
a color antitriplet $\bar{\textbf{3}}$ and a color sextet 6. In a similar way,
combining two antiquarks reduces to
$|\bar{q}\bar{q}\rangle:\bar{\textbf{3}}\bigotimes\bar{\textbf{3}}=\textbf{3}\bigoplus\bar{\textbf{6}}$,
a triplet 3 and antisextet $\bar{\textbf{6}}$. Accordingly, combining an
antitriplet diquark and a triplet antidiquark reduce to
$|\left[qq\right]-\left[\bar{q}\bar{q}\right]\rangle:\textbf{3}\bigotimes\bar{\textbf{3}}=\textbf{1}\bigoplus\textbf{8}$,
and form a color singlet for which the one-gluon exchange potential is
attractive. The antitriplet state has a color factor $\kappa=-2/3$ which is
attractive whereas the sextet state has a color factor $\kappa=+1/3$ which is
repulsive. Therefore we will only consider diquarks in the antitriplet color
state. This conclusion was achieved, for example by Ref. [76], in which it is
shown that for single-flavor tetraquarks, only the antitriplet diquarks can
build pure states.

This difference in color structure of the quark-antiquark and quark-quark
systems make possible to extend the quark-antiquark model to be valid in the
model of quark-quark system, just by changing the color factor $\kappa$ and
the string tension $b$. Changing color factor $\kappa=-4/3$ (for quark-
antiquark system in color singlet state) to $\kappa=-2/3$ (quark-quark system
in the antitriplet color state) is equivalent of introducing a factor of 1/2
in the Coulomb part of the Cornell potential for the conventional quark-
antiquark system. This factor should be taken as a global factor since it
comes from the color structure of the wave function. Therefore the string
tension should also be divided by a factor of 2. The general rule for diquark
potential from quark-antiquark potential is making $V_{qq}=V_{q\bar{q}}/2$.
This conclusion was done in different tetraquark models [17, 77, 78, 79]. So
we also divide our potential (Eq. 1) by a factor of 2 for obtaining diquark
spectra.

###  III.1  Diquarks

We now present our numerical calculations for masses of diquarks composed of
charm ($c$) and bottom ($b$) quarks, which are also equivalent for antidiquark
cases in this work. We assume that information in the spin-dependent
interaction described in Eq. (7) is inherited when changing from quark-
antiquark system to quark-quark system. In other words, spin-spin interaction
is encoded in the diquarks. Similar assumptions were made in [17, 79, 80] and
gave reliable results. In the present work, we choose the attractive color
antitriplet state which is antisymmetric in the color wavefunction. In order
to maintain Pauli exclusion principle, the diquark total spin must be 1,
$S=1$. So the total wavefunction of the diquark will be antisymmetric.

Based on the previous arguments, the results for the ground and first radially
excited diquark masses are presented in Table 2.

Table 2: Results for the four diquark masses for $cc$ and $bb$. Corresponding diquark masses are same as $cc$ and $bb$, respectively. All results are in MeV unit.  Diquark | $N^{2S+1}L_{J}$ | Mass
---|---|---
$cc$ |  |
| $1^{3}S_{1}$ | 3114
| $2^{3}S_{1}$ | 3443
$bb$ |  |
| $1^{3}S_{1}$ | 9792
| $2^{3}S_{1}$ | 10011

Comparison of ground and first radially excited states for the $cc$ diquark
masses of earlier works are shown in Table 3.

Table 3: Comparison of ground and first radially excited states for $cc$ diquark masses. All results are in GeV unit.  State | This work | [17] | [81]
---|---|---|---
$1S$ | 3.114 | 3.133 | 3.13
$2S$ | 3.443 | 3.456 | 3.47

As can be seen from Table 3, masses for ground and first radially excited $cc$
diquark states are in good agreement. We note that Ref. [17] used a Cornell
potential with nonrelativistic framework and Ref. [81] used a nonrelativistic
QCD motivated potential. We also compare ground states masses of $cc$ diquark
with the available theoretical studies in Table 4.

Table 4: Comparison of $cc$ diquark masses of this work and results from other works. All diquarks are considered to be in the ground state $N^{2S+1}L_{J}=1^{3}S_{1}$. All results are in MeV. Reference | Mass
---|---
This work | 3114
[17] | 3133
[26] | 3226
[41] | 3130
[43] | 3204
[79] | 3128
[81] | 3130
[82] | 3329
[83] | 3144
[84] | 3510 $\pm$ 350
[85] [Mod. Ia] | 3400
[85] [Mod. Ib] | 3370
[85] [Mod. IIa] | 3420
[85] [Mod. IIb] | 3370

As can be seen from Table 4, $cc$ diquark masses lie in the range of
$3.1-3.5~{}\text{GeV}$. Differences are due to the models of the references.
In the present work, diquark masses depend on the parameters of the potential
model, which were given in Section II. Compared with the values of the
reference works, the results deviate with at most 400 MeV. Our result is in
good agreement with the results in Refs. [17, 41, 79, 81] which used
nonrelativistic quark model. Ref. [26] used diquark-antidiquark picture in the
framework of the relativistic quark model based on the quasipotential
approach. In their model, quark-quark and diquark-antidiquark interactions in
the potential are constructed similar to the mesons and baryons. Their $cc$
diquark mass is approximately 100 MeV heavier than our result. In Ref. [43],
diquark masses are calculated by a framework based on existing empirical
information about mesons and baryons. The predicted value is approximately 100
MeV higher than our result. In Ref. [82], diquark masses are calculated by
taking into account spin-spin, spin-orbit, and tensor interactions in quark
model. The predicted value of Ref. [82] is higher by about 200 MeV compared to
our result. In Ref. [83], Schwinger-Dyson and Bethe-Salpeter equations are
solved in order to obtain diquark masses by $\text{vector}\times\text{vector}$
interaction model (denoted as CI in reference paper) in which the obtained
mass agree well compared to our result. A QCD Sum Rule study was done in Ref.
[84] yielding a mass value which is consistently bigger by about 400 MeV than
our result. Ref. [85] made an exploratory study with doubly heavy baryon which
is composed of a heavy diquark and a light quark and calculated diquark masses
in the Bethe-Salpeter (BS) formalism. They used two different parameters sets
with confining parameter $\kappa^{\prime}$ and coupling strength $\alpha_{s}$.
They also took care of heavy quark limit. Mod. Ia refers $m_{Q}\to\infty$,
i.e., taking the heavy quark limit while Mod. Ib refers
$m_{Q}\to\text{finite}$, i.e., without heavy quark limit with appropriate
values of $\kappa^{\prime}$ and $\alpha_{s}$. The heavy quark limits are the
same but parameters are different for Mod. IIa and Mod. IIb. Our results are
at most 300 MeV lower than the results of Ref. [85].

In Table 5, we compare $bb$ diquark masses of ground and first radially
excited states of different works.

Table 5: Comparison of ground and first radially excited states for $bb$ diquark. All results are in GeV unit.  State | This work | [81] | [86]
---|---|---|---
$1S$ | 9.792 | 9.72 | 9.74
$2S$ | 10.011 | 10.01 | 10.02

It can be seen from Table 5 that ground and first excited $bb$ diquark masses
are in good agreement. We note that Ref. [86] used the QCD potential of
Buchmüller and Tye under nonrelativistic quark model framework.

As happened before in the $cc$ diquark case, we compare ground state masses of
$bb$ diquark with the available theoretical studies in Table 6.

Table 6: Comparison of $bb$ diquark masses of this work and results from other works. All diquarks are considered to be in the ground state $N^{2S+1}L_{J}=1^{3}S_{1}$. All results are in MeV. Reference | Mass
---|---
This work | 9792
[26] | 9778
[41] | 9720
[43] | 9719
[79] | 9643
[81] | 9720
[82] | 9845
[83] | 9491
[84] | 8670 $\pm$ 690
[85] [Mod. Ia] | 10070
[85] [Mod. Ib] | 10050
[85] [Mod. IIa] | 10080
[85] [Mod. IIb] | 10050
[86] | 9740
[87] | 9850

The $bb$ diquark masses lie in the range of $8.6-10.0~{}\text{GeV}$. This gap
is bigger with respect to $cc$ diquark mass which was around 400 MeV. The
reason for this could be the constituent mass of $b$ quark, which is roughly
three or four times $c$ quark mass in constituent quark models.
Nonrelativistic quark model results agree well with our result [17, 41, 79,
81, 86]. Ref. [87] calculated $bb\bar{b}\bar{b}$ tetraquark mass in a
nonrelativistic effective field theory and in a relativized diquark model
characterized by one-gluon-exchange (OGE) plus a confining potential. In order
to do this, they estimated $bb$ diquark mass by binding the OGE plus a
confining potential. Their result agree well with our result.

The differences of our results in $cc$ and $bb$ diquark masses apart from
nonrelativistic formalisms could be as a result of the relativistic effects,
spin-orbit and tensor interactions, model parameters or method uncertainties
(there exist a typical 10% uncertainty in the results of QCD Sum Rules).

###  III.2 Tetraquarks

In tetraquark systems, we treated fully-heavy tetraquark as a two body
($QQ-\bar{Q}\bar{Q}$) system with equal masses, $m_{QQ}=m_{\bar{Q}\bar{Q}}$.
We also treated them as axial-vector states composed of diquarks and
antidiquarks and we assumed that the interaction between the diquarks and
antidiquarks mimickes the ordinary quarkonia interaction. It should be also
noted that, tetraquarks made up of quarks of the same flavor can only be made
up of axial-vector diquark-antidiquark pairs because of Pauli principle.

One of the important advantages of quark model is that, it allows to study
systems with various combinations of principal quantum number, orbital angular
momentum, spin and total angular momentum. Coupling of spin 1 diquark and spin
1 antidiquark produces tetraquark of total spin 0, 1 and 2 which we denote as
$S_{T}=0,1,2$. We transferred quantum mechanic couplings into this
nonrelativistic approach where total spin $S_{T}$ and orbital angular momentum
$L_{T}$ couple into the total angular momentum $J_{T}$. Obtaining charge and
parity quantum numbers of tetraquark states is quite different than the
conventional method. We obtain charge and parity quantum numbers of
tetraquarks as explained in [17, 80, 88]. We will use the following notation

$|T_{QQ\bar{Q}\bar{Q}}\rangle=|S_{d},S_{\bar{d}},S_{T},L_{T},J_{T}\rangle,$
(14)

where $S_{d}$ is the total spin of the diquark, $S_{\bar{d}}$ is the total
spin of the antidiquark, $S_{T}$ is the total spin of the tetraquark, $L_{T}$
is the orbital angular momentum relative to the diquark-antidiquark system (in
the two-body scheme), $J_{T}$ is the total angular momentum of the tetraquark
comes from the coupling of $S_{T}$ and $L_{T}$. Then, charge and parity
quantum numbers can be written as

$\displaystyle C_{T}$ $\displaystyle=$ $\displaystyle(-1)^{L_{T}+S_{T}},$ (15)
$\displaystyle P_{T}$ $\displaystyle=$ $\displaystyle(-1)^{L_{T}}.$ (16)

Since we use diquarks/antidiquarks with spin 1 in the antitriplet color
configuration, the resulting tetraquark states in $S-$wave can have the
following possibilities:

$\displaystyle|0^{++}\rangle$ $\displaystyle=$
$\displaystyle|S_{d}=1,S_{\bar{d}}=1,S_{T}=0,L_{T}=0,J_{T}=0\rangle,$ (17)
$\displaystyle|1^{+-}\rangle$ $\displaystyle=$
$\displaystyle|S_{d}=1,S_{\bar{d}}=1,S_{T}=1,L_{T}=0,J_{T}=1\rangle,$ (18)
$\displaystyle|2^{++}\rangle$ $\displaystyle=$
$\displaystyle|S_{d}=1,S_{\bar{d}}=1,S_{T}=2,L_{T}=0,J_{T}=2\rangle.$ (19)

Considering diquarks and antidiquarks as consituents of tetraquarks, masses of
the tetraquarks can be calculated using the diquark masses with
$N^{2S+1}L_{J}=1^{3}S$ in Table 2 for $cc$ and $bb$ diquarks. We also show the
values of $\Delta=M_{\text{tetra}}-M_{\text{threshold}}$, where
$M_{\text{tetra}}$ is the tetraquark mass and $M_{\text{threshold}}$ is the
mass of its lowest meson-meson threshold. A negative $\Delta$ means that the
tetraquark lies below the threshold of the fall-apart decay into two mesons
and thus should be a narrow state. Besides that a state with small positive
$\Delta$ value could also be observed as a resonance since its decay rate will
be suppressed by the phase space. All other states with large positive
$\Delta$ values are expected to be broad and difficult to observe in the
experiments. The results are presented in Table 7.

Table 7: Masses of ground, first, and second radially excited fully-heavy $cc\bar{c}\bar{c}$ and $bb\bar{b}\bar{b}$ tetraquark states. $E_{th}$ is the threshold mass of two heavy $Q\bar{Q}$ mesons. All results are given in MeV. Configuration | $N^{2S_{T}+1}L_{T_{J_{T}}}$ | $J^{PC}$ | Mass | Threshold | $E_{th}$ | $\Delta$
---|---|---|---|---|---|---
$cc\bar{c}\bar{c}$ | $1^{1}S_{0}$ | $0^{++}$ | 6322 | $\eta_{c}(1S)\eta_{c}(1S)$ | 5968 | 354
| $1^{1}S_{0}$ | $0^{++}$ | 6322 | $J/\psi(1S)J/\psi(1S)$ | 6194 | 128
| $1^{3}S_{1}$ | $1^{+-}$ | 6354 | $\eta_{c}(1S)J/\psi(1S)$ | 6081 | 273
| $1^{5}S_{2}$ | $2^{++}$ | 6385 | $J/\psi(1S)J/\psi(1S)$ | 6194 | 191
| $2^{1}S_{0}$ | $0^{++}$ | 6575 | $\eta_{c}(1S)\eta_{c}(1S)$ | 5968 | 607
| $2^{1}S_{0}$ | $0^{++}$ | 6575 | $J/\psi(1S)J/\psi(1S)$ | 6194 | 381
| $2^{3}S_{1}$ | $1^{+-}$ | 6609 | $\eta_{c}(1S)J/\psi(1S)$ | 6081 | 528
| $2^{5}S_{2}$ | $2^{++}$ | 6639 | $J/\psi(1S)J/\psi(1S)$ | 6194 | 445
| $2^{1}S_{0}$ | $0^{++}$ | 6575 | $\eta_{c}(2S)\eta_{c}(2S)$ | 7276 | -701
| $2^{1}S_{0}$ | $0^{++}$ | 6575 | $\psi(2S)\psi(2S)$ | 7372 | -797
| $2^{3}S_{1}$ | $1^{+-}$ | 6609 | $\eta_{c}(2S)\psi(2S)$ | 7324 | -715
| $2^{5}S_{2}$ | $2^{++}$ | 6639 | $\psi(2S)\psi(2S)$ | 7372 | -733
| $3^{1}S_{0}$ | $0^{++}$ | 6782 | $\cdots$ | $\cdots$ | $\cdots$
| $3^{3}S_{1}$ | $1^{+-}$ | 6814 | $\cdots$ | $\cdots$ | $\cdots$
| $3^{5}S_{2}$ | $2^{++}$ | 6842 | $\cdots$ | $\cdots$ | $\cdots$
$bb\bar{b}\bar{b}$ | $1^{1}S_{0}$ | $0^{++}$ | 19666 | $\eta_{b}(1S)\eta_{b}(1S)$ | 18798 | 868
| $1^{1}S_{0}$ | $0^{++}$ | 19666 | $\Upsilon(1S)\Upsilon(1S)$ | 18920 | 746
| $1^{3}S_{1}$ | $1^{+-}$ | 19673 | $\eta_{b}(1S)\Upsilon(1S)$ | 18859 | 814
| $1^{5}S_{2}$ | $2^{++}$ | 19680 | $\Upsilon(1S)\Upsilon(1S)$ | 18920 | 760
| $2^{1}S_{0}$ | $0^{++}$ | 19841 | $\eta_{b}(1S)\eta_{b}(1S)$ | 18798 | 1043
| $2^{1}S_{0}$ | $0^{++}$ | 19841 | $\Upsilon(1S)\Upsilon(1S)$ | 18920 | 921
| $2^{3}S_{1}$ | $1^{+-}$ | 19849 | $\eta_{b}(1S)\Upsilon(1S)$ | 18859 | 990
| $2^{5}S_{2}$ | $2^{++}$ | 19855 | $\Upsilon(1S)\Upsilon(1S)$ | 18920 | 935
| $2^{1}S_{0}$ | $0^{++}$ | 19841 | $\eta_{b}(2S)\eta_{b}(2S)$ | 19998 | -157
| $2^{1}S_{0}$ | $0^{++}$ | 19841 | $\Upsilon(2S)\Upsilon(2S)$ | 20046 | -205
| $2^{3}S_{1}$ | $1^{+-}$ | 19849 | $\eta_{b}(2S)\Upsilon(2S)$ | 20022 | -173
| $2^{5}S_{2}$ | $2^{++}$ | 19855 | $\Upsilon(2S)\Upsilon(2S)$ | 20046 | -191
| $3^{1}S_{0}$ | $0^{++}$ | 20001 | $\cdots$ | $\cdots$ | $\cdots$
| $3^{3}S_{1}$ | $1^{+-}$ | 20012 | $\cdots$ | $\cdots$ | $\cdots$
| $3^{5}S_{2}$ | $2^{++}$ | 20021 | $\cdots$ | $\cdots$ | $\cdots$

We see from this table that the predicted masses of ground and first radially
excited $cc\bar{c}\bar{c}$ and $bb\bar{b}\bar{b}$ tetraquarks lie
significantly higher than the corresponding lowest (1S-1S) meson-meson
thresholds. However masses of first radially excited $cc\bar{c}\bar{c}$ and
$bb\bar{b}\bar{b}$ tetraquarks are below than their corresponding (2S-2S)
meson-meson thresholds. 1S $cc\bar{c}\bar{c}$ tetraquarks are located at
$\sim$ 6.3 GeV and the mass splitting is about 30 MeV. $J^{PC}=0^{++}$ 1S
$cc\bar{c}\bar{c}$ state is about $\sim$ 300 MeV above than $\eta_{c}\eta_{c}$
mass threshold while $\sim$ 100 MeV above than $J/\psi J/\psi$ mass threshold.
This suggests that $J^{PC}=0^{++}$ $cc\bar{c}\bar{c}$ state might easily decay
into $\eta_{c}\eta_{c}$ and $J/\psi J/\psi$ final states through quark
rearrangements, respectively. $J^{PC}=1^{+-}$ 1S $cc\bar{c}\bar{c}$ state lies
about 270 MeV above than the $\eta_{c}J/\psi$ mass threshold. It might easily
decay into $\eta_{c}J/\psi$ final state through quark rearrangements.
$J^{PC}=2^{++}$ 1S $cc\bar{c}\bar{c}$ tetraquark is about 200 MeV above than
the mass threshold of $J/\psi J/\psi$, which suggests that it might easily
decay into $J/\psi J/\psi$ final state through quark rearrangements. The
decays of 1S $cc\bar{c}\bar{c}$ tetraquarks are not suppressed dynamically or
kinematically. Therefore these states should be broad and hard to be observed
in experimental studies.

The previous discussion can be easily adapted to 1S $bb\bar{b}\bar{b}$
tetraquarks. The $J^{PC}=0^{++},1^{+-}$, and $J^{PC}=2^{++}$ 1S
$bb\bar{b}\bar{b}$ tetraquarks are located at $\sim$ 19.6 GeV and the mass
splitting is about 10 MeV. This result is expected since spin-spin interaction
term depend on mass as $\sim 1/m^{2}$. As the mass increases the contribution
of spin-spin interaction term becomes small. $J^{PC}=0^{++}$ 1S
$bb\bar{b}\bar{b}$ tetraquark is about $\sim$ 800 MeV above than
$\eta_{b}\eta_{b}$ mass threshold and $\sim$ 700 MeV above than
$\Upsilon\Upsilon$ mass threshold, respectively. This significant amount
suggest that it might easily decay into $\eta_{b}\eta_{b}$ and
$\Upsilon\Upsilon$ final states through quark rearrangements, respectively.
$J^{PC}=1^{+-}$ 1S $bb\bar{b}\bar{b}$ state lies about 800 MeV above than the
$\eta_{b}\Upsilon$ mass threshold. This significant amount suggests that it
might easily decay into $\eta_{b}\Upsilon$ final state through quark
rearrangements. $J^{PC}=2^{++}$ 1S $bb\bar{b}\bar{b}$ tetraquark is about 700
MeV above than the mass threshold of $\Upsilon\Upsilon$, and might easily
decay into $\Upsilon\Upsilon$ final state through quark rearrangements. The
decays of 1S $bb\bar{b}\bar{b}$ tetraquarks are not suppressed dynamically or
kinematically. Therefore these states should be broad and difficult to be
observed in experimental studies. This conclusion about 1S $bb\bar{b}\bar{b}$
tetraquarks agree with current experimental data of LHCb [39] and CMS [89]
Collaborations. No narrow beautiful tetraquarks in the $\Upsilon(1S)$-pair
production were observed in these experiments. Indeed, a lattice QCD study was
performed to investigate the existence of a $bb\bar{b}\bar{b}$ tetraquark
bound state by using the first-principles lattice nonrelativistic QCD
methodology [90]. They studied $S$-wave $bb\bar{b}\bar{b}$ tetraquark in three
channels $0^{++}$, $1^{+-}$ and $2^{++}$, and found no evidence of a QCD bound
tetraquark below the lowest noninteracting thresholds in the channels. A
similar conclusion was made in [47] by using a new compact tetraquark model.
They found that the decay modes of fully-$b$ with $J^{PC}=0^{++}$ tetraquarks
may be difficult to access in experiments. In Ref. [91], a quantitative
analysis of the $bb\bar{b}\bar{b}$ tetraquark decays into 4 muons,
$B^{+}B^{-}$, $B^{0}\bar{B}^{0}$, and $B^{0}_{s}\bar{B}^{0}_{s}$ channels was
made. They concluded that the 4 muon (4$\mu$) signal for the $J^{PC}=0^{++}$
ground state of $bb\bar{b}\bar{b}$ is likely to be too small for the upgraded
LHCb, but there is a hope for the $J^{PC}=2^{++}$ $bb\bar{b}\bar{b}$
tetraquark. A similar study with considering different channels was done in
[92]. They studied the production of the fully-heavy tetraquarks by focusing
on the $\gamma\gamma\rightarrow{\cal{Q}}{\cal{Q}}$
(${\cal{Q}}=J/\psi,\,\Upsilon$) subprocess mediated by the $T_{4Q}$ resonance
in the $s$-channel. They found that experimental study of this process is
feasible and can be used to explore the existence and properties of the
$T_{4c}(6900)$ and $T_{4b}(19000)$ states.

The obtained mass values of 2S $cc\bar{c}\bar{c}$ and $bb\bar{b}\bar{b}$
tetraquarks are significantly above than the fall-apart decays to the lowest
allowed (1S-1S) two-mesons. These states should be broad and difficult to be
observed experimentally.

On the other hand, the predicted masses of radially excited 2S
$cc\bar{c}\bar{c}$ and $bb\bar{b}\bar{b}$ tetraquarks are lower than their
corresponding thresholds. The $J^{PC}=0^{++},1^{+-}$, and $J^{PC}=2^{++}$ 2S
$cc\bar{c}\bar{c}$ tetraquarks are located at $\sim$ 6.6 GeV and mass
splitting is about 30 MeV. The $J^{PC}=0^{++},1^{+-}$, and $J^{PC}=2^{++}$ 2S
$bb\bar{b}\bar{b}$ tetraquarks are located at $\sim$ 19.8 GeV and mass
splitting is about 10 MeV. Very recently, the LCHb Collaboration observed a
broad structure at around $6.2$-$6.8$ GeV which is very close to twice the
$J/\psi$ mass threshold, a narrow structure around 6.9 GeV, and a vague
structure around 7.2 GeV [16]. We have several candidates of excited 2S
$cc\bar{c}\bar{c}$ tetraquarks in the $6.2$-$6.8$ mass range for the broad
structure observed by LHCb Collaboration [16]. Narrow structure around 6.9 GeV
might be orbital excitations of 2S or 3S $cc\bar{c}\bar{c}$ tetraquarks. A
recent study about LHCb of $J/\psi$ pairs was performed by using coupled
channel methods in two channels $J/\psi J/\psi$ and $\psi(2S)J/\psi$ and in
three channels $J/\psi J/\psi$, $\psi(2S)J/\psi$, and $\psi(3770)J/\psi$ [93].
They found that in both channels, there is a hint of an existence of a near-
threshold in the $J/\psi J/\psi$ system with $J^{PC}=0^{++}$ or
$J^{PC}=2^{++}$ quantum numbers having a mass around 6.2 GeV which they refer
this state as $X(6200)$. Our results of $J^{PC}=0^{++}$ and $J^{PC}=2^{++}$
$cc\bar{c}\bar{c}$ tetraquarks are $\sim$ 100 MeV above than their prediction.
In Ref. [94], it was suggested that the structure around 7.2 GeV which was
named as $X(7200)$ can be related to $\chi_{c1}(3872)$ (formerly known as
$X(3872)$) via heavy anti-quark di-quark symmetry. Orbital angular momentum
addition to our calculation for 3S $cc\bar{c}\bar{c}$ tetraquarks could yield
a candidate for $X(7200)$.

The obtained mass values of 2S $cc\bar{c}\bar{c}$ and $bb\bar{b}\bar{b}$
tetraquarks are significantly above than the thresholds of decays to the
lowest (1S-1S) two-mesons. These states should be broad and are difficult to
observe experimentally.

We compare our predictions for the masses of fully-heavy tetraquarks with the
results of previous calculations. At first step we compare results of some
works which considered excited fully-heavy tetraquarks. The results can be
seen in Table 8.

Table 8: Comparison of calculated fully-heavy tetraquark $cc\bar{c}\bar{c}$ and $cc\bar{c}\bar{c}$ $bb\bar{b}\bar{b}$ masses for ground (1S) and radial excited states (2S, 3S) with quantum numbers $J^{PC}=0^{++},1^{+-}$ and $2^{++}$. All results are in MeV. $J^{PC}$ | $0^{++}$ | $1^{+-}$ | $2^{++}$
---|---|---|---
State | $1S$ | $2S$ | $3S$ | $1S$ | $2S$ | $3S$ | $1S$ | $2S$ | $3S$
$cc\bar{c}\bar{c}$ |  |  |  |  |  |  |  |  |
This work | 6322 | 6575 | 6782 | 6354 | 6609 | 6814 | 6385 | 6639 | 6842
Ref. [95] | 6476 | 6908 | 7296 | 6441 | 6896 | 7300 | 6475 | 6921 | 7320
Ref. [96] | 5883 | 6573 | 6948 | 6120 | 6669 | 7016 | 6246 | 6739 | 7071
$bb\bar{b}\bar{b}$ |  |  |  |  |  |  |  |  |
This work | 19666 | 19841 | 20001 | 19673 | 19849 | 20012 | 19680 | 19855 | 20021
Ref. [95] | 19226 | 19583 | 19887 | 19214 | 19582 | 19989 | 19232 | 19594 | 19898
Ref. [96] | 18748 | 19335 | 19644 | 18828 | 19366 | 19665 | 18900 | 19398 | 19688

In Table 8, Ref. [95] calculated the masses and wave functions of the exotic
hadron states $cc\bar{c}\bar{c}$ and $bb\bar{b}\bar{b}$ solving the four-body
Schrödinger equation by including Cornell potential and spin-spin interaction
in vacuum and strongly interacting matter in quark model. Ref. [96] used a
relativized diquark model Hamiltonian built with relativistic kinetic
energies, a one gluon exchange potential and linear confinement, in order to
obtain masses of $cc\bar{c}\bar{c}$ and $bb\bar{b}\bar{b}$, and some other
systems. Our results are compatible compared to reference studies. The reason
of deviances in the results should be the different models used in the
studies. At the second step, we compare our results for the ground states of
fully-heavy tetraquarks of available studies. The results can be seen in Table
9.

Table 9: Comparison of theoretical predictions for the ground state masses of the $cc\bar{c}\bar{c}$ and $bb\bar{b}\bar{b}$ tetraquarks. All results are in MeV. Reference | $cc\bar{c}\bar{c}$ | $bb\bar{b}\bar{b}$
---|---|---
| $0^{++}$ | $1^{+-}$ | $2^{++}$ | $0^{++}$ | $1^{+-}$ | $2^{++}$
This work | 6322 | 6354 | 6385 | 19666 | 19673 | 19680
[17] | 5969 | 6021 | 6115 | $\cdots$ | $\cdots$ | $\cdots$
[25] | 6271 | 6231 | 6287 | 18981 | 18969 | 19000
[26] | 6190 | 6271 | 6367 | 19314 | 19320 | 19330
[36] | 6542 | 6515 | 6543 | 19255 | 19251 | 19262
[37] | 6314 | 6375 | 6407 | 19237 | 19264 | 19279
[41] | 5966 | 6051 | 6223 | 18754 | 18808 | 18916
[42] | 6440–6470 | 6370–6510 | 6370–6510 | 18460–18490 | 18320–18540 | 18320–18530
[43, 29] | $6192\pm 25$ | $\cdots$ | $6429\pm 25$ | $18826\pm 25$ | $\cdots$ | $18956\pm 25$
[101, 46] | 5990$\pm$ 80 | 6050$\pm$ 80 | 6090$\pm$ 80 | 18840 $\pm$ 90 | 18840$\pm$ 90 | 18850$\pm$ 90
[76] | 7016 | 6899 | 6956 | 20275 | 20212 | 20243
[79]-Mod. I | 5960 | 6009 | 6100 | 18723 | 18738 | 20243
[79]-Mod. II | 6198 | 6246 | 6323 | 18754 | 18768 | 18797
[95] | 6476 | 6441 | 6475 | 19226 | 19214 | 19232
[96] | 5883 | 6120 | 6246 | 18748 | 18828 | 18900
[97] | 6487 | 6500 | 6524 | 19322 | 19329 | 19341
[98] | 6425 | 6425 | 6432 | 19247 | 19247 | 19249
[99] | 6407 | 6463 | 6486 | 19329 | 19373 | 19387
[100] | $\cdots$ | $\cdots$ | $\cdots$ | 19178 | 19226 | 19236
[102] | 6695 | 6528 | 6573 | $\cdots$ | $\cdots$ | $\cdots$
[103] | 6480 | 6508 | 6565 | $\cdots$ | $\cdots$ | $\cdots$

Ref. [97] used a nonrelativistic quark model within a potential model by
including the linear confining potential, Coulomb potential, and spin-spin
interactions in diquark-antidiquark picture. Ref. [98] studied the mass
spectra of the $S-$wave fully-heavy tetraquarks in two nonrelativistic quark
models. Color-magnetic interaction model (CMIM), a traditional constituent
quark model (CQM) and a multiquark color flux-tube model (MCFTM) were used to
investigate the fully-heavy tetraquarks in diquark-antidiquark scenario [99].
Ref. [100] used nonrelativistic chiral quark model using the Gaussian
expansion method for ground state of $bb\bar{b}\bar{b}$ tetraquark. Refs. [46,
101] made QCD Sum Rule studies on fully-heavy tetraquarks. In another work,
only fully-charmed tetraquarks were considered using a parameterized
Hamiltonian with a sufficiently large but finite oscillator basis [102]. Ref.
[103] calculated tetraquark masses in a constituent quark model, where the
Cornell-like potential and one-gluon exchange spin-spin coupling were taken
into account.

It can be seen from Tables 8 and 9 that there are significant disagreements
between the masses of fully-heavy $cc\bar{c}\bar{c}$ and $bb\bar{b}\bar{b}$
tetraquarks. Especially from Table 8, in the $cc\bar{c}\bar{c}$ sector, masses
of $J^{PC}=0^{++},1^{+-}$, and $J^{PC}=2^{++}$ lie in the range of $5.8-7.0$
GeV, $6.0-6.8$ GeV and $6.1-6.9$ GeV, respectively. In the $bb\bar{b}\bar{b}$
sector, masses of $J^{PC}=0^{++},1^{+-}$, and $J^{PC}=2^{++}$ lie in the range
of $18.7-20.0$ GeV, $18.3-20.2$ GeV and $18.3-20.2$ GeV, respectively. Refs.
[17, 29, 43, 46, 101, 79, 96] predict ground state fully-heavy tetraquark
masses below or slightly above the corresponding two meson threshold masses
which means these states could be observed in the experiments. On the other
side, present model of this work and other approaches predict ground state
fully-heavy tetraquark masses significantly above the corresponding
thresholds, and thus it is difficult to observe these states in the
experiments. 2S fully-heavy tetraquark masses are significantly above then the
corresponding lowest (1S-1S) two-meson threshold masses. Concerning excited 2S
fully-heavy tetraquark masses, our results are significantly below than the
corresponding two-meson thresholds and there can be an open window for
observing these structures in the experiments.

It would be good to analyze the response of the obtained results to minor
variations of the numerical input of the potential model. The potential model
has eight parameters. Therefore we divided into two group. In the first group
we have changed mass $m$ and $\alpha_{s}$ values and fixed the rest of
parameters. In the second group we have changed coupling strength $h$, range
$b$, $c$ and $V_{0}$ values and fixed the rest of parameters. The parameters
of group 1 is shown in Table 10.

Table 10:  Group 1 potential parameters.  Input | Set 1 | Set 2 | Set 3 | Set 4
---|---|---|---|---
$m_{c}$ | $1.41~{}\text{GeV}$ | $1.39~{}\text{GeV}$ | $1.4~{}\text{GeV}$ | $1.4~{}\text{GeV}$
$m_{b}$ | $4.822~{}\text{GeV}$ | $4.802~{}\text{GeV}$ | $4.812~{}\text{GeV}$ | $4.812~{}\text{GeV}$
$\alpha_{s}(m_{c}^{2})$ | 0.37 | 0.37 | 0.38 | 0.36
$\alpha_{s}(m_{b}^{2})$ | 0.26 | 0.26 | 0.27 | 0.25

As can be seen from Table 10, in Sets 1 and 2 we fixed $\alpha_{s}$ values
constant and changed the mass values while in Sets 3 and 4, we fixed the mass
values and changed $\alpha_{s}$ values. The results are given in Table 11 for
$cc\bar{c}\bar{c}$ tetraquark and in Table 12 for $bb\bar{b}\bar{b}$
tetraquark.

Table 11:  Responses in $cc\bar{c}\bar{c}$ tetraquark mass values to potential parameter input in group 1. PM stands for results of potential model used in this work. All results are in MeV. State | Set 1 | Set 2 | Set 3 | Set 4 | PM
---|---|---|---|---|---
$1^{1}S_{0}$ | 6331 | 6314 | 6320 | 6321 | 6322
$1^{3}S_{1}$ | 6364 | 6342 | 6351 | 6353 | 6354
$1^{5}S_{2}$ | 6393 | 6371 | 6383 | 6384 | 6385
Table 12:  Responses in $bb\bar{b}\bar{b}$ tetraquark mass values to potential parameter input in group 1. PM stands for results of potential model used in this work. All results are in MeV. State | Set 1 | Set 2 | Set 3 | Set 4 | PM
---|---|---|---|---|---
$1^{1}S_{0}$ | 19687 | 19646 | 19665 | 19666 | 19666
$1^{3}S_{1}$ | 19698 | 19654 | 19672 | 19673 | 19673
$1^{5}S_{2}$ | 19708 | 19759 | 19680 | 19680 | 19680

As can be seen from Tables 11 and 12, the responses in the numerical input
changes of group 1 did not significantly effect the results. The parameters of
group 2 are given in Table 13.

Table 13:  Group 2 potential parameters.  Input | Set 5 | Set 6 | Set 7 | Set 8
---|---|---|---|---
$h$ | $0.21~{}\text{GeV}$ | $0.19~{}\text{GeV}$ | $0.20~{}\text{GeV}$ | $0.20~{}\text{GeV}$
$b$ | $0.5~{}\text{GeV}$ | $0.3~{}\text{GeV}$ | $0.4~{}\text{GeV}$ | $0.4~{}\text{GeV}$
$c$ | $0.193~{}\text{GeV}^{2}$ | $0.193~{}\text{GeV}^{2}$ | $0.203~{}\text{GeV}^{2}$ | $0.183~{}\text{GeV}^{2}$
$V_{0}$ | $-0.223~{}\text{GeV}$ | $-0.223~{}\text{GeV}$ | $-0.233~{}\text{GeV}$ | $0.213~{}\text{GeV}$

As can be seen from Table 13, in Sets 5 and 6 we fixed $c$ and $V_{0}$ values
constant and changed $h$ and $b$ values while in Sets 7 and 8, we fixed $h$
and $b$ values and changed $c$ and $V_{0}$ values. The results are given in
Table 14 for $cc\bar{c}\bar{c}$ tetraquark and in Table 15 for
$bb\bar{b}\bar{b}$ tetraquark.

Table 14:  Responses in $cc\bar{c}\bar{c}$ tetraquark mass values to potential parameter input in group 2. PM stands for results of potential model used in this work. All results are in MeV. State | Set 1 | Set 2 | Set 3 | Set 4 | PM
---|---|---|---|---|---
$1^{1}S_{0}$ | 6314 | 6329 | 6330 | 6314 | 6322
$1^{3}S_{1}$ | 6345 | 6372 | 6362 | 6344 | 6354
$1^{5}S_{2}$ | 6376 | 6401 | 6394 | 6375 | 6385
Table 15:  Responses in $bb\bar{b}\bar{b}$ tetraquark mass values to potential parameter input in group 2. PM stands for results of potential model used in this work. All results are in MeV. State | Set 1 | Set 2 | Set 3 | Set 4 | PM
---|---|---|---|---|---
$1^{1}S_{0}$ | 19676 | 19680 | 19667 | 19664 | 19666
$1^{3}S_{1}$ | 19682 | 19686 | 19674 | 19672 | 19673
$1^{5}S_{2}$ | 19689 | 19694 | 19680 | 19768 | 19680

As can be seen from Tables 14 and 15, the responses in the numerical input
changes of group 2 did not significantly effect the results.

We believe that, the results of this work should serve as a benchmark for
other studies, which are necessary to identify and develop theoretical models.
This effort is important since predictions with respect to different models or
approaches should be in well agreement and depend mildly to the corresponding
model parameters.

##  IV Summary and Concluding Remarks

In this work, adapting a perspective in which tetraquarks are assumed to be
compact and consist of axial-vector diquarks and antidiquarks, masses of
fully-heavy tetraquark states were calculated. We used a nonrelativistic
potential quark model with the linear confining potential, Hulthen potential,
and spin-spin interactions. We predicted ground, first and second radially
excited mass spectra of fully-heavy tetraquarks, $cc\bar{c}\bar{c}$ and
$bb\bar{b}\bar{b}$. In order to do this, we adopted a model which was
formulated for quarkonia systems and made use of this model for diquark and
tetraquark systems.

We have factorized four-body system into two two-body systems: quarks form
diquark and antiquarks form antidiquarks, then diquark and antidiquark form
tetraquark. By doing this, the complex structure of four-body problem became
more accessible. This flowchart allows us to study both ground and excited
tetraquark states. Furthermore, various spin configurations were conducted
easily in this manner.

At first step, we obtained diquark/antidiquark masses. Since diquark masses
equal antidiquark masses, in the second step we computed tetraquark masses. In
comparison with other models, the obtained $cc$ and $bb$ diquark masses are in
good agreement with earlier results. This situation can be seen in Tables 4
and 6. Indeed, radial excited states of $cc$ and $bb$ diquark masses agree
well with available results in the literature which can be seen in Tables 3
and 5. After obtaining diquark masses, we computed ground, first, and second
radially excited $cc\bar{c}\bar{c}$ and $bb\bar{b}\bar{b}$ tetraquark masses.
The results are presented in Table 7. It was found that both ground states of
$cc\bar{c}\bar{c}$ and $bb\bar{b}\bar{b}$ tetraquarks masses are considerably
higher than their corresponding lowest (1S-1S) two-meson decay threshold
masses. Therefore they might easily decay into the lowest two-meson allowed
states and these decays proceed through quark rearrangements. Such decays are
not suppressed dynamically or kinematically. These states should be broad and
are difficult to observe experimentally. We observed that first excited states
of fully-charmed and fully-bottomed tetraquarks are lower than their
corresponding (2S-2S) two-meson thresholds. These tetraquarks might be
observed by obtaining enough center of mass energy and luminosity in the
experiments. We also observed that the masses of 2S $cc\bar{c}\bar{c}$ and
$bb\bar{b}\bar{b}$ tetraquarks are significantly higher than the corresponding
lowest (1S-1S) two-meson threshold masses. First and second radially excited
masses of fully-heavy tetraquarks are compatible with the reference studies
given in Table 8.

The comparison of ground state masses of $cc\bar{c}\bar{c}$ and
$bb\bar{b}\bar{b}$ tetraquarks are presented in Table 9. As argued before, the
results deviate up to 1.0 GeV for $cc\bar{c}\bar{c}$ tetraquark and 2.0 GeV
for $bb\bar{b}\bar{b}$ tetraquark which are considered as high differences.
Significant deviations between diquark and tetraquark masses may present a
deeper perspective. In quark model, quarks are point-like objects. The
dynamics, interactions and bound states for two and three quarks systems are
described well. The diquark-antidiquark models require additional assumptions.
In this perspective most diquark-antidiquark models assume diquarks and
antidiquarks are point-like objects so that the interactions are effectively
the same as for ordinary quarkonia systems. This was argued, for example in
Ref. [26], in which diquarks have a structure (size) and this effect the
results. In addition to this, chromomagnetic interaction inside the diquarks
or between diquarks are hypothesized [80, 104].

To sum up, we believe that our predictions for fully-charmed and fully-
bottomed ground and excited states may help to experimental studies for a
possible tetraquark signal in corresponding channels.

## References

  * [1] M. Gell-Mann, Phys. Lett. 8 (1964), 214-215
  * [2] R. L. Jaffe, Phys. Rev. D 15 (1977), 267.
  * [3] R. L. Jaffe, Phys. Rev. D 15 (1977), 281.
  * [4] J. P. Ader, J. M. Richard and P. Taxil, Phys. Rev. D 25 (1982), 2370.
  * [5] S. K. Choi et al. [Belle], Phys. Rev. Lett. 91 (2003), 262001, [arXiv:hep-ex/0309032 [hep-ex]].
  * [6] P. del Amo Sanchez et al. [BaBar], Phys. Rev. D 82 (2010), 011101, [arXiv:1005.5190 [hep-ex]].
  * [7] B. Aubert et al. [BaBar], Phys. Rev. Lett. 102 (2009), 132001, [arXiv:0809.0042 [hep-ex]].
  * [8] R. Aaij et al. [LHCb], Nucl. Phys. B 886 (2014), 665-680 ,[arXiv:1404.0275 [hep-ex]].
  * [9] B. Aubert et al. [BaBar], Phys. Rev. Lett. 95 (2005), 142001, [arXiv:hep-ex/0506081 [hep-ex]].
  * [10] M. Ablikim et al. [BESIII], Phys. Rev. Lett. 110 (2013), 252001, [arXiv:1303.5949 [hep-ex]].
  * [11] Z. Q. Liu et al. [Belle], Phys. Rev. Lett. 110 (2013), 252002 [erratum: Phys. Rev. Lett. 111 (2013), 019901], [arXiv:1304.0121 [hep-ex]].
  * [12] R. Aaij et al. [LHCb], Phys. Rev. Lett. 115 (2015), 072001, [arXiv:1507.03414 [hep-ex]].
  * [13] Z. Ghalenovi and M. M. Sorkhi, Eur. Phys. J. Plus 135 (2020) no.5, 399.
  * [14] F. K. Guo, C. Hanhart, U. G. Meißner, Q. Wang, Q. Zhao and B. S. Zou, Rev. Mod. Phys. 90 (2018) no.1, 015004, [arXiv:1705.00141 [hep-ph]].
  * [15] A. Ali, L. Maiani and A. D. Polosa, “Multiquark Hadrons,” (Cambridge University Press, Cambridge, 2019).
  * [16] R. Aaij et al. [LHCb], Sci. Bull. 65 (2020), 1983-1993, [arXiv:2006.16957 [hep-ex]].
  * [17] V. R. Debastiani and F. S. Navarra, Chin. Phys. C 43 (2019) no.1, 013105, [arXiv:1706.07553 [hep-ph]].
  * [18] B. C. Yang, L. Tang and C. F. Qiao, [arXiv:2012.04463 [hep-ph]].
  * [19] B. D. Wan and C. F. Qiao, [arXiv:2012.00454 [hep-ph]].
  * [20] C. Gong, M. C. Du, B. Zhou, Q. Zhao and X. H. Zhong, [arXiv:2011.11374 [hep-ph]].
  * [21] J. W. Zhu, X. D. Guo, R. Y. Zhang, W. G. Ma and X. Q. Li, [arXiv:2011.07799 [hep-ph]].
  * [22] Q. F. Cao, H. Chen, H. R. Qi and H. Q. Zheng, [arXiv:2011.04347 [hep-ph]].
  * [23] Z. H. Guo and J. A. Oller, [arXiv:2011.00978 [hep-ph]].
  * [24] J. R. Zhang, [arXiv:2010.07719 [hep-ph]].
  * [25] X. Z. Weng, X. L. Chen, W. Z. Deng and S. L. Zhu, [arXiv:2010.05163 [hep-ph]].
  * [26] R. N. Faustov, V. O. Galkin and E. M. Savchenko, Phys. Rev. D 102 (2020), 114030, [arXiv:2009.13237 [hep-ph]].
  * [27] Y. Q. Ma and H. F. Zhang, [arXiv:2009.08376 [hep-ph]].
  * [28] Z. G. Wang, [arXiv:2009.05371 [hep-ph]].
  * [29] M. Karliner and J. L. Rosner, Phys. Rev. D 102 (2020) no.11, 114039, [arXiv:2009.04429 [hep-ph]].
  * [30] R. Maciuła, W. Schäfer and A. Szczurek, Phys. Lett. B 812 (2021), 136010, [arXiv:2009.02100 [hep-ph]].
  * [31] G. Eichmann, C. S. Fischer, W. Heupel, N. Santowsky and P. C. Wallbott, Few Body Syst. 61 (2020) no.4, 38, [arXiv:2008.10240 [hep-ph]].
  * [32] J. Z. Wang, D. Y. Chen, X. Liu and T. Matsuki, [arXiv:2008.07430 [hep-ph]].
  * [33] K. T. Chao and S. L. Zhu, Sci. Bull. 65 (2020), 1952-1953 [arXiv:2008.07670 [hep-ph]].
  * [34] J. M. Richard, [arXiv:2008.01962 [hep-ph]].
  * [35] C. Becchi, J. Ferretti, A. Giachino, L. Maiani and E. Santopinto, Phys. Lett. B 811 (2020), 135952, [arXiv:2006.14388 [hep-ph]].
  * [36] Q. F. Lü, D. Y. Chen and Y. B. Dong, Eur. Phys. J. C 80 (2020) no.9, 871, [arXiv:2006.14445 [hep-ph]].
  * [37] X. Jin, Y. Xue, H. Huang and J. Ping, Eur. Phys. J. C 80 (2020) no.11, 1083, [arXiv:2006.13745 [hep-ph]].
  * [38] S. Durgut, Search for Exotic Mesons at CMS, 2018 (http:// meetings.aps.org/Meeting/APR18/Session/U09.6)
  * [39] R. Aaij et al. [LHCb], JHEP 10 (2018), 086, [arXiv:1806.09707 [hep-ex]].
  * [40] L. Heller and J. A. Tjon, Phys. Rev. D 32 (1985), 755.
  * [41] A. V. Berezhnoy, A. V. Luchinsky and A. A. Novoselov, Phys. Rev. D 86 (2012), 034004, [arXiv:1111.1867 [hep-ph]].
  * [42] W. Chen, H. X. Chen, X. Liu, T. G. Steele and S. L. Zhu, Phys. Lett. B 773 (2017), 247-251, [arXiv:1605.01647 [hep-ph]].
  * [43] M. Karliner, S. Nussinov and J. L. Rosner, Phys. Rev. D 95 (2017) no.3, 034011, [arXiv:1611.00348 [hep-ph]].
  * [44] Y. Bai, S. Lu and J. Osborne, Phys. Lett. B 798 (2019), 134930, [arXiv:1612.00012 [hep-ph]].
  * [45] E. Eichten and Z. Liu, [arXiv:1709.09605 [hep-ph]].
  * [46] Z. G. Wang, Eur. Phys. J. C 77 (2017) no.7, 432, [arXiv:1701.04285 [hep-ph]].
  * [47] A. Esposito and A. D. Polosa, Eur. Phys. J. C 78 (2018) no.9, 782, [arXiv:1807.06040 [hep-ph]].
  * [48] H. X. Chen, W. Chen, X. Liu and S. L. Zhu, Sci. Bull. 65 (2020), 1994-2000, [arXiv:2006.16027 [hep-ph]].
  * [49] R. Aaij et al. [LHCb], Phys. Rev. Lett. 119 (2017) no.11, 112001, [arXiv:1707.01621 [hep-ex]].
  * [50] G. Malli, J. Oreg, Chem. Phys. Lett. 69, 313–314 (1980).
  * [51] J. Lindhard, P.G. Hansen, Phys. Rev. Lett. 57, 965–967 (1986).
  * [52] C.S. Jia, J.Y.Wang, S. He, L.T. Sun, J. Phys. A Math. Gen. 33, 6993–6998 (2000).
  * [53] N. Candemir, O. Bayrak, J. Math. Phys. 54, 042104 (2013).
  * [54] M. K. Bahar, A. Soylu and A. Poszwa, IEEE Transactions on Plasma Science, 44, no. 10, 2297-2306 (2016).
  * [55] H. H. Karayer, A. D. Demirhan, F. Büyükkılıç, Turk J. Phys. 41 551-559 (2017).
  * [56] P. Pyykko, J. Jokisaari, Chem. Phys. 10, 293–301 (1975)
  * [57] S. M. Ikhdair, R. Sever, J. Math. Chem. 42, 3, 461-471 (2007).
  * [58] A. K. Roy, Int. J. Quantum Chem. 115, 937–947 (2015).
  * [59] C.S. Lam, Y.P. Varshni, Phys. Rev. A 4, 1875 (1971).
  * [60] R.L. Hall, Phys. Rev. A 32, 14 (1985).
  * [61] U. Laha and J. Bhoi, Phys. Rev. C 91 no.3, 034614 (2015).
  * [62] J. Bhoi and U. Laha, J. Phys. G 40, 045107 (2013).
  * [63] M. Hosseinpour, F. M. Andrade, E. O. Silva and H. Hassanabadi, Eur. Phys. J. C 77 no.5, 270 (2017) [erratum: Eur. Phys. J. C 77 no.6, 373 (2017)], [arXiv:1608.03558 [hep-th]].
  * [64] M. Aktas, EPL, 121 10005, (2018).
  * [65] Bhaghyesh, K. B. Vijaya Kumar and A. P. Monteiro, J. Phys. G 38 (2011), 085001.
  * [66] John R. Taylor, Scattering Theory: The Quantum Theory of Nonrelativistic Collisions, 2006, Dover Publications.
  * [67] W. Lucha and F. F. Schöberl, Int. J. Mod. Phys. A 29 (2014) no.29, 1450181, [arXiv:1408.4957 [hep-ph]].
  * [68] A. De Rujula, H. Georgi and S. L. Glashow, Phys. Rev. D 12 (1975), 147-162.
  * [69] E. Eichten and F. Feinberg, Phys. Rev. D 23 (1981), 2724 doi:10.1103/PhysRevD.23.2724
  * [70] W. Lucha, F. F. Schoberl and D. Gromes, Phys. Rept. 200 (1991), 127-240.
  * [71] W. Lucha and F. F. Schoberl, in International Summer School for Students on Development in Nuclear Theory and Particle Physics (1995), [arXiv:hep-ph/9601263 [hep-ph]].
  * [72] M. B. Voloshin, Prog. Part. Nucl. Phys. 61 (2008), 455-511, [arXiv:0711.4556 [hep-ph]].
  * [73] Bhaghyesh, K. B. Vijaya Kumar and Y. L. Ma, Int. J. Mod. Phys. A 27 (2012), 1250011.
  * [74] S. S. Gershtein, V. V. Kiselev, A. K. Likhoded and A. V. Tkabladze, Phys. Rev. D 51 (1995), 3613-3627, [arXiv:hep-ph/9406339 [hep-ph]].
  * [75] W. Lucha and F. F. Schoberl, Int. J. Mod. Phys. C 10, 607-620 (1999), [arXiv:hep-ph/9811453 [hep-ph]].
  * [76] J. Wu, Y. R. Liu, K. Chen, X. Liu and S. L. Zhu, Phys. Rev. D 97 (2018) no.9, 094015, [arXiv:1605.01134 [hep-ph]].
  * [77] D. Ebert, R. N. Faustov, V. O. Galkin and W. Lucha, Phys. Rev. D 76 (2007), 114015, [arXiv:0706.3853 [hep-ph]].
  * [78] Q. F. Lü and Y. B. Dong, Phys. Rev. D 94 (2016) no.9, 094041, [arXiv:1603.06417 [hep-ph]].
  * [79] P. Lundhammar and T. Ohlsson, Phys. Rev. D 102 (2020) no.5, 054018, [arXiv:2006.09393 [hep-ph]].
  * [80] L. Maiani, F. Piccinini, A. D. Polosa and V. Riquer, Phys. Rev. D 89 (2014), 114010, [arXiv:1405.1551 [hep-ph]].
  * [81] V. V. Kiselev, A. K. Likhoded, O. N. Pakhomova and V. A. Saleev, Phys. Rev. D 66 (2002), 034030, [arXiv:hep-ph/0206140 [hep-ph]].
  * [82] J. Ferretti, Few Body Syst. 60 (2019) no.1, 17.
  * [83] M. A. Bedolla, Few Body Syst. 60 (2019) no.2, 24.
  * [84] S. Esau, A. Palameta, R. T. Kleiv, D. Harnett and T. G. Steele, Phys. Rev. D 100 (2019), 074025, [arXiv:1905.12803 [hep-ph]].
  * [85] Q. X. Yu and X. H. Guo, Nucl. Phys. B 947 (2019), 114727, [arXiv:1810.00437 [hep-ph]].
  * [86] S. S. Gershtein, V. V. Kiselev, A. K. Likhoded and A. I. Onishchenko, Phys. Rev. D 62 (2000), 054021.
  * [87] M. N. Anwar, J. Ferretti, F. K. Guo, E. Santopinto and B. S. Zou, Eur. Phys. J. C 78 (2018) no.8, 647, [arXiv:1710.02540 [hep-ph]].
  * [88] L. Maiani, F. Piccinini, A. D. Polosa and V. Riquer, Phys. Rev. D 71 (2005), 014028, [arXiv:hep-ph/0412098 [hep-ph]].
  * [89] A. M. Sirunyan et al. [CMS], Phys. Lett. B 808 (2020), 135578, [arXiv:2002.06393 [hep-ex]].
  * [90] C. Hughes, E. Eichten and C. T. H. Davies, Phys. Rev. D 97 (2018) no.5, 054505, [arXiv:1710.03236 [hep-lat]].
  * [91] C. Becchi, A. Giachino, L. Maiani and E. Santopinto, Phys. Lett. B 806 (2020), 135495, [arXiv:2002.11077 [hep-ph]].
  * [92] V. P. Goncalves and B. D. Moreira, [arXiv:2101.03798 [hep-ph]].
  * [93] X. K. Dong, V. Baru, F. K. Guo, C. Hanhart and A. Nefediev, [arXiv:2009.07795 [hep-ph]].
  * [94] M. Z. Liu and L. S. Geng, [arXiv:2012.05096 [hep-ph]].
  * [95] J. Zhao, S. Shi and P. Zhuang, Phys. Rev. D 102 (2020) no.11, 114001, [arXiv:2009.10319 [hep-ph]].
  * [96] M. A. Bedolla, J. Ferretti, C. D. Roberts and E. Santopinto, Eur. Phys. J. C 80 (2020) no.11, 1004, [arXiv:1911.00960 [hep-ph]].
  * [97] M. S. Liu, Q. F. Lü, X. H. Zhong and Q. Zhao, Phys. Rev. D 100 (2019) no.1, 016006, [arXiv:1901.02564 [hep-ph]].
  * [98] G. J. Wang, L. Meng and S. L. Zhu, Phys. Rev. D 100 (2019) no.9, 096013, [arXiv:1907.05177 [hep-ph]].
  * [99] C. Deng, H. Chen and J. Ping, Phys. Rev. D 103 (2021) no.1, 014001, [arXiv:2003.05154 [hep-ph]].
  * [100] X. Chen, Eur. Phys. J. A 55 (2019) no.7, 106, [arXiv:1902.00008 [hep-ph]].
  * [101] Z. G. Wang and Z. Y. Di, Acta Phys. Polon. B 50 (2019), 1335, [arXiv:1807.08520 [hep-ph]].
  * [102] R. J. Lloyd and J. P. Vary, Phys. Rev. D 70 (2004), 014009, [arXiv:hep-ph/0311179 [hep-ph]].
  * [103] Z. Zhao, K. Xu, A. Kaewsnod, X. Liu, A. Limphirat and Y. Yan, [arXiv:2012.15554 [hep-ph]].
  * [104] J. M. Richard, A. Valcarce and J. Vijande, Phys. Rev. D 95 (2017) no.5, 054019, [arXiv:1703.00783 [hep-ph]].