Spectrum of fully-heavy tetraquarks from a diquark+antidiquark perspective

Using a relativized quark model Hamiltonian, we explore the possibility that fully-heavy tetraquarks can be formed as bound-states of elementary colour-antitriplet diquarks and colour-triplet antidiquarks. Regarding ground-states in the $J^{PC}=0^{++}$ channel, the analysis reveals that narrow resonance-like structures exist near the lowest meson+meson thresholds in the following systems: $b s \bar b \bar s$, $bb \bar n \bar n$ ($n=u, d$), $bb \bar s \bar s$, $cc\bar c \bar c$, $b b \bar b \bar b$, $b c\bar b \bar c$, $b b \bar c \bar c$. We also compute extensive spectra for the fully-heavy quark flavour combinations. A reliable reaction model must be developed before a clear structural picture of any such states can be formed.


Introduction
Until the current millennium, the spectrum of known hadrons was limited to systems that fit simply into the patterns typical of constituent-quark models [1,2], i.e. quark-antiquark (qq) mesons and three-quark (qqq) baryons. Notwithstanding this, Refs. [1,2] also raised the possibility of complicated hadrons, e.g. qqqq and qqqqq. Today, a large amount of data, obtained at both e + e − and hadron colliders, has provided evidence for the possible existence of such exotic hadrons.
The first exotic discovered was the electric-charge neutral X(3872), now named χ c1 (3872) [3]. Potentially a QqQq system, where Q denotes a heavy quark, it was seen in the decay B ± → K ± X (X → J/ψπ + π − ) by the Belle Collaboration [4]. Regarding QQqqq systems, states identified as pentaquarks -P c (4312), P c (4440), P c (4457) -have recently been reported by the LHCb Collaboration in studies of the decay Λ 0 b → J/ΨK − p [5,6]. More recently, data supporting discovery of a doublecharm baryon [7]: Ξ ++ cc (3621), has focused further attention on the prospects for heavy-quark systems to reveal novel features of the Standard Model.
Owing to the large masses of the valence degrees of freedom, the possible existence of fully-heavy QQQQ bound-states (Q = c, b) and similar, mixed systems (ccbb,ccbb ∼ ccbb) can reasonably be explored using nonrelativistic tools for QCD phenomenology and theory. Here, in contrast to systems involving light-quarks, for which both light-meson and gluon exchange may play a role in tetraquark formation, binding in fullyheavy systems is very probably dominated by gluon-exchange forces because the typical gluon mass-scale (m g ∼ 0.5 GeV [64]) is much lighter than that of any necessarily-heavy meson that could be exchanged between two subsystems within the tetraquark composite. It is thus natural to suppose that the favoured structural configuration for a fully-heavy tetraquark bound-state is diquark+antidiquark.
It has been argued [65] that if stable cccc and/or bbbb tetraquarks exist, they should be observable at the Large Hadron Collider (LHC). However, the only search to date, focusing on the Υ(1S ) µ + µ − invariant-mass distribution obtained from high-energy pp collisions, was unsuccessful [66]; possibly because the width of a bbbb state is too small [67]. Experimental searches continue, motivated by theoretical analyses which predict the existence of a bbbb bound-state with mass near the η b η b (ΥΥ) threshold, e.g. Refs. [47,[68][69][70][71][72][73]. Plainly, if a stable bbbb ground-state exists, one may expect at least a few radial and orbital excitations, i.e. a spectrum of bbbb excited states.
With these motivations, we compute the spectra of cccc, ccbb,ccbb ∼ ccbb, bbbb tetraquarks from the diquark+antidiquark perspective, using a potential model characterised by linear confinement and one-gluon exchange. The Hamiltonian eigenvalue problem is solved by means of a numerical variational method based on harmonic oscillator trial wave functions, employed elsewhere for calculations of meson and baryon spectra [47,49,[74][75][76]. Using the same approach and assuming isospin symmetry, we also calculate the groundstate masses of similarly viewed bqbq, bbqq systems (q = u, s). In these cases, the justification for a diquark+antidiquark picture is weaker, but comparison of the computed masses with those of accessible colour-singlet final states can still provide hints about the possible stability of such systems.

Relativized Diquark Model
We assume that the putative tetraquark states are colourantitriplet (3 c ) diquark + colour-triplet (3 c ) antidiquark (DD) bound-states. Furthermore, the constituent D,D are each treated as being inert against internal spatial excitations [74][75][76][77]. This should be a fair approximation for fully-heavy systems owing to the suppression of quark exchange between the diquark subclusters in this case [78]. Consequently, dynamics within the DD system can be described by a single relative coordinate r rel , with conjugate momentum q rel .
To describe the internal dynamics of a D aDb system, we choose the Hamiltonian constrained elsewhere for kindred bound-states [47,49]: with the interaction being the sum of a linear-confinement term and a one-gluon exchange (OGE) potential: [25,49,79,80], where the Coulomb-like piece is [25,80] G(r rel ) = − 4α s (r rel ) 3r rel = − k 4α k 3r rel Erf(τ D aDb k r rel ) . Here, Erf is the error function and [25,80]: The parameters defining our Hamiltonian are listed in Ta 4), the values of the parameters α k and γ k (k = 1, 2, 3), σ 0 and s are drawn from Refs. [25,80]. This leaves the diquark masses; and they are all determined by using a Hamiltonian like that in Eq. (1) to solve for the mass of the given (q 1 q 2 ) sc,ax system, {q 1 , q 2 = n, s, c, b}, n = u = d, using the same constituent-quark masses employed for mesons (in GeV) [25]: M n = 0.22, M s = 0.419, M c = 1.628, M b = 4.977. Hence, the results we subsequently report are parameter-free predictions.

bbqq and bqbq ground-state masses
As an exploratory exercise, we first compute the masses of J = 0 ++ heavy-light tetraquarksbqbq, bbqq systems (q = n, s) -and compare the results with the closest meson+meson thresholds in order to obtain an indication of the possible stability of each such system.
Using the Hamiltonian specified by Eq. (1) and the parameters in Table 1, we obtain the following ground-state masses: where the energies of the two possible DD configurations are both shown, viz. scalar-scalar and axial-vector-axial-vector.
Evidently, when combining3 c and 3 c constituents, the OGE colour-hyperfine interaction favours a lighter av-av combination. This is because the spin-spin interaction in Eq. (2) is attractive. (The nature of our uncertainty estimate is discussed in Appendix A.) To gauge the possibility of stability for these systems we compare our calculated masses with prospective two-meson thresholds. 1 The lightest available final states are [3]: Hence, plausibly, bsbs tetraquark configurations may be stable. Turning now to bbqq systems, we obtain Despite a positive experiment-model mass-balance, the comparison between Eqs. (8) and (10) indicates that, in each case, our model produces a two-body final state that is lighter than the initial tetraquark; hence, bbnn and bbss tetraquarks are probably unstable.

cccc, bbbb, bcbc, bbcc ground states
In QQQQ systems treated as bound-states of colour tripletantitriplet pairs, fermion statistics also precludes a role for scalar diquarks. Consequently, the ground-state cccc is an av-av combination; and using Eq. (1) we find This value is below the empirical η c η c threshold (5.968 GeV); but a comparison with our computed value (5.82 (12) GeV) is less favourable. We conclude, therefore, that the probability of a stable cccc bound-state constituted as (cc)¯3 c (cc) 3 c is marginal. Table 2 lists our prediction for the mass of this system alongside a sample of values obtained elsewhere [47,[68][69][70][71][72][73]. Our conclusion is supported by the fact that these other analyses produce masses larger than ours. 1 Experimental masses are used here because quark models are typically not appropriate for QCD's Nambu-Goldstone bosons, especially the η-η sector.  Using the same framework, the calculated mass of the analogous bbbb system is Once again, this value is below the empirical η b η b threshold (18.797 GeV), but lies above our computed result (18.66 (13) GeV). Notably, too, our tetraquark mass is lighter than that obtained in most other analyses. It follows that one cannot confidently predict existence of a stable J = 0 ++ bbbb tetraquark. 2 Considering the bcbc case, both sc-sc and av-av configurations can exist; and we find M gs bcbc = 12.52 (08) GeV (sc-sc configuration) 12.37 (09) GeV (av-av configuration) .
The pattern observed above is repeated here. The mass of the lighter av-av configuration is (slightly) below the empirical η b η c threshold (12.383 GeV), but it lies above our computed value (12.24 (12) GeV). Given, too, that our predicted masses lie below those obtained elsewhere (see Table 2), a stable (bc)¯3 c (bc) 3 c system appears unlikely. One can also imagine J = 0 ++ bbcc (bbcc) configurations. In this case, only the av-av (bb)¯3 c (cc) 3 c configuration is possible and its ground-state mass is M gs bbcc = 12.45 (11) GeV .
The now standard pattern is evident here. Namely, our predicted mass lies below the empirical B cBc threshold (12.55 GeV) but above the model-consistent calculated value (12.36 (17) GeV). Again, therefore, a stable tetraquark in this configuration is unlikely.

Complete Tetraquark Spectra
In the preceding subsections we showed that the internally consistent application of Eq. (1) does not support stable   Table 3: Spectra obtained by solving the eigenvalue problem defined by Eq. (1). Following Appendix A, the model uncertainty in each result is 2 %. The states are labelled thus: N is the radial quantum number (N = 1 is the ground state); S D , SD are the spin of the diquark and antidiquark, respectively, coupled to the total spin of the meson, S ; the latter is coupled to the orbital angular momentum, L, to get the total angular momentum of the tetraquark, J. Degenerate states are orthogonal combinations of diquark+antidiquark spin vectors [49]. diquark¯3 c +antidiquark 3 c J = 0 ++ tetraquark systems. Notwithstanding that, these states might exist as narrow resonance-like structures above the lightest breakup threshold, but development of a reliable reaction model for tetraquark production and decays would be necessary before the character of such systems could be elucidated. We remark on this problem in Sec. 4. Neglecting decays, the Hamiltonian in Eq. (1) predicts a rich spectrum; and in Tables 3, 4 we report the lightest states in the spectra of cccc, bbbb, bcbc and bbcc systems. The typical level-ordering is illustrated using the bbbb system in Fig. 1. These results should serve as useful benchmarks for other analyses, which are necessary in order to identify model-dependent artefacts and develop a perspective on those predictions which might only be weakly sensitive to model details. Moreover, given that the decay modes of J = 0 ++ tetraquarks may be difficult to access experimentally [67], our predictions for orbitallyexcited and J 0 tetraquarks may serve useful in guiding new experimental searches for fully-heavy four-quark states.
As we have already highlighted, one source of uncertainty in our results is the choice of model Hamiltonian: Appendix A explains how we have attempted to estimate the size of this sensitivity. Another lies in the approximations used to simplify the tetraquark wave function. Within the diquark+antidiquark framework, this uncertainty arises because one can produce an overall colour singlet from both3 c × 3 c and 6 c ×6 c . Consequently, to obtain the "physical" tetraquark colour wave function, a mixing angle should be introduced: α 2 + β 2 = 1, as described, e.g. in Refs. [53,71]. Here   where quarks and antiquarks in the fundamental representations 3 c and3 c , respectively, are combined to obtain diquark (antidiquark) colour wave functions3 c , 6 c (3 c ,6 c ); and, finally, these diquark and antidiquark colour wave functions are combined into a colour singlet tetraquark configuration. If one considers systems with only a single diquark (antidiquark) as, e.g. when describing baryons as quark+diquark bound-states [82][83][84][85][86], the 6 (6 c ) is ignored because one-gluon exchange is repulsive in this channel [87]. Additionally, with diquarks (antidiquarks) treated as elementary degrees-offreedom, it is not possible to use a typical two-body Hamiltonian to determine the relative weights of the |Φ 1,3 c 3 c and |Φ 1,6 c6c components in the wave function. The mixing angle is then a free parameter, which may only be determined once substantial, reliable data becomes available. A similar problem is manifest in the spectroscopy of meson-meson molecular states where are both admissible components of the wave function. Given these issues, herein, as in other analyses, e.g. Refs. [38,47,49,67], we have only considered β = 0 in Eq. (15).

Possible Tetraquark Decay Modes
In considering the prospects for tetraquark discovery, it is important to discuss the likely decay modes. We begin with the 0 ++ fully-heavy systems. Plainly, no open-beauty final states exist for bbbb. Leptonic decays are possible, e.g. bbbb → Υ µ + µ − , and readily accessible experimentally, but estimates suggest the widths are small [67]. Regarding bbcc tetraquarks, numerous weak decays are possible and a few might be measurable [88]. Moreover, as noted above, resonance-like QQQQ systems, Q = c or b, can decay into purely hadronic final states [47,70], perhaps with an appreciable phase space. Again, however, experimental detection would likely be challenging [67]. One can also imagine the possibility of openflavour baryonic decays: QQQQ → QQq +QQq transitions, where q = u, d, s. Observation of such decay products would be a fairly unambiguous signal of a four-quark initial state; but unless one considers radial excitations of the tetraquark, the baryon-antibaryon threshold will be too high. (Computed spectra of doubly-heavy baryons are reported elsewhere, e.g. Refs. [78,[89][90][91].) With 0 ++ systems difficult to observe, it may be better to search for the J PC 0 ++ states listed in Tables 3, 4. A prime example is presented by the 1 −− systems. Possessing the same Poincaré-invariant quantum numbers as the photon, such states would be accessible via photoproduction or using e − e + colliders. Moreover, since even the lightest such states lie 300 MeV above the lowest open heavy pseudoscalar meson thresholds, there is likely sufficient phase space to enable detection. The decays could proceed as illustrated in Fig. 2; but a reliable picture of the internal structure of fully-heavy diquarks must be developed before predictions for the widths become possible. Notably, since heavy-quark exchange/rearrangement is kinematically suppressed, both the production and decay of   Table 3, obtained by solving the eigenvalue problem defined by Eq. (1).
fully-heavy tetraquarks will be difficult to observe.

Summary and Perspective
Adopting a perspective in which tetraquarks are viewed as bound-states of elementary colour-antitriplet diquarks and colour-triplet antidiquarks and using a well-constrained model Hamiltonian, built with relativistic kinetic energies, a onegluon exchange potential and linear confinement [Sec. 2], we computed the masses of ground-state bbqq, bbqq tetraquarks, q = n, s (n = u = d) and extensive spectra for cccc, bbbb, bcbc, bbcc states [Sec. 3 and Tables 3, 4]. The eigenvalue problems were solved using a numerical variational procedure in concert with harmonic-oscillator trial wave functions.
In each channel, comparing our prediction for the mass of the J PC = 0 ++ ground-state with the experimental value of the lowest meson-meson threshold, we found tetraquarks marginally stable against strong decays in almost all channels, viz. S = {bsbs, bbnn, bbss, cccc, bbbb, bcbc, bbcc}. The bnbn system lies above the η b η threshold. On the other hand, when compared with meson thresholds computed using the same Hamiltonian, all ground-state tetraquarks are marginally unstable. We therefore judge that narrow resonance-like tetraquark structures might exist near the lowest meson+meson thresholds in those channels contained in S.
Our analysis can be improved, most notably by forgoing the elementary diquark approximation and solving a four-body problem in which the internal structure of diquark correlations is resolved, e.g. using methods such as those in Refs. [53,92]. One might also tackle tetraquark systems using few body methods in quantum field theory, following Refs. [93,94]. It is perhaps most important, however, to emphasise that no clear picture of putative heavy-tetraquark states can be drawn before a reliable reaction model is developed to describe their production and decay. There is a pressing need for progress in this direction, which can yield estimates of production cross-sections and principal decay modes. We use a model Hamiltonian, Eq. (1), to compute tetraquark masses. Although constrained by an array of applications, it is still a model; hence, there is a model uncertainty. In order to provide an estimate of its size, we also computed tetraquark masses using the relativised quark model (RQM) Hamiltonian introduced in Ref. [25]. Only a few obvious changes are necessary because this Hamiltonian was also constructed to bind a colour triplet-antitriplet pair into a colour-singlet system.
When forming a S -wave system from two axial-vector constituents, the only contribution from spin-dependent interactions in the RQM is that produced by the contact term, V cont : (A.1) Contrarily, in the case of tensor, V tens , and spin-orbit, V so , interactions, one obtains the matrix elements S L J | V tens (r) |S L J = 0 0 0| V tens (r) |0 0 0 where Y (2) is a L = 2 spherical harmonic [95], and S L J | V so (r) |S L J = 0 0 0| V so (r) |0 0 0 ∝ √ L(L + 1)(2L + 1) = 0 . (A. 3) The smearing function coefficient employed in Ref. [25], σ C 1 C 2 , with C 1,2 denoting the constituents, is the same as that we use, given by Eq. (4).
As an illustrative example, consider the fully-b J = 0 tetraquark. Our prediction for the ground-state mass is reported in Eq. (12). Using the RQM Hamiltonian and a computed value of σ (bb)(bb) = 77.9 fm −1 , one finds Our mass prediction cannot be judged more accurate than the difference between this result and that in Eq. (12), viz. 74 MeV.
We therefore list this value as the uncertainty in Eq. (12).