Doubly heavy tetraquarks in the chiral quark soliton model

Recent discovery of a doubly charmed tetraquark ${\cal T}_{cc}^{+}$ with a mass of $\sim 3875$ MeV, approximately $\sim 300$ keV below the $D^{*+}D^{0}$ threshold, by the LHCb Collaboration [1, 2] triggered a number of theoretical studies of exotic heavy-light states. A comprehensive review of multiquark states, both experimental and theoretical, before ${\cal T}_{cc}^{+}$ discovery can be found in Ref. [3] and more recently after the discovery of ${\cal T}_{cc}^{+}$ in Ref. [4] and references therein. The up to date compilation of theoretical results is best illustrated in Fig. 42 of Ref. [4].

The existence of heavy tetraquarks has been anticipated theoretically already many yers ago [5]. In 1993, Manohar and Wise [6] showed using heavy quark symmetry [7] that $QQq_{1}q_{2}$ tetraquarks are bound in the limit $m_{Q}\rightarrow\infty$ (see also [8, 9]). This has been also pointed out more recently in Ref. [10]. To the best of our knowledge, the first estimate of a doubly heavy tetraquark mass is from Lipkin in 1986 [11] (although the fourfold heavy tetraquarks were discussed even earlier in 1982 [12]). We have reviewed the variational approach of Ref.  [11] in Ref. [13] adding new information coming from the discovery of $\Xi_{cc}^{++}(3621)$ [14] and showing that the upper bound on a ${\cal T}_{cc}^{+}$ mass is approximately 60 MeV above the $DD^{*}$ threshold. On the contrary, the bound on a ${\cal T}_{bb}$ mass was 224 MeV below the threshold. In the same paper, we advocated the possibility of using the chiral quark soliton model ($\chi$QSM) to estimate the ${\cal T}_{QQ}$ mass.

A mean field description of heavy baryons as a light quark-soliton and a heavy quark has been introduced and developed in Refs. [15, 16, 17, 18]. This approach is a modification of the $\chi$QSM used previously to describe light baryons (see [19] and Refs. [20, 21, 22] for review) where the soliton is constructed from $N_{c}$ light quarks. To describe heavy baryons one has to remove one light quark from the valence level and add a heavy quark instead. In the large $N_{c}$ limit this replacement hardly changes the mean fields of the soliton.

Support for such a treatment can be inferred from Ref. [23] where the authors studied soliton behavior in the limit where the current quark masses are $m\rightarrow\infty$. Although such a limit may at first sight be in contradiction with the chiral symmetry, which is the main theoretical basis of the model, it gave very good phenomenological results when compared to lattice data at finite $m_{\pi}$. At sufficiently large $m$, the soliton ceases to exist, and the correct heavy quark limit is achieved.

In the $\chi$QSM, the soliton mass is given as a sum over the energies of the valence quarks and the sea quark energies computed with respect to the vacuum and appropriately regularized [23],

In the present context, Eq. (1) takes the following form:

As has been argued in Ref. [23], for large $m_{Q}$, the sum over the sea quarks in the second line of Eq. (2) vanishes, and $E_{\rm val}(m_{Q})\approx m_{Q}$. One copy of the soliton ceases to exist; however, the remaining $N_{c}-1$ quarks still form a stable soliton.

Such a soliton does not carry any quantum numbers except for the baryon number resulting from the valence quarks. Spin and isospin appear when the soliton rotations in space and flavor are quantized. This procedure results in a collective Hamiltonian analogous to the one of a quantum mechanical symmetric top; however, due to the Wess-Zumino-Witten term [24, 25], the allowed Hilbert space is truncated to the representations that contain states of hypercharge $Y^{\prime}=N_{\rm val}/3$. For $N_{\rm val}=N_{c}=3$ (1), these are octet and decuplet of ground state baryons. For $N_{\rm val}=N_{c}-1=2$ (2), we have antitriplet of spin 0 and sextet of spin 1. It is therefore convenient to label heavy quark baryons (and tetraquarks as well) by the SU(3)_flavor representation of the light subsystem.

From this perspective, the soliton is reminiscent of a diquark, and the quantization rule $Y^{\prime}=(N_{c}-1)/3$ selects SU(3)_flavor representations identical as the ones of the quark model. Given the success of the $\chi$QSM in reproducing the data [15, 16, 17, 18], we propose here to use the same strategy to describe the doubly heavy tetraquarks replacing heavy quark $Q$ by a heavy (anti)diquark $\bar{Q}\bar{Q}$.

We observe that two heavy quarks of the same flavor (say $cc$ or $bb$) can form a color antitriplet (antisymmetric in color) provided they are symmetric in spin [26]. Therefore, they form a tight object of spin 1. Hence, two heavy antiquarks are in color $\mathbf{3}$ and spin 1, behaving as a spin 1 heavy quark. A tetraquark can be therefore viewed as being composed of a heavy (anti)diquark of spin 1 and a ($N_{c}-1$)-quark soliton.¹¹1In what follows, we use term diquark referring both to $QQ$ and $\bar{Q}\bar{Q}$ states.

There are three main lessons that we have learned from our previous studies of heavy baryons [15, 16, 17, 18]:

• •

  the soliton properties do not depend on the mass of the heavy quark,

• •

  neither do they depend on the spin coupling between a soliton and a heavy quark,

• •

  hyperfine splittings scale like $1/m_{Q}$.

This is discussed in detail in Sec. II.

Therefore, a very simple and predictive picture of a soliton + heavy object (that is $\boldsymbol{\overline{3}}$ in color) bound state emerges, where the mass is simply given as a sum of the soliton mass (including $m_{s}$ and rotational splittings), mass of a heavy object (quark or a diquark), and the hyperfine splitting. This picture is very reminiscent to the one of Ref. [10]. Mass formulas for such states are therefore identical to the ones of heavy baryons, with some modification due to the spin 1 character of the heavy diquark; this is elaborated on in detail in Sec. III. So the main problem is to estimate the diquark mass. Here, we propose to use the Cornell potential as described later in Sec. IV.

We find that only bottom antitriplet tetraquarks, both nonstrange and strange, are bound by approximately 140 and 60 MeV, respectively. We present numerical results for antitriplet and sextet tetraquarks in Sec. V and conclude in Sec. VI.

Let us first recall how baryon masses are calculated in the present model. We quantize the soliton as if it were constructed from $N_{c}-1$ rather than $N_{c}$ light quarks. Then, in the chiral limit the soliton energy is given as

Here, $M_{\mathrm{sol}}$ is a classical soliton mass, $I_{1,2}$ denote moments of inertia, $C_{2}(p,q)$ is the SU(3)_flavor Casimir, and $J$ corresponds to spin. In our case $Y^{\prime}=2/3$ and the allowed SU(3) representations correspond to $\boldsymbol{\overline{3}}$ with spin $J=0$ and $\boldsymbol{6}$ with spin $J=1$ [15].

SU(3) splittings are given by the operator

where constants $\alpha$, $\beta$ and $\gamma$ can be expressed through generalized moments of inertia (see e.g. Eq. (4) in Ref [15]) and can be computed ab initio in some specific versions of the model. In the most simple case with the pseudoscalar fields only, the numerical values can be found e.g. in Ref. [27], and in the context of heavy baryons in Ref. [28]. In both cases they lead to reasonable phenomenology. However, in reality one should take into account all possible chiral fields: scalar, pseudoscalar, vector, axial and tensor [29] for which full numerical analysis has not been performed. Here, the explicit forms of $\alpha,~{}\beta$, and $\gamma$ are not needed as we treat them as free parameters.

Heavy baryon masses are calculated by adding the mass of the heavy quark to the soliton mass and by taking into account the hyperfine splitting given by the following Hamiltonian:

where ${\boldsymbol{J}}$ and ${\boldsymbol{S}}_{Q}$ stand for the soliton and heavy quark spin, respectively. We have assumed here that the possible $m_{Q}$ dependence of $\varkappa$, due to the presence of the wave function squared in (5), can be ignored. Since the spin of the $\boldsymbol{\overline{3}}$ representation is zero, there is no hyperfine splitting in this case. In the case of $\boldsymbol{6}$ we have two sets of heavy baryons with spins $1/2$ and $3/2.$ This pattern is seen in the experimental data [30].

Mass formulas for heavy baryons read therefore as follows [13]:

Here splitting parameters $\delta_{\boldsymbol{\overline{3}}}$ and $\delta_{\boldsymbol{6}}$ are known functions of parameters $\alpha$, $\beta$ and $\gamma$ (see Eq. (9) in Ref. [15]), and $Y_{B}$ stands for a hypercharge of a given baryon.

Let us examine the consequences of the mass formulas (9). First of all, as has been already observed in [15], Eqs. (9) admit one parameter independent sum rule in the sextet

which for charm is satisfied at the level $1.4$ MeV. We use (10) to estimate $M_{\Omega_{b}^{\ast}}=6076.37$ MeV when we compute average sextet masses in the $b$ sector.

To get rid of the hyperfine splittings, we average out spin dependence in sextets by defining

Average masses $M_{B,\boldsymbol{6}}$ and masses in triplets should be equally spaced with $Y_{B}$. independently of the heavy quark.²²2We neglect small isospin violation. For $\bar{\bf 3}$ we have (in MeV)

which is satisfied with 2 % accuracy. In the case of $\boldsymbol{6}$, we have more relations (in MeV),

which are satisfied with 4 % accuracy.³³3Note that for numerical analysis in the present paper, we have used most recent version of PDG [30], and therefore, there are small numerical differences with respect to Ref. [15].

We can also form differences of average multiplet masses between the $b$ and $c$ sectors to compute the heavy quark mass difference (in MeV),

Furthermore, we can extract the hyperfine splitting parameter testing our assumptions concerning the Hamiltonian (5),

(in MeV). From these estimates, we get

with the average value of $0.283$. The PDG values of the $\overline{\rm MS}$ heavy quark masses lead to $m_{c}/m_{b}=0.3$ where both masses are evaluated at the renormalization point $\mu=m_{Q}$ [30]. Of course, heavy quark masses in the effective models, like the one considered in this paper, may differ from the QCD masses. It is therefore encouraging that we get a mass ratio close to the ratio of the QCD masses. Nevertheless, quark masses extracted from Eqs. (15) and (18),

are a bit higher (especially for $m_{b}$) than those quoted by PDG [30]. For $m_{c}/m_{b}=0.283$, we get $m_{c}=1314.1$ MeV and $m_{b}=4641.5$ MeV, which are still lower than the effective values used, e.g., in Ref. [31].

Finally, to test heavy quark dependence of the mass formulas (9), we can compute the nonstrange moment of inertia from the sextet-${\overline{\boldsymbol{3}}}$ average mass differences,

in MeV. We see that indeed heavy quark masses cancel with very high precision. This, together with Eq. (15), suggests that possible nonlinear in $m_{Q}$ binding effects are very small if not vanishing. We can therefore safely assume that formulas (15) are valid for any heavy object replacing $Q$. We pursue this possibility in the next section.

In the present case, instead of a heavy quark, we add to the soliton a heavy diquark $\bar{Q}\bar{Q}$ of spin 1. Assuming that the soliton is not changed by this replacement we arrive at the following mass formulas for tetraquarks:

where $C_{s}$ is a spin factor arising from the fact that both the sextet soliton and the diquark have spin 1:

Here $M_{B,\overline{\boldsymbol{3}}}$ is a heavy baryon mass in SU(3)_flavor $\overline{{\boldsymbol{3}}}$ and $M_{B,\boldsymbol{6}}$ is a spin averaged mass of a sextet baryon (11).

Sextet splittings satisfy the following relation

Before proceeding to numerical calculations let us discuss strong decay thresholds. Since the ground state $\boldsymbol{\overline{3}}$ tetraquarks have $J^{P}=1^{+}$, they can decay to $D+D^{\ast}$ or $B+B^{\ast}$. The corresponding thresholds are listed in the second rows of Tables 1 and 2 for nonstrange and strange tetraquarks, respectively. In the latter case, $D_{s}D^{*}$ and $B_{s}B^{*}$ thresholds are lighter than $DD_{s}^{\ast}$ or $BB_{s}^{\ast}$.

In the case of the sextet tetraquarks, we have three families of spin $0$, 1, and $2$ of nonstrange, strange and doubly strange tetraquarks. Pertinent thresholds (averaged over isospin) are listed in Tables 1–3.

Mass formulas (22) relate tetraquark masses directly to heavy baryon masses and therefore are fairly model independent. They are analogous to the masses given in Eq.(1) of Ref. [10]. The spin part has been discussed in [10] and in [32]; however, the hyperfine coupling has not been specified. Here, we know the value of $\varkappa/m_{c,b}$ (17), so in order to estimate tetraquark masses, we only need heavy diquark mass $m_{\bar{Q}\bar{Q}}$ for $m_{Q}$ in the range (19).

The main problem in predicting heavy tetraquark masses in the present model is to have a reliable estimate of the heavy diquark mass, as it is beyond the large $N_{c}$ effective theory that we have used for the light sector. To this end, we propose to apply a nonrelativistic Schrödinger equation with the Cornell potential [33]

with $\kappa=C_{F}\alpha_{\rm s}$, which has been successfully used to describe heavy $Q\bar{Q}$ spectra (see, e.g., Ref. [34]).

There are two practical reasons to use the Cornell potential in the present context. The first one is that in order to compute $QQ$ (or $\bar{Q}\bar{Q}$) masses one has to rescale model parameters by a factor of 2. This follows from the fact that the color charge $\left<{\boldsymbol{\lambda}}\cdot{\boldsymbol{\lambda}}\right>$ is factor 2 smaller when quark color charges are in an (anti)triplet than in a singlet (see, e.g., Table III in Ref. [31]). As this is quite obvious for the Coulomb term, lattice calculations suggest the same behavior of the confining part [35]. Also the chromomagnetic spin interaction, which we neglect in the following, scales in the same way.

The second reason is that the Coulomb part in potential (24) can be in fact considered as a perturbation to the linear potential, for which solutions in terms of the Airy functions are known semianalytically. We have checked that it is enough to consider only the first order perturbation theory.

We are interested in the $S$ states only, so we put $l=0$ in the pertinent Schrödinger equation. The reduced mass of the equal mass system entering the Schrödinger equation is $\mu=m_{Q}/2$. So we are looking for a solution in terms of a $u_{n}$ function defined as follows:

It is convenient to introduce a dimensionless variable $\rho$,

and rescaled dimensionless parameters $\lambda$ and $\zeta$,

With these substitutions, the Schrödinger equation takes a very simple form,

For $\lambda=0$, Eq. (28) reduces to the Airy equation, and the unperturbed energies are given in terms of the zeros $z_{n}$ of the Airy function ${\rm Ai}(\rho-\zeta)$. This follows from the boundary condition $u_{n}(0)=0$. Therefore, we have energy quantization,

Note that these zeros are negative, so the energy $\zeta_{n}^{(0)}$ is positive. The normalized solution is

First order perturbative correction is linear in $\lambda$, so the full energy reads

We need energies for two first levels only, for which $a_{1}=0.835$ and $a_{2}=0.582$. Masses of the $Q\bar{Q}$ states read

where we have introduced two new parameters,

For a given $m_{Q}$ from the range covering (19), we have computed parameters $\varepsilon_{Q}$ and $\tilde{\kappa}$ from the two lowest $Q\bar{Q}$ states.⁴⁴4Parameters $\varepsilon_{Q}$ and $\tilde{\kappa}$ must be positive. It turns out that there are no such solutions for too low $m_{Q}$. Since we need to estimate the mass of a spin 1 diquark, we have chosen as inputs $J/\psi(3096.6)$ and $\psi_{2S}(3686.1)$ for charm and $\Upsilon_{1S}(9399.0)$ and $\Upsilon_{2S}(10023.3)$ for bottom. We have checked that the original parameters $\kappa$ and $\sigma$ obtained that way are in qualitative agreement with numerical results of Ref. [34].

Having $\varepsilon_{Q}$ and $\tilde{\kappa}$ fixed, we can easily compute diquark masses in color (anti)triplet by rescaling $\kappa\rightarrow\kappa/2$ and $\sigma\rightarrow\sigma/2$, leading to $\varepsilon_{Q}\rightarrow\varepsilon_{Q}/4^{1/3}$ and $\tilde{\kappa}\rightarrow\tilde{\kappa}/4$. It is important to realize that the two terms in Eq. (32) scale differently with this change of parameters. The confining positive part is reduced by a factor $(1/4)^{1/3}\simeq 0.63$ while the Coulomb negative part is reduced by $(1/4)^{2}/3\simeq 0.4$. This delicate balance can make the diquark mass higher than the $Q\bar{Q}$ ground state. This happens, however, only at sufficiently high $m_{Q}$ where the first order perturbation theory breaks down.

The diquark masses for charm and bottom are plotted in Fig. 1. One can see that at sufficiently large mass the Coulomb term becomes equal to the confining term, and the diquark mass becomes lighter than $2m_{Q}$ signaling the break down of the first order perturbation theory. However, in the range of model masses (19), the linear confining term dominates, and the first order perturbation theory is sufficient. In Ref. [13], we have naively approximated $m_{QQ}\simeq 2m_{Q}$, whereas for the Cornell potential, we get $m_{QQ}\simeq(2.1\div 2.3)\,m_{Q}$ in the mass range (19). This seemingly small difference led to the overbinding observed in [13].

It is of course legitimate to ask how the diquark masses depend on the potential that one chooses do describe heavy quark dynamics. One could try, for example, a harmonic oscillator potential, which for $m_{c}=1400$ MeV and $\omega=590$ MeV reproduces masses of $J/\psi$ and $\psi_{2S}$. After rescaling $\omega^{2}\rightarrow\omega^{2}/2$ (which is a naive implementation of the rescaling valid for the Cornell potential), one obtains $m_{cc}=3008$ MeV, 60 MeV below the mass following from the Cornell potential for $m_{c}=1400$ MeV. Nevertheless, as we shall shortly see, this reduction does not lead to a bound tetraquark state.

Motivated by the success of the chiral quark soliton model in describing the heavy baryon spectra, we have constructed mass formulas for heavy tetraquarks with two heavy quarks of the same flavor. We first discussed baryon phenomenology to conclude that the properties of the light sector do not depend on the heavy quark properties. This is quite expected on the grounds of heavy quark symmetry. It is therefore legitimate to replace heavy quark $Q$ in color $\boldsymbol{3}$ by a heavy anti-diquark, that differs from $Q$ by mass and spin. Mass formulas (9) relate tetraquark masses to the masses of heavy baryons, and the only model parameter borrowed from the baryon phenomenology is the hyperfine splitting parameter (17) $\varkappa/m_{Q}$. In this sense, our approach, although derived from the $\chi$QSM, is fairly model independent. This is why formulas (9) are identical to the ones derived in the heavy quark limit from QCD in Ref. [10].

The only unknown ingredient of the present approach is the heavy diquark mass. To this end, we have used the Cornell potential, first to fit potential parameters to reproduce lowest spin 1 onia, both in charm and bottom sectors, and then, after rescaling these parameters, to compute the spin 1 diquark masses. We find that only bottom tetraquarks in flavor antitriplet are bound, while the charm ones are above the threshold. This is true also in the case when the structure of a heavy diquark can be resolved by the light quarks and repulsive color ${\bf 6}$ channel is included [36].

Numerical results presented in Tables 4 and 5 are in a very good agreement with the results of Ref. [10] where all necessary parameters have been extracted from data, including the mass of $\Xi_{cc}^{++}$ [14]. No model calculations have been performed in Ref. [10], and in turn, we did not use any input from doubly charmed $\Xi$. This is a strong argument in favor of our approach to the heavy diquark mass.

It is interesting to observe that our model has a completely different $N_{c}$ counting than typical models discussed, e.g., in Refs. [26] or [36], where tetraquarks are composed from four quarks for any $N_{c}$. In our case, the soliton for large $N_{c}$ belongs to a color representation $\cal{R}$ corresponding to an antisymmetric product on $N_{c}-1$ quarks. This is because we have to take one light quark from the soliton and add one heavy quark to construct a heavy baryon. For $N_{c}=2$, this is $\cal{R}={\overline{\boldsymbol{3}}}$. In order to construct a tetraquark, we need to put heavy antiquarks in a complex conjugate representation $\overline{\cal{R}}$ corresponding to $N_{c}-1$ antisymmetrized antiquarks. For $N_{c}=2$, this is ${\boldsymbol{3}}$. So for arbitrary $N_{c}$, our tetraquark consits of $N_{c}-1$ light quarks and $N_{c}-1$ heavy antiquarks, see Fig. 4. Such a configuration has been briefly discussed in Ref. [26]. A system composed of $N_{c}-1$ heavy (anti)quarks is amenable to semiclassical treatment. It would be interesting to pursue this possibility in constructing a model for a diquark.

Finally, we have to confront the LHCb tetraquark [1, 2] which is just below the $DD^{*}$ threshold. Here, two possibilities exist. Either our model is not accurate enough to deal with dynamics which gives binding energies of the order of hundreds keV, or the LHCb tetraquark corresponds to a different configuration that is out of reach for the soliton models. Obviously, charm quark mass is far from infinity and $1/m_{c}$ corrections might finally lower our predictions. In the present approach, however, we have no systematic scheme that would allow one to include such effects. Also the diquark model can be responsible for overshooting the physical mass. Nevertheless, given the very good accuracy of the $\chi$QSM predictions for heavy baryon masses and very good agreement with the phenomenological analysis of Ref. [10], one is perhaps more inclined towards the second possibility. Indeed, the LHCb [2] estimated the size of ${\cal T}_{cc}^{+}$ to be of the order of $7$ fm, significantly larger than the typical size of heavy flavor hadrons. This suggests a molecular structure of ${\cal T}_{cc}^{+}$ [37, 38].

In order to compute the space (or momentum) structure of tetraquarks in the present model, one should resort to a dynamical description of the soliton in terms of quark degrees of freedom. Some studies in this direction within $\chi$QSM have been undertaken in the case of singly heavy baryons. In Ref. [39] electromagnetic form factors and in Refs. [40, 41] gravitational form factors have been studied. It follows that heavy baryons are more compact than the proton. That conclusion should also apply to the present case, as the heavy quark or diquark is treated here merely as a static color source. The internal structure of heavy tetraquarks certainly deserves detailed studies, it is, however, beyond scope of the present paper.

This research has been supported by the Polish National Science Centre Grants No. 2017/27/B/ST2/01314 and No. 2018/31/B/ST2/01022. The author thanks the Institute for Nuclear Theory at the University of Washington for its kind hospitality, stimulating research environment and partial support by the INT’s U.S. Department of Energy Grant No. DE-FG02-00ER41132.