Spectroscopy and decays of the fully-heavy tetraquarks

We discuss the possible existence of the fully-heavy tetraquarks. We calculate the ground-state energy of the $bb \bar b \bar b$ bound state, where $b$ stands for the bottom quark, in a nonrelativistic effective field theory framework with one-gluon-exchange (OGE) color Coulomb interaction, and in a relativized diquark model characterized by OGE plus a confining potential. Our analysis advocates the existence of uni-flavor heavy four-quark bound states. The ground state $bb\bar b\bar b$ tetraquark mass is predicted to be $(18.72\pm0.02)$~GeV. Mass inequality relations among the lowest $QQ\bar{Q}\bar{Q}$ state, where $Q\in \{c, b\}$, and the corresponding heavy quarkonia are presented, which give the upper limit on the mass of ground state $QQ\bar{Q}\bar{Q}$. The possible decays of the lowest $bb\bar{b}\bar{b}$ are highlighted, which might provide useful references in the search for them in ongoing LHC experiments, and its width is estimated to be a few tens of MeV.


I. INTRODUCTION
Heavy quarkonium spectroscopy is the arena of quark potential models, in which the full interaction incorporates the short-range one-gluon-exchange (OGE) and long-range confining potential [1][2][3][4][5]. To put heavy quarkonium spectroscopy in a model-independent framework, efforts have been made to develop a nonrelativistic effective field theory formalism [6][7][8].
Higher-order perturbative calculations show that the lowest-lying heavy quarkonia can be assigned as weakly coupled states [9][10][11]. They are characterized by a momentum scale mv Λ QCD , where m and v are the heavy quark mass and velocity, respectively, and Λ QCD is the nonperturbative scale of quantum chromodynamics (QCD), and their dynamics is mainly dominated by the short-distance interaction. In particular, there is some evidence that the first few bottomonia and B c mesons can be identified as weakly coupled systems [12][13][14], while the J/ψ lies on the borderline of this classification [8]. The lowest-lying doubly-and triply-heavy baryons have also been studied as weakly coupled bound states [15,16] but, due to the lack of experimental incentive, 1 this area is not as well explored as that of heavy quarkonia.
Fully-heavy tetraquarks, such as bbbb and the charm analog cccc, are of considerable interest, since they are free of light degrees of freedom and might be used as an elegant probe to investigate the interplay between both aspects of QCD: perturbative and nonperturbative. Owing to the lack of experimental incentive, only a few theoretical studies have been conducted in this field, and whether Nature allows the bound states of fully-heavy tetraquarks or not is still an open question. If it does, what are their masses and decay properties? Recently, this debate has given rise to several studies within different approaches, including the constituent quark and diquark model [39][40][41][42], chiral quark model [43], chromoelectric potential model [44], and QCD sum rules [45,46]. All of these predict the existence of fully-heavy (bottom) tetraquarks, except for Ref. [44], in which the authors argue that the stability of a fully-heavy system of quarks should rely on more subtle effects that are not included in the simple picture of constituent quarks. The fullyheavy tetraquark case is similar to that of polyelectrons Ps 2 , the bound state of two electrons and two positrons (e − e − e + e + ), discussed long ago by Wheeler [47]. Just after the prediction of Ps 2 , a debate on their stability started [48,49], but it took half a century to get experimental evidence of them [50].
In this paper, we compute the bbbb ground-state energy. The calculations are performed within the formal-ism of a nonrelativistic effective field theory (NREFT) at the leading order (LO), which neglects all spin-dependent and confining interactions, and in a relativized diquarkantidiquark model, both of which are characterized by the OGE potential. We also make a rough estimate of the bbbb total width, in order to guide experimental searches.
The paper is organized as follows. In section II, we describe the NREFT approach to solve the system considered. At the end of this section, we include our diquark model results. Section III is devoted to discussing the results, and we compare our predictions of the mass of the X bbbb ground-state with those of previous studies. In Section IV, we derive a set of mass inequalities, and provide the upper limits of the ground state masses of several fully-heavy tetraquarks. In Sec. V, we highlight some possible decay modes of the X bbbb and give the ballpark estimate of its total decay width. Finally, we provide a short summary.

II. FORMALISM
In this section, we describe an NREFT approach to the fully-heavy tetraquark spectroscopy, characterized by the OGE interaction. This is used to estimate the mass of the bbbb ground state. As all of the four quarks are very heavy, the tetraquark system is weakly coupled with a small size of the order 1/ 1 is the heavy quark velocity. In this case, the dynamics of the system is dominated by the OGE which provides a color Coulomb potential at the LO, and the long-distance confining potential and spin-dependent interactions become perturbations.
In addition, we also give a tetraquark (diquarkantidiquark) model prediction for the bbbb ground-state mass.

A. Nonrelativistic Hamiltonian
The nonrelativistic Hamiltonian describing the 2Q−2Q system takes the following form at the LO where T i = m i + p 2 i /(2m i ), r ij ≡ |r i − r j |, and the second term is the spin-independent OGE color Coulomb potential, where λ i and λ j are color matrices and α s is the strong coupling strength. The quark mass m b is to be determined from the ground state bottomonium masses, and α s is taken at the scale of the typical momentum transfer.

B. Solving the four-body problem
The four-body problem is notoriously delicate. One should be careful about the choice of the wave functions, because a crude adoption may give rise to misleading conclusions, as first illustrated by Ore in 1946 in Ref. [48] in the case of polyelectrons (e − e − e + e + bound states). One year later the author, in collaboration with Hylleraas [49], came up with an elegant prescription to handle four-body systems. This is what we are willing to use to calculate the bbbb ground-state energy, provided that we make the substitution e − → b, e + →b.
Since the ground state is non-degenerate and its symmetries are governed by those of the Hamiltonian H NR , the spatial wave function should be symmetric under the exchange of Q 1 Q 2 ↔Q 3Q4 . This leads to ψ spatial (Q 1 Q 2 ) = ψ spatial (Q 3Q4 ). This symmetry [49] helps to reduce the number of integration variables and simplify the four-body problem. To describe the quark relative motion, we define the following Jacobi coordinates, which are shown in Fig. 1. The "physical" tetraquark color wave function can be written as a superposition of different color configura-tions, 2 |ψ c = α |3 12 3 34 ; 1 1234 + β |6 12634 ; 1 1234 , (5) where the coefficients α and β can be obtained by diagonalizing the model Hamiltonian on the tetraquark wave function. When the two identical heavy quarks are in thē 3 color representation and in an S-wave, the spin must be equal to 1 so as to satisfy the Pauli principle. For the two heavy quarks being in the 6 color representation and in an S-wave, the spin needs to be 0. Therefore, the mixing between3 and 6 color representations requires a flip of heavy quark spin, and thus is of higher order in NREFT (suppressed by v 2 ∼ 0.1). Since the OGE gives an attractive interaction for3 and a repulsive interaction for 6, at LO the ground-state color configuration should be |3 12 3 34 ; 1 1234 . Such a ground-state configuration has been suggested [51,52], and more recently in Refs. [40,43,53]. It is also supported by recent lattice QCD studies on fully-heavy four quark systems [54].
On the basis of the |3 12 3 34 ; 1 1234 color configuration, the kinetic energy matrix elements can be written as where m 1 , m 2 and m 3 are the reduced masses of Q 1 Q 2 , Q 3Q4 and Q 1 Q 2 −Q 3Q4 , respectively. The spatial trial wave function is written in terms of the Jacobi coordinates of Eq. (4) and Fig. 1, where (ξ 1 , ξ 2 , ξ 3 ) ≡ (σ, ρ, λ), β i are the oscillatory (variational) parameters, and is the overall normalization constant. It has been recently shown that the above Gaussian variational basis is a powerful tool to obtain the ground state energy of heavy tetraquarks [53]. It is worth noting that the color configurations |1 13 1 24 ; 1 1234 and |8 13 8 24 ; 1 1234 are linear superpositions of the two color representations of Eq. (5).

C. bbbb ground-state
The eigenvalue problem is solved by means of a numerical variational approach with Gaussian trial wave functions [57,58]. The quark mass m b is extracted via Eq. (3), where one also has to include the spin-averaged mass of the 1S (ground-state) bottomonia, The other parameter, i.e., the running coupling constant α s (µ), is taken at the scale of the typical momentum transfer in the ground-state bottomonium (µ = 1.5 GeV). The extracted values are The spin-averaged binding energy of the 1S bottomonium is We can estimate the energy of the X bbbb ground-state by means of the variational method with Gaussian trial wave functions previously used in the bottomonium case and the model parameter values of Eq. (10). We get where the uncertainty is given by the product of the binding energy and α s (1.5 GeV). The binding energy of the lowest four b-quark bound state is obtained as The extracted optimal values of the variational parameters β i are β 1 = β 2 = 0.77 GeV and β 3 = 0.60 GeV, where β 1 = β 2 is due to the symmetry of the spatial wave function [49], as discussed in the last subsection.

D. bbbb ground state in a relativized diquark model
In this subsection, we provide an estimate of the fullybottom tetraquark ground-state energy by treating the X bbbb configuration as a diquark, bb, and antidiquark, bb, bound state in a relativized diquark model. In this model, the diquark is an effective degree of freedom which describes two strongly correlated quarks with no internal spatial excitations. The lowest-energy diquark configurations, scalar and axial-vector, are both color anti-triplets but have different spin and flavor quantum numbers. The former, i.e., the scalar diquark, has spin-0 and its flavor wave function is antisymmetric. The latter has spin-1 and is flavor-symmetric. As mentioned in the last section, because of the Pauli principle, a (point-like) diquark made up of quarks of the same flavor can only be of the axial-vector type [51,59].
To calculate the energy of the bbbb ground state in diquark anti-diquark configuration, we need an evaluation of the bb axial-vector diquark mass, M av,bb . The diquark (antidiquark) mass is estimated by binding a bb (bb) pair via the OGE plus a confining potential [2], and we have 4 M av,bb = 9.85 GeV. (14) After doing that, we can compute the 1S, 0 ++ tetraquark ground-state energy in the relativized diquark model. The model Hamiltonian is characterized by the OGE plus confining potential. For more details on the model and the values of the model parameters, we refer to Ref. [60].
With the use of the previous value of M av,bb , we get which is consistent with the result obtained in the previous subsection. See also Table I, where our relativized diquark model prediction is compared to those of previous studies. The theoretical uncertainty in the above prediction can be estimated by using the typical error of the quark model calculations. The intrinsic error in the quark model predictions is of the order of 30 ∼ 50 MeV [38]; therefore, this diquark model prediction has the cumulative uncertainty of O(50 MeV).

III. RESULTS AND DISCUSSION
In Table I, our NREFT and relativized diquark model results for the mass of the bbbb ground-state are compared with those of previous studies. The results are strongly model-dependent and vary in the range of (18.7 ± 0.2) GeV, approximately. The differences among the results mainly stem from different choices of effective Hamiltonians, the model-parameter fitting procedures, and the use of distinct approximations in the tetraquark wave function. This last is related to the possible ways of combining the quark color representations to obtain a color singlet wave function for the tetraquark. In the last column of the table, we also report the differences between the predicted masses and the η b η b threshold, (18798 ± 4) MeV [18] whose central value is used. It is worthwhile to remind some previous studies. For example, Karliner et al. [41] made a phenomenological estimation of the QQQQ ground-state energies, with Q = b and c. The authors utilized a relation between meson and baryon masses to extrapolate the binding energy of QQQQ systems in the diquark model. They made the ballpark estimates of B bb (1S) and the binding energy of the tetraquark with respect to the two-bb-diquark threshold, and obtained a value of about 1/2 for the ratio B bb (1S)/B bbbb (1S).
Bai et al. calculated the bbbb ground-state energy by means of a phenomenological potential, whose parameters were fitted to the bb spectrum and verified by lattice simulation [39]. The bbbb ground-state energy was obtained by solving the Schrödinger equation numerically; spin-dependent corrections and a linear confining potential were included. Their result agrees with our LO NREFT result very well, indicating that neglecting the spin-dependent and confining parts provides a very good approximation to the system under study.
A diquark model prediction of the bbbb and cccc ground-state energies was given in Ref. [42]. There, the authors assumed tetraquarks to be made up of two almost point-like diquarks in the color-triplet configuration. The model parameters were fitted by solving a non-relativistic Schrödinger equation for the charmonium and bottomonium spectrum. The bbbb ground-state mass agrees with our diquark model well.
In Ref. [61], the authors computed the cccc groundstate energy within a nonrelativistic potential model. The authors concluded that the possibility of obtaining bound states depends on the assumptions made on the quark dynamics and the flavor configurations. There are also early predictions for the mass [62] and width [63] of the fully-charm tetraquark, and the ground-state energies of Q 2Q2 systems, with Q ∈ {s, c, b, t} [64]. For the full mass spectrum of all-charm tetraquarks, we refer to Ref. [59]. All-charm four-quark bound state has also been studied in the Bethe-Salpeter approach and the ground state was reported to be deeply bound (650 MeV) below the 2η c threshold [65].
A study of the importance of mixing effects between tetraquark and molecular-type components is also worthwhile to be carried out, especially for those states which are close to meson-meson thresholds. The probability of a fully-heavy tetraquark to be a (QQ) − (QQ) hadronic molecule might be calculable without any ambiguity. It is worthwhile to mention that the possibility of an η bη b bound state has been studied by computing the QCD van der Vaals force in the framework of potential nonrelativistic QCD in Ref. [66]. Given that the recent lattice nonrelativistic QCD calculation did not find any state below the noninteracting η b η b threshold [67], such a study seems necessary to understand what happens.
It would also be very interesting to test the possibility of a bbcc tetraquark that remains stable against strong decays, 5 but unfortunately there is no experimental evidence yet. The stability of heavy-light tetraquark systems is still an open question. The QQqq was shown to be stable against strong decays by Lipkin long ago 6 [68]. Very recently, bbqq was shown to be stable against strong decays but not its charm counterpart ccqq, nor the mixed (beauty+charm) bcqq state [69,70]. For a detailed discussion on the stability of such systems, we refer to the recent study in Ref. [71].

IV. TETRAQUARK-MESON MASS INEQUALITIES
One can estimate the lower bounds of the ground state energy levels of the fully-heavy four-quark systems by using the variational approach, as suggested long ago by Nussinov [72] and by Bertlmann and Martin [73]. This approach was tested by the lattice QCD calculations of Weingarten [74] and also by the rigorous calculations in the vectorlike gauge theories (QCD) by Witten [75]. We attempt to extend Nussinov's approach to the fully-heavy four-quark bound states and obtain inequalities between the tetraquark states and the corresponding heavy quarkonia in the following. 5 All the other possible fully-heavy tetraquarks can decay strongly by annihilating at least a pair of quarks and antiquarks of the same flavor. 6 Before Lipkin, the authors of Ref. [61] also argued that X QQqq is stable if the quark mass ratio m Q /mq is large enough. For reasonable values of the quark mass ratio, 5 ∼ 20, they predict a X QQqq above but not very far from the corresponding threshold. Note that the ratios m b /mu ≈ 16, mc/mu ≈ 5, and m b /mc ≈ 3, if one uses typical values for the constituent quark masses. Let us consider the general Hamiltonian of four heavy quarks with pair-wise interactions: The color-antitriplet potential V (3) QQ between any quark (or antiquark) pair can be related to the color-singlet quark-antiquark potential, V (1) is the overall ratio of the color coefficient in the leading SU(N c ) group. The above relation can also be verified by using the eigenvalues of the Casimir invariants from the Table II, where as before we consider only the color-antitriplet for Q i Q j and triplet forQ iQj , and we have written the colorsinglet quark-antiquark (Q iQj ) potential from applying Eq. (17) to Q i Q j as V (1) ij . To work out the last term in the above equation, one has to calculate the eigenvalues of the color matrices, viz. the Casimir invariants whereĈ i,j = λ i,j /2. The eigenvalues of these Casimir invariants are explicitly given in Table II. We need to know in which color representation the Q iQj pair is. For instance, by writing down explicitly the color wave function (see, e.g., Ref. [77]), we can get for i = 1, 2; j = 3, 4. This result enables us to write the pair-wise potential of four-quark state in terms of the quark-antiquark color-singlet potential, viz., the quarkonium potential, 12 +V (1) 34 14 +V (1)
(22) This means that under the approximation of one-gluon exchange and that the two quarks are in color anti-triplet, the four-quark Hamiltonian [(Eq. (18)] can be expressed in the following form, where H ij = T i +T j +V (1) ij (r ij ) is the quarkonium Hamiltonian. Taking quarkonium wave functions as the trial wave function, and applying the variational principle [72], we obtain an upper bound on the ground state energy of four-heavy-quark system by computing the expectation value of the Hamiltonian of the subsystems, namely where ψ ij are the corresponding ground state wave functions of the ij subsystems. Here, we use instead of ≤ because we have made the approximation that the two quarks are in color anti-triplet and the two anti-quarks are in color triplet, though it is expected to work rather well since the mixing with the sextet-antisextet configuration is expected to be suppressed by v 2 ∼ 0.1 for fullybottom and ∼ 0.3 for fully-charm four-quark systems. In this sense, the inequalities derived for baryons by Nussinov [72] are more rigorous because the two quarks in a baryon must be in an anti-triplet without any approximation.
The inequalities for all possible fully-heavy tetraquarks are listed in Table III. The corresponding numerical values are obtained using the spin-averaged quarkonium masses M bb (1S) = 9.445 GeV, M cc (1S) = 3.069 GeV and M cb (1S) = 6.324 GeV. Since the B * c (1 3 S 1 ) meson has not been observed, for M cb (1S) we used the average value of the theoretical predictions of Godfrey-Isgur Model [78] and the Cornell potential model from Ref. [79].
If we compare the so-obtained approximate upper bound in the fully-bottom sector with the earlier theoretical predictions listed in Table I, we find this value is indeed larger than all the model predictions except for that of Ref. [46], which is calculated by means of QCD sum rules and has a large uncertainty.

V. POSSIBLE DECAYS OF X bbbb
The main difficulties in the experimental observation of fully-bottom tetraquark mesons are related to the production mechanisms and observation of the main decay modes. A very recent study discussed the production of the bbbb ground-state at the LHC, concluding that the bbbb is supposed to be very narrow and likely to be discovered [40]. Another recent study [80] also discussed the possible production of a narrow scalar resonance around 18 ∼ 19 GeV at the LHC. 7 In the following, we will discuss the decays of such a state and provide a rough estimate of its width. We will argue that the width is at least a few tens of MeV.
As shown in Table I, most of the models predict the ground-state fully-bottom tetraquark to be below the η b η b threshold. Because of this, one expects its width to be almost saturated by the following decay modes: 1) decays into final states containing a pair of bottom and anti-bottom quarks; 2) decays into hadrons of lighter flavors. The former is dominated by one-gluon exchange, while the latter should be dominated by two gluons. Figure 2 shows an example of the decay into a pair of openbottom mesons. Because the gluons are at the scale of m b , one expects the first type of decays to dominate over the second. As it is not an easy task to calculate partial decay widths into given exclusive decay modes, in the following we provide a rough estimate of the inclusive decay width of the X bbbb ground-state based on the first decay mode.
The inclusive width of the X bbbb into final states with a pair of b andb quarks can be described as a two-step process triggered by a transition operator T as where h 1 and h 2 indicate possible bottom hadrons allowed by the spin-parity quantum numbers and the available phase space. Here T 1 is operator for the coupling of a heavy quark-antiquark pair to a gluon, given by √ 4πα sQ 1 2 λ a γ µ Q a µ , and T 2 is responsible for the transition from bbg to the final hadronic states. We need to sum up all the possible final states generated at the second vertex of Fig. 2, which includes not only two-body but also many-body final states. As a result, only the second factor (T 1 ) in the above equation matters.
TheQQ annihilation at short distances requires a direct dependence on the zero-point wave function of the color-octetQQ inside the four-body bound state, R QQ (8) (0). Hence, It is well-known that the η b decay width is saturated by two-gluon exchange, and the η b inclusive decay width [81]. The wave function at origin of a color-octet bb (8) pair might be larger than that of the asymptotic color-singlet bb (1) , namely a bottomonium [82]. However, at the present stage there is no need for a precise calculation of the width. For an order-of-magnitude estimate, one may simply assume |R bb (8) (0)| 2 ∼ |R bb (1) (0)| 2 . Therefore, we can estimate the inclusive width of the X bbbb decays into hadrons containing b andb as The η b dominantly decays into hadrons, hence, . Taking α s (m b ) = 0.22 [18], and neglecting all other possible decay modes which should be subdominant, we get The previous width (of the order of a few tens of MeV) is large enough to make the resonance observable. Our estimate of Eq. (28) is of the same order of magnitude as an earlier prediction for similar systems [63], while it is much larger than the estimate of Ref. [41], 1.2 MeV. Similarly, we expect the width of the fully-charm tetraquark, if below the η c η c threshold, to be larger, The fully-bottom tetraquark states could be searched for in final states including a pair of bottom hadrons, such as BB, Λ bΛb , Ξ bΞb , Σ bΣb and Ω bΩb . They can also decay into a fully leptonic final state via an intermediate Υ(1S)X state as where l can be τ , µ or e, and X could be the off-shell lowest vector bottomonium, X ≡ Υ(1S) * . This decay involves the annihilation of two bb pairs into virtual photons, so the branching fraction is expected to be small, O(10 −4 ∼ 10 −8 ), as estimated in [41]. 8 Despite of this, the multi-lepton final states are expected to provide a clean signal with a low background. The ideal place to look for the decays of Eq. (30) is the Large Hadron Collider experiments, where the Higgs boson cross section was measured by reconstructing a four-lepton final state [83]. Because of this, a scan at relatively lower energies, of the order of 2M η b (1S) , should be almost straightforward at LHC.

VI. SUMMARY
We calculated the bbbb ground-state energy in terms of two different approximations for the tetraquark wave function. They were used to simplify the solution of the four-body problem.
In the first case, we provided an evaluation of the bbbb ground-state energy in an NREFT at the LO, where the potential is an OGE-induced color Coulomb potential. A nice feature of our approach is that the color wave function is completely given by |3 12 3 34 ; 1 1234 at the LO, and its mixing with |6 12634 ; 1 1234 only occurs at higher orders. The solution of the four-body problem was simplified by making use of the symmetries introduced by Hylleraas and Ore in their study of polyelectrons, namely the bound states of two electrons and positrons. In our specific case, one of the above mentioned symmetries could be expressed as ψ spatial (bb) = ψ spatial (bb), where ψ spatial (bb) and ψ spatial (bb) are the spatial wave functions of the bb andbb systems, respectively. Thanks to this, the four-body problem was simplified by reducing the number of integration variables. In the second case, we calculated the bbbb ground-state energy in a relativized diquark model. This model is characterized by OGE plus a confining potential. In the diquark model, the effective degree of freedom of the diquark, describing two strongly correlated quarks with no internal spatial excitations, is introduced. Tetraquark mesons are then obtained as twobody diquark-antidiquark bound states. Our results in both approaches-NREFT at LO (18.72±0.02 GeV) and relativized diquark model (18.75 GeV)-only differ by a few tens of MeV, and suggests the existence of a bbbb bound-state below the η b η b threshold.
We also derived a set of approximate inequalities for the binding energies of the fully-heavy tetraquarks in terms of those of various heavy quarkonia. Instead of giv-ing the values of the ground-state energies of the states of interest, the inequalities provide upper bounds on them. As expected, our LO NREFT and relativized diquark model results on the X bbbb ground-state mass satisfy the corresponding inequality. Finally, we discussed the possible decay modes of the lowest fully-heavy tetraquarks and estimated the decay width of the ground state bbbb to be of O(50 MeV). We hope our results might provide useful references in the search for fully-heavy tetraquarks in ongoing LHC experiments. In particular, we suggest searching for the lowest X bbbb in the relevant center-of-mass energy region around 18.7 GeV in the final states of four leptons or a pair of bottom hadrons such as BB, Λ bΛb , Ξ bΞb , Σ bΣb and Ω bΩb ,

NOTE ADDED
During the revision process of this manuscript, a preliminary CMS analysis (not yet approved by the collaboration) was presented in APS April meeting [84], hinting at a potential four-lepton excess in the mass range of 18 ∼ 19 GeV.