Masses of Doubly Heavy Tetraquarks with Error Bars

In the heavy-quark limit, the two heavy quarks in a doubly heavy baryon or a doubly heavy tetraquark are bound by their color-Coulomb potential into a compact diquark. The doubly heavy hadrons are related by the approximate heavy-quark--diquark symmetry of QCD to the heavy hadrons obtained by replacing the heavy diquark by a heavy antiquark. Effective field theories can be used to expand the masses of singly heavy hadrons and doubly heavy hadrons in inverse powers of the heavy quark masses. The coefficients in the expansions for doubly heavy tetraquarks can be determined from those for heavy mesons, heavy baryons, and doubly heavy baryons using heavy-quark--diquark symmetry. We predict the masses of the ground-state doubly heavy tetraquarks with error bars using as inputs the masses of heavy mesons and heavy baryons measured in experiments and the masses of doubly heavy baryons calculated using lattice QCD. The only doubly heavy tetraquarks predicted to be stable with respect to strong decays are $bb$ tetraquarks with light flavor $\bar u \bar d$, $\bar s \bar u$ and $\bar s \bar d$.


I. INTRODUCTION
One of the most basic properties of a quantum field theory is its particle spectrum. The stable particles are particularly important, because they can exist as asymptotic states. The spectrum of quantum chromodynamics (QCD) consists of color-singlet clusters of quarks and gluons, including 2-quark (qq) mesons and 3-quark (qqq) baryons. The flavor symmetries of QCD guarantee that the lightest 2-quark meson with any net flavor and the lightest 3-quark baryon with any flavor are stable with respect to QCD interactions, so they can decay only through electromagnetic or weak interactions. The spectrum of QCD may also include exotic hadrons with additional constituents, such as tetraquark (qqqq) mesons and pentaquark (qqqqq) baryons. Since the constituents of such an exotic hadron can be rearranged into two color-singlet clusters, it can decay into two hadrons unless the mass of the exotic hadron is below the strong-decay threshold. Whether there are any such exotic hadrons that are stable is a dynamical question that depends on the parameters of QCD and specifically on the quark masses.
The possibility that the spectrum of QCD includes stable tetraquark mesons containing two heavy quarks (QQqq) was first studied by Richard and collaborators using quark potential models [1,2]. Manohar and Wise pointed out that QCD predicts that there must be stable QQqq tetraquark mesons in the limit of infinite heavy-quark mass [3]. In this limit, the attractive color-Coulomb potential between the two heavy quarks decreases the energy of the QQqq meson to below the energy of two separated Qq mesons. The color-Coulomb potential binds the QQ pair into a compact diquark whose effect on the lighter QCD fields is the same as that of a single heavy antiquark. The resulting approximate symmetry of QCD, which is called heavy-quark-diquark symmetry, relates doubly heavy QQq baryons to heavyQq mesons [4][5][6] and doubly heavy QQqq tetraquarks to heavyQqq anti-baryons [7].
Whether QCD predicts the existence of stable doubly heavy tetraquarks containing bottom or charm quarks is a dynamical question that can only be answered using quantitative information from QCD. Bicudo et al. used lattice QCD calculations of static potentials for two heavy quarks together with the Born-Oppenheimer approximation to present evidence for the existence of a stable tetraquark with flavor bbūd [8,9]. A convincing case that QCD predicts the existence of a stable bbūd tetraquark was made by Karliner and Rosner [10] and by Eichten and Quigg [11]. Karliner and Rosner also predicted that there are bc and cc tetraquarks with masses near their strong-decay thresholds. Eichten and Quigg also predicted that there are stable bbsū and bbsd tetraquarks, but that all bc and cc tetraquarks have masses well above their strong-decay thresholds.
Doubly heavy tetraquarks have been investigated directly using lattice QCD. Francis et al. presented strong evidence for the existence of deeply bound tetraquarks with flavors bbūd, bbsū, and bbsd [12,13]. The Hadron Spectrum Collaboration studied cc tetraquarks in lattice QCD with large unphysical u and d masses and found no evidence for bound states [14]. Junnarkar, Mathur and Padmanath verified the existence of deeply bound tetraquarks with flavors bbūd, bbsū, and bbsd [15]. Leskovec et al. calculated the mass of the ground-state bbūd tetraquark with quantum numbers 1 + with all the major systematic errors quantified [16].
Lattice QCD should eventually be able to provide definitive predictions for the masses of all the doubly heavy tetraquarks. However other methods based on QCD can still provide useful insights into the spectrum of doubly heavy tetraquarks, especially if the errors in their predictions can be quantified. In this paper, we quantify the errors in the approach used by Eichten and Quigg to demonstrate the existence of stable bb tetraquarks [11]. The weakest point in their analysis was in their values for the masses of doubly heavy baryons. To determine the masses of cc tetraquarks, they used the LHCb measurement of the mass of the double-charm baryon Ξ ++ cc , which has a well defined error bar [17]. However to determine the masses of bc and bb tetraquarks, they used predictions for the masses of bc and bb baryons from a quark-diquark potential model [18]. We instead use the masses of doubly heavy baryons calculated using lattice QCD, which have well-defined error bars. This allows us to give predictions for the masses of doubly heavy tetraquarks in the heavy-diquark limit with error bars. We also improve on the analysis of Ref. [11] in several other ways, including correcting some errors in the heavy-quark mass dependence of the coefficients in an effective Hamiltonian for doubly heavy hadrons.
In Section II, we use the measured masses of the ground-state heavy mesons and the ground-state heavy baryons to determine coefficients in the expansions of their masses in inverse powers of the heavy-quark mass. We also use the masses of ground-state heavy baryons calculated using lattice QCD to determine coefficients in the expansions. An effective field theory for doubly heavy hadrons [5,6] can be used to organize their masses into expansions in inverse powers of the heavy quark masses. In Section III, we use the masses of ground-state doubly heavy baryons calculated using lattice QCD to determine coefficients in the expansions. We then use heavy-quark-diquark symmetry to determine coefficients in the expansions of the masses for doubly heavy tetraquarks. Finally we use these expansions to predict the masses of doubly heavy tetraquarks with error bars. The only doubly heavy tetraquarks with masses below the strong-decay thresholds have flavors bbūd, bbsū, and bbsd. We conclude in Section IV by summarizing our results and discussing the prospects for more accurate predictions of the masses for doubly heavy tetraquarks.

II. SINGLY HEAVY HADRONS
In this section, we consider mesons and baryons that contain a single heavy quark. Heavy Quark Effective Theory (HQET) can be used to organize their masses into expansions in inverse powers of the heavy quark mass m Q . We determine coefficients in the expansions of the masses of heavy mesons and heavy baryons using the masses measured in experiments. We also determine coefficients in the expansions of the masses of heavy baryons using masses calculated using lattice QCD.

A. Heavy Mesons
A heavy hadron contains a single heavy quark Q and light QCD fields. The flavor of the heavy quark can be charm (c) or bottom (b). It has spin quantum number 1 2 , and it is in a color-triplet (3) state. In the presence of a heavy quark Q, the light QCD fields can have finite energy only if they can combine with the heavy quark to form a color singlet. This constrains the flavor of the light QCD fields. The simplest possibility is the flavor of a light antiquarkq. The corresponding hadron is a heavy meson with flavor Qq. In the presence of the heavy quark, the lowest-energy states of the light QCD fields are discrete states ℓ that can be specified by several quantum numbers: • three possible flavorsq. They can be arranged into an isospin doublet (ū,d) and an isospin singlets.
• a principal quantum number ν = 1, 2, . . . that labels successive states with the same j P .
In a heavy meson, the light QCD states with the lowest principal quantum number can be labeled by ℓ =q, j P . HQET is an effective field theory for the sector of QCD with a single heavy hadron [19]. The Lagrangian for HQET is organized into an expansion in powers of the inverse heavy quark mass 1/m Q . At first order in 1/m Q , the terms in the Lagrangian are a kinetic-energy term and a spin-dependent term that depends on the spin S of the heavy quark. From the form of the HQET Lagrangian, we can infer the form of the Hamiltonian H Q ℓ for a heavy meson with light QCD fields in the state ℓ. Through first order in 1/m Q , the only operators in the Hamiltonian are the heavy-quark spin S and the total angular momentum j ℓ of the light QCD fields. The quantum number for j ℓ is the half-odd integer j. The Hamiltonian H Q ℓ for a heavy meson through first order in 1/m Q is The only dependence on the heavy flavor Q is through the mass m Q . The coefficients E ℓ , K ℓ , and S ℓ depend on the state ℓ of the light QCD fields. The total angular momentum of the heavy meson is J = S + j ℓ . The quantum number for J is the spin J of the heavy meson. The heavy mesons with light QCD state ℓ =q, j P form a spin doublet consisting of two states with spins J = j − 1 2 , j + 1 2 whose mass splitting is proportional to 1/m Q . That spin splitting is referred to as a hyperfine splitting.
The lowest-energy states ℓ of the light QCD fields can be deduced from the masses of the observed heavy mesons containing a single charm or bottom quark. The ground state of the light QCD fields for a given light flavorq has j P = 1 2 − . It can be interpreted as a constituent antiquark in an S-wave orbital. The ground-state heavy meson with light flavorq is therefore a spin doublet consisting of two states with J P = 0 − , 1 − . The charmmeson doublets with light flavorsū,d, ands are (D 0 , D * 0 ), (D + , D * + ), and (D + s , D * + s ). The bottom-meson doublets with light flavorsū,d, ands are (B − , B * − ), (B 0 ,B * 0 ), and (B 0 s ,B * 0 s ).

B. Heavy Baryons
The second simplest possibility for the flavor of the light QCD fields in the presence of a heavy quark Q is the flavor qq ′ of two light quarks. The corresponding hadron is a heavy baryon with flavor Qqq ′ . In the presence of the heavy quark, the lowest-energy states of the light QCD fields are discrete states ℓ that can be specified by several quantum numbers: • the angular-momentum/parity quantum numbers j P , where j = 0, 1, 2, . . . is an integer and P = +, −, • a principal quantum number ν = 1, 2, . . . that labels successive states with the same j P .
In a heavy baryon, the light QCD states ℓ with the lowest principal quantum number can be labeled by [q q ′ ], j P or {q q ′ }, j P . From the form of the HQET Lagrangian, we can infer the form of the Hamiltonian H Q ℓ for a heavy baryon with light QCD fields in the state ℓ. Through first order in 1/m Q , the only operators in the Hamiltonian are the heavy-quark spin S and the total angular momentum j ℓ of the light QCD fields. The quantum number for j ℓ is the integer j. The Hamiltonian H Q ℓ for a heavy baryon through first order in 1/m Q has the same form as Eq. (1). The total angular momentum of the heavy baryon is J = S + j ℓ . The quantum number for J is the spin J of the heavy baryon. If j = 0, the heavy baryon is a spin singlet with spin J = 1 2 . If j ≥ 1, the heavy baryons form a spin doublet consisting of two states with spins J = j − 1 2 , j + 1 2 whose mass splitting is proportional to 1/m Q . That spin splitting is referred to as a hyperfine splitting.
The lowest-energy states ℓ of the light QCD fields can be deduced from the masses of the observed heavy baryons containing a single charm or bottom quark. For an antisymmetric light flavor [q q ′ ], the ground state of the light QCD fields has j P = 0 + . It can be interpreted as two constituent quarks that are antisymmetric in their colors, antisymmetric in their spins, and in S-wave orbitals. The ground-state heavy baryon with antisymmetric light flavor [q q ′ ] is a spin singlet with J P = 1 2 + . For heavy flavor Q, the baryon with the isospinsinglet light flavor [u d] is Λ Q . The baryons with the isospin-doublet light flavors ([d s], [u s]) are Ξ Q , with a superscript that specifies the electric charge. For a symmetric light flavor {q q ′ }, the ground state of the light QCD fields has j P = 1 + . It can be interpreted as two constituent quarks that are antisymmetric in their colors, symmetric in their spins, and in S-wave orbitals. The ground-state heavy baryons with symmetric flavor {q q ′ } are spin doublets consisting of two states with J P = 1 2 + , 3 2 + . For a heavy flavor Q, the spin doublets with the isospin-triplet light flavors (dd, {u d}, uu) are (Σ Q , Σ * Q ). The spin doublets with the isospin-doublet light flavors ({d s}, {u s}) are (Ξ ′ Q , Ξ * Q ). The spin doublet with the light flavor ss is (Ω Q , Ω * Q ). For each spin doublet of baryons, the members of an isospin multiplet are distinguished by giving their electric charge as a superscript.

C. Heavy Quark Masses
In Ref. [11], Eichten and Quigg presented an updated analysis of the 1/m Q expansion for the masses of ground-state heavy hadrons. Their prescription for the heavy quark mass m Q was half the mass of the most deeply bound quarkonium state with J P C = 1 −− , which is J/ψ for charmonium and Υ for bottomonium. The resulting masses for the charm and bottom quarks are m c = 1.55 GeV and m b = 4.73 GeV. This prescription is not ideal, because it is biased by the large binding energy of Υ and by the large spin splitting between J/ψ and η c . These biases can be avoided by determining m Q from the masses of heavy mesons and heavy baryons. An alternative prescription for m Q is the difference between the sum of the spin centroids of the masses for the doublets of ground-state mesons with flavors Qū and Qd and the mass of the ground-state baryon with flavor Q[u d]. The resulting masses for the charm and bottom quarks are The masses of B * 0 and B * + have not been measured separately, so they have been set equal to the measured mass M B * . Using the spin centroids of the heavy meson masses reduces the bias from the larger hyperfine splittings of the charm mesons. Contributions to the energies of the light QCD fields that can be interpreted as constituent masses of theū andd in the heavy mesons and the u and d in the heavy baryon are canceled by the subtractions in Eq. (2). If the heavy hadron masses are expanded to first order in 1/m Q using Eq. (1), our prescription for m Q implies the constraints Thus the energies E ℓ are not strictly independent of the heavy quark mass, but include some resummation of higher orders in 1/m Q . For the lowest-energy states ℓ of the light QCD fields, the coefficients E ℓ , K ℓ , and S ℓ in the Hamiltonian in Eq. (1) can be deduced from the masses of observed heavy hadrons containing a single charm or bottom quark. The Hamiltonian in Eq. (1) does not take into account the contributions to the masses from electromagnetic interactions between quarks. The electromagnetic interactions and the u−d mass difference both contribute to the isospin splitting between two heavy hadrons that differ only by the replacement of a u quark by a d quark. The resulting isospin splittings are at most about 5 MeV. In our analysis, we will not take isospin splittings into account in the Hamiltonian in Eq. (1). We do however take them into account through theoretical errors on hadron masses. We will try to take into account all contributions to hadron masses that may be larger than the isospin splittings, but we will ignore contributions that are much smaller than 10 MeV.
Hyperfine splittings between heavy hadrons come from the S · j ℓ term in the Hamiltonian in Eq. (1). The hyperfine splittings between charm mesons are roughly 140 MeV, and the hyperfine splittings between bottom mesons are smaller by about a factor of 3. A correction to the coefficient S ℓ of order 1/m Q could therefore give contributions to charm hadron masses that are greater than 10 MeV, and it should therefore be taken into account. These corrections can be taken into account by replacing S ℓ by a different coefficient S ℓ,Q for the two heavy flavors Q = c, b. It is also convenient to absorb the kinetic energy term in the Hamiltonian in Eq. (1) into the coefficient E ℓ by allowing it to be different for the two heavy flavors Q = c, b: With these two changes, our Hamiltonian for heavy hadrons with heavy quark Q and light QCD fields in the state ℓ reduces to An analysis based on this Hamiltonian should be able to predict the masses of heavy hadrons with an accuracy of about 5 MeV.
If the coefficients in the Hamiltonian in Eq. (5) have been determined for Q = c and b, they can be extrapolated in m Q to obtain estimates for the coefficients in the Hamiltonian in Eq. (1). The asymptotic coefficient of the kinetic energy term is The asymptotic spin-splitting coefficient is The difference between the coefficients E ℓ,Q in Eq. (5) for states ℓ whose light flavor differs only by the replacement of u or d by s has a contribution from the mass of the strange quark. We denote that difference by m s,Q , with the dependence on the heavy flavor explicit but the dependence on the other quantum numbers in ℓ suppressed. For states ℓ that only have flavors lighter than s, such asū, [ud], or uu, we sometimes also suppress the dependence of E ℓ,Q on the other quantum numbers, denoting it by E u/d,Q . We can interpret m s,Q as a constituent mass for the strange quark, although it also includes contributions from the difference between the kinetic energies of s and u or d. The mass splittings between strange and nonstrange hadrons are roughly 100 MeV for both charm hadrons and bottom hadrons. The difference between m s,c and m s,b could be larger than 10 MeV, and it should therefore be taken into account. The dependence of m s,c and m s,b on the suppressed quantum numbers for the light QCD state ℓ may also be significant.

D. Coefficients for Heavy Mesons from Experiment
We proceed to apply the Hamiltonian in Eq. (5) to the ground-state heavy mesons. We use the masses of heavy mesons in the 2020 Review of Particle Physics [20]. The masses have been measured for all the charm mesons in the spin doublets (D 0 , D * 0 ), (D + , D * + ), and (D + s , D * + s ) and for all the bottom mesons in the spin doublets (B − , B * − ), (B 0 ,B * 0 ), and (B 0 s ,B * 0 s ), except that the masses of B * − andB * 0 are not resolved separately. We ignore the effects of electromagnetism and the u − d mass difference in the Hamiltonian in Eq. (5). We instead treat their contributions to the masses as theoretical errors. We use isospin splittings to estimate those theoretical errors. We obtain the theoretical error in the charm-meson masses from the measured D + − D 0 mass difference and from the difference between the measured D * + and D * 0 masses. Their weighted average is 4.8 MeV. We take the theoretical error in the bottom-meson masses to be theB 0 − B − mass difference, which is 0.3 MeV. The much smaller isospin splittings for bottom mesons can be attributed to a near cancellation between the contributions to the mass from the Coulomb energy between the quark and antiquark and from the d − u mass difference. For charm mesons, these two contributions have the same sign. We obtain the total error in each heavy meson mass by adding the experimental error and the theoretical error linearly.
In a ground-state heavy meson, the light QCD fields are in the state ℓ =q, 1 2 − . We determine the coefficients E u/d,c , m s,c , and S c in Eq. (5) by minimizing the χ 2 for the masses of the 6 ground-state charm mesons. We determine the coefficients E u/d,b , m s,b , and S b by minimizing the χ 2 for the 5 measured masses of the ground-state bottom mesons. The results are given in Table I Table I are those for the Hamiltonian in Eq. (5). Extrapolations in 1/m Q can be used to estimate the coefficients in the Hamiltonian in Eq. (1) for asymptotically large m Q . The coefficient K ℓ of the kinetic energy term can be estimated by inserting E u/d,c and E u/d,b into Eq. (6). The resulting estimate for the kinetic energy K ℓ /(2m c ) in charm mesons with the lightest flavorsū andd is about 11 MeV. By replacing E ℓ,c − E ℓ,b in Eq. (6) by m s,c − m s,b , the additional kinetic energy in charm mesons with the light flavors is estimated to be about 22 MeV. These kinetic energies are not small compared to isospin splittings, so it is important to take them into account. By inserting S c and S b into Eq. (7), the asymptotic coefficient S ℓ from Eq. (7) is estimated to be (0.447 ± 0.008) GeV 2 . The contribution to the hyperfine splittings of charm mesons from the difference between S ℓ,c and S ℓ is (S ℓ,c − S ℓ )/2m c , which is about 8 MeV. This is not small compared to the isospin splittings, so it is important to take it into account.

E. Coefficients for Heavy Baryons from Experiment
We proceed to apply the Hamiltonian in Eq. (5) to the ground-state heavy baryons. We use the masses of heavy baryons in the 2020 Review of Particle Physics [20]. The masses for heavy baryons with antisymmetric light flavor have been measured for the charm baryons Λ + c , Ξ + c , and Ξ 0 c and the corresponding bottom baryons Λ 0 b , Ξ 0 b , and Ξ − b . The masses for heavy baryons with symmetric light flavor have been measured for both members of all six charmbaryon doublets: . They have been measured for both members of three of the six bottom-baryon doublets: . They have also been measured for Ξ * 0 b and Ω − b . We ignore the effects of electromagnetism and the u − d mass difference in the Hamiltonian in Eq. (5). We instead treat their contributions to the masses as theoretical errors. We use isospin splittings to estimate those theoretical errors. There are some isospin splittings for heavy baryons that have been measured directly with a better precision than the differences between the measured masses. There are 7 pairs of ground-state charm baryons that differ by the replacement of a u quark by a d quark. We take the theoretical error in the masses for the charm baryons to be the weighted average of the absolute values of these 7 isospin splittings, which is 1.7 MeV. The isospin splittings for the ground-state bottom baryons that are known are those  by a d quark, and those for Σ which differ by the replacement of two u quarks by two d quarks. We take the theoretical error in the masses for the bottom baryons to be the weighted average of these 4 isospin splittings, which is 4.9 MeV. We obtain the total error in each heavy baryon mass by adding the theoretical error and the experimental error linearly.
In a ground-state heavy baryon with antisymmetric light flavor [q q ′ ], the light QCD fields are in the state ℓ = [q q ′ ], 0 + . We determine the coefficients E u/d,c and m s,c in Eq. (5) by minimizing the χ 2 for the masses of the 3 charm baryons Λ + c , Ξ + c , and Ξ 0 c . We determine the coefficients E u/d,b and m s,b by minimizing the χ 2 for the masses of the 3 bottom baryons The results are given in Table II. In a ground-state heavy baryon with symmetric light flavor {q q ′ }, the light QCD fields are in the state ℓ = {q q ′ }, 1 + . We determine the coefficients E u/d,c , m s,c , and S c by minimizing the χ 2 for the masses of the 12 baryons in the six charm-baryon doublets. We determine the coefficients E u/d,b , m s,b , and S b by minimizing the χ 2 for the masses of the 8 baryons in the six bottom-baryon doublets whose masses have been measured. The results are given in Table II. The energies E u/d,Q in Table II for heavy baryons with antisymmetric light flavor and symmetric light flavor are larger than those for heavy mesons in Table I by about 310 MeV and 520 MeV, respectively. The strange-quark masses m s,Q in Table II for heavy baryons with antisymmetric light flavor and symmetric light flavor are larger than those for heavy mesons in Table I by about 85 MeV and 25 MeV, respectively. The spin-splitting coefficients S Q for heavy baryons with symmetric light flavor in Table II are smaller than those for heavy mesons in Table I by about a factor of 1/3. Some of the differences between the coefficients for Q = c and b in Table II correspond to energy differences that are larger than isospin splittings. The energies E u/d,c are larger than E u/d,b by about 15 MeV. The strange-quark masses m s,c are larger than m s,b by about 7 MeV. The spin-splitting coefficients S c and S b are equal to within the errors.

F. Coefficients for Heavy Baryons from Lattice QCD
The coefficients in the Hamiltonians for heavy hadrons in Eq. (5) can also be determined from the masses of hadrons calculated using lattice QCD. There have been several calculations of the masses of baryons containing heavy quarks using lattice QCD with all the major sources of systematic uncertainties under control. The systematic errors include those from the extrapolations to zero lattice spacing, to infinite volume, and to the physical light quark masses. The b quark is always treated using lattice NRQCD, which introduces additional systematic errors. Briceno et al. [21], Alexandrou et al. [22], and Brown et al. [23] have calculated the masses of ground-state baryons containing one, two or three c quarks with dynamical light quarks and with the systematic errors quantified. There are also calculations of these masses with only a single lattice spacing [24][25][26][27][28]. Brown et al. have calculated the masses of the ground-state baryons containing one, two, or three b quarks, including those containing bc, bcc, and bbc, with dynamical light quarks and with the systematic errors quantified [23]. There is a previous calculation of the masses of ground-state baryons containing one or two b quarks at a single lattice spacing [29]. There is also a recent calculation of the masses of ground-state baryons containing a single b quark, including those containing bc and bcc, with only statistical errors [30]. In all these lattice QCD calculations, the effects of electromagnetism and the u − d mass difference were ignored. We proceed to apply the Hamiltonian in Eq. (5) to the ground-state heavy baryons using masses calculated using lattice QCD. We include only those calculations in which all the important systematic errors have been quantified. For charm baryons, we use the masses calculated in Refs. [21][22][23]. The charm baryons with antisymmetric light flavor are the two spin singlets Λ + c and Ξ c . The charm baryons with symmetric light flavor are the three spin doublets (Σ c , Σ * c ), (Ξ ′ c , Ξ * c ), and (Ω 0 c , Ω * 0 c ). We have suppressed the superscript with the electric charge on members of isospin doublets and isospin triplets. For bottom baryons, we use the masses calculated by Brown et al. [23]. Some mass differences, including all the hyperfine splittings for the 3 charm-baryon doublets and the 3 bottom-baryon doublets, were also calculated in Ref. [23]. This is particularly important for bottom baryons, because the errors in the mass differences are smaller than the errors obtained by subtracting the masses by about a factor of 10. For some charm baryons, the errors in the mass differences are smaller by about a factor of 4. We obtain a single error for each mass or mass difference calculated using lattice QCD by adding the statistical and systematic errors linearly.
The only term in the Hamiltonian in Eq. (5) that contributes to the hyperfine splittings within a heavy-baryon spin multiplet is the spin-dependent term with coefficient S ℓ,Q . We determine S ℓ,Q for ℓ = {q q ′ }, 1 + by using the hyperfine splittings in heavy-baryon doublets calculated in Ref. [23]. The values of S c and S b are determined by minimizing the χ 2 for the hyperfine splittings in the 3 charm-baryon doublets and the 3 bottom-baryon doublets, respectively. The results are given in Table III.
We determine the coefficients E u/d,Q and m s,Q in the Hamiltonian in Eq. (5) for heavy baryons by minimizing the χ 2 for the masses of heavy baryons calculated using lattice QCD. For charm baryons with antisymmetric light flavor, we determine E u/d,c and m s,c by minimizing the χ 2 for the 6 values of the masses for Λ + c and Ξ c calculated in Refs. [21][22][23]. For bottom baryons with antisymmetric light flavor, we determine E u/d,b and m s,b by minimizing the χ 2 for the masses of Λ 0 b and Ξ b calculated in Ref. [23]. The results are given in Table III. For heavy baryons with symmetric light flavor, we replace the coefficient S ℓ,Q in Eq. (5) by the central value of S c or S b in Table III. We ignore the errors in S c or S b , because the errors in the hyperfine splittings used to determined these coefficients are small compared to the errors in the masses. For charm baryons with symmetric light flavor, we determine E u/d,c and m s,c by minimizing the χ 2 for the 18 values of the 6 masses for charm baryons in the spin doublets calculated in Refs. [21][22][23]. For bottom baryons with symmetric light flavor, we determine E u/d,b and m s,b by minimizing the χ 2 for the masses of the 6 bottom baryons in the spin doublets calculated in Ref. [23]. The results are given in Table III.
The coefficients E u/d,Q , m s,Q , and S Q in Table III are the same for Q = c and b to within the errors. The errors for c are smaller than those for b, because the χ 2 is calculated using many more charm baryon masses. The coefficients for c Table III are consistent within their errors with those in Table II that were determined from measured baryon masses. For c {q q ′ }, the differences between the coefficients in Table III and Table II are more significant. The coefficient E u/d,c in Table III is larger by 1.8 error bars. The coefficient m s,c in Table III is smaller by 2.1 error bars. The coefficient S c in Table III is larger by 1.2 error bars.

III. DOUBLY HEAVY HADRONS
In this section, we consider baryons and mesons that contain two heavy quarks. An effective field theory for doubly heavy hadrons [5,6] can be used to organize their masses into expansions in inverse powers of the heavy quark masses. We determine the coefficients in the expansions of the masses of doubly heavy baryons using masses calculated using lattice QCD. We then use the heavy-quark-diquark symmetry to determine coefficients in the expansions of the masses for doubly heavy tetraquarks. Finally we use these expansions to predict the masses of doubly heavy tetraquarks. This section follows the strategy of Eichten and Quigg in Ref. [11], but with some corrections and improvements.

A. Heavy Diquark
In a doubly heavy hadron, the heavy diquark has the flavor Q 1 Q 2 of two heavy quarks. Since the two heavy quarks are color triplets (3), the color state of the heavy diquark is a linear combination of anti-triplet (3 * ) and sextet (6) states. If the heavy quarks are separated by a small distance r, they have a large associated potential energy. If their color state is 3 * , the color-Coulomb potential at short distances is −2α s /3r, which is attractive. If their color state is 6, the color-Coulomb potential at short distances is +α s /3r, which is repulsive. We assume the attractive color-Coulomb potential binds the two heavy quarks into a compact diquark in the anti-triplet color state. We ignore the color-sextet component of the heavy diquark.
Given that the color state of the two heavy quarks is 3 * , they are antisymmetric in color. If the two heavy quarks have the same flavor Q, the combination of their spatial and spin states must be symmetric. The total spin S = S 1 + S 2 of the diquark is the sum of the spins of the two heavy quarks. The spin state can be antisymmetric with spin quantum number S = 0 or symmetric with S = 1. If the two heavy quarks have orbital-angular-momentum quantum number L, their spatial state is symmetric if L is even and antisymmetric if L is odd. Thus L + S must be odd. We will consider only the L = 0 state of the heavy diquark, which presumably has the lowest energy. In this case, we must have S = 1.
If the two heavy quarks have distinct flavors Q 1 and Q 2 , the combination of their flavor, spatial, and spin states must be symmetric. In the flavor-symmetric state {Q 1 , Q 2 }, L + S must be odd. In the flavor-antisymmetric state [Q 1 , Q 2 ], L+S must be even. We will consider only the L = 0 state of the heavy diquark, which presumably has the lowest energy. In this case, the flavor-antisymmetric state [Q 1 , Q 2 ] must have S = 0 and the flavor-symmetric state {Q 1 , Q 2 } must have S = 1. The flavor-antisymmetric state [Q 1 , Q 2 ] is expected to be lower in energy than the flavor-symmetric state {Q 1 , Q 2 }.

B. Doubly Heavy Baryon
In the presence of a diquark consisting of two heavy quarks Q 1 Q 2 , the light QCD fields can have finite energy only if they can combine with the two heavy quarks to form a color singlet. This constrains the flavor of the light QCD fields. If the color state of the heavy diquark is 3 * , the simplest possibility is the flavor of a light quark q. The resulting hadron is a doubly heavy baryon with flavor Q 1 Q 2 q. The lowest-energy states of the light QCD fields are discrete states ℓ that are the charge conjugates of the discrete statesl in a heavy meson with flavor Qq. The states ℓ can be specified by the light flavor q = u, d, s, the angular-momentum/parity quantum numbers j P , where j is a half-integer, and a principal quantum number. The light QCD states with the lowest principal quantum number in a doubly heavy baryon can be labeled by ℓ = q, j P .
An effective field theory called potential NRQCD (pNRQCD) for doubly heavy baryons has been developed by Brambilla, Vairo, and Rosch (BVR) [5]. This effective field theory is analogous to pNRQCD for heavy quarkonium [31]. BVR considered the general case where the two heavy quarks could have different flavors Q 1 and Q 2 and different masses m Q 1 and m Q 2 . The effective field theory is formulated in terms of a triplet field T for the heavy diquark in a 3 * color state and a sextet field Σ for the heavy diquark in a 6 color state together with the usual QCD fields for gluons and light quarks. The Lagrangian can be expanded in powers of the inverse heavy quark mass 1/m Q and the radius r of the multipole expansion. In the case of distinct heavy quarks Q 1 and Q 2 , the expansion in powers of the inverse heavy quark mass is actually an expansion in powers of 1/m Q 1 , 1/m Q 2 , and 1/(m Q 1 + m Q 2 ), which we will represent collectively by 1/m Q . BVR wrote down the terms in the Lagrangian explicitly at leading order, at first order in 1/m Q , and at first order in r. They determined the coefficients in the effective Lagrangian by matching Green functions with those of nonrelativistic QCD (NRQCD). There are three terms in the Lagrangian at first order in 1/m Q that depend on the triplet field T but not on the sextet field Σ: is the reduced mass of the two heavy quarks, • a spin term T † (S 1 /m Q 1 + S 2 /m Q 2 ) · B T that depends on the spins S 1 and S 2 of the two heavy quarks and on the chromomagnetic field B, • an orbital-angular-momentum term T † L · B T that depends on the relative orbital angular momentum L of the two heavy quarks and has a coefficient proportional to ]. The effective field theory pNRQCD for doubly heavy baryons was also considered by Fleming and Mehen (FM) [6]. FM considered only the case of identical heavy flavors QQ. They derived terms in the Lagrangian for the triplet field T by starting from NRQCD and making a sequence of field transformations. FM determined explicitly the T † (S 1 + S 2 ) · B T /m Q term in the Lagrangian. Their result for its coefficient agrees with that of BVR, but it differs by a factor of 2 from some previous results.
The effective field theory pNRQCD for doubly heavy baryons developed by BVR [5] can be used to organize the masses of doubly heavy baryons into expansions in powers of inverse heavy quark masses. We assume the terms in the pNRQCD Lagrangian involving the sextet field Σ can be ignored. From the form of the pNRQCD Lagrangian for the triplet field T and the light QCD fields, which includes the three terms itemized above, we can infer the form of the Hamiltonian H Q 1 Q 2 ℓ for doubly heavy baryons with light QCD fields in the state ℓ. Through first order in 1/m Q , the only operators in the Hamiltonian are the heavy-quark spins S 1 and S 2 , the relative orbital angular momentum L of the two heavy quarks, and the total angular momentum j ℓ of the light QCD fields. The quantum number for j ℓ is the half-integer j. Since we only consider the L = 0 state of the heavy diquark, we can drop the term in the Hamiltonian involving L. The resulting Hamiltonian for doubly heavy baryons to first order in 1/m Q is By heavy-quark-diquark symmetry, the coefficients E ℓ , K ℓ , and S ℓ in Eq. (8) are identical to those in the Hamiltonian H Q ℓ in Eq. (1) for heavy mesons. The Hamiltonian in Eq. (8) depends on the heavy flavor Q 1 Q 2 through the masses m Q 1 and m Q 2 and the diquark energy E Q 1 Q 2 . In Ref. [11], Eichten and Quigg assumed the kinetic term could be obtained from that of a heavy hadron by replacing the heavy quark mass m Q by the total mass m Q 1 + m Q 2 of the diquark. The denominator of their kinetic term was therefore 2(m Q 1 + m Q 2 ) instead of 2µ Q 1 Q 2 . In the case of equal-mass heavy quarks, the correct kinetic energy term is larger by a factor of 4. In the case of the unequal quark masses for b and c in Eqs. (2), the correct kinetic energy term is larger by a factor of 5.35. The prefactor of 1/2 in the spin-dependent term in Eq. (8) takes into account that the two heavy quarks are in an anti-triplet (3 * ) color state. This factor was first derived using effective field theory in Refs. [5,6]. The total angular momentum of the doubly heavy baryon is J = S + j ℓ , where S = S 1 + S 2 . The quantum number for J is the spin J of the doubly heavy baryon. If the quantum number S for the total spin of the heavy diquark is 0, the doubly heavy baryon is a spin singlet with spin J = j. If S is 1, the doubly heavy baryons form a spin doublet with spins J = 1 2 , 3 2 if j = 1 2 , and they form a spin triplet with spins J = j −1, j, j +1 if j ≥ 3 2 . We proceed to simplify the spin-dependent term in the Hamiltonian in Eq. (8). If the spins of the heavy quark are in a state with definite total spin quantum number S, the Wigner-Eckart theorem can be used to reduce the sum of the two interaction terms to a single term proportional to S · j ℓ . The result from the Wigner-Eckart theorem is most conveniently expressed in terms of expectation values of operators in a state of nonzero quantum number S for the total spin of the heavy quarks: In the last step, we used S i · S = 1 2 S 2 . In Ref. [11], Eichten and Quigg assumed the spindependent term could be obtained from that of a heavy hadron in Eq. (1) by replacing the heavy quark mass m Q in the denominator by the total mass m Q 1 + m Q 2 of the diquark. They did not take into account the prefactor of 1/2 in the spin-dependent term in Eq. (8), which comes from the two heavy quarks being in the 3 * color state. The denominator of their spin-dependent term was therefore 2(m Q 1 + m Q 2 ) instead of 8µ Q 1 Q 2 as in Eq. (9). There is a fortuitous agreement between these two denominators in the case of equal-mass heavy quarks. In the case of the unequal quark masses for b and c in Eqs. (2), the correct spin-dependent energy is larger by a factor of 1.34.
The spin-dependent term in the Hamiltonian for heavy hadrons is large enough that corrections to the spin-splitting coefficient S ℓ of order 1/m Q can give contributions to masses that are larger than isospin splittings. We can allow for dependence of the spin-splitting coefficient for a doubly heavy baryon on the heavy flavors Q 1 Q 2 by replacing S ℓ in Eq. (9) by S ℓ,Q 1 Q 2 . The dependence on the heavy flavors includes that from matching between pNRQCD and NRQCD at higher orders in the 1/m Q expansion. The pNRQCD Lagrangian at order 1/m Q includes the terms T † S 1 ·B T /m Q 1 and T † S 2 ·B T /m Q 2 with equal coefficients. At order 1/m 2 Q , there are corrections to the coefficients of these two terms proportional to 1/m Q 1 and 1/m Q 2 , respectively. These corrections are the same as for the corresponding term in HQET for a singly heavy hadron. They can therefore be taken into account by replacing S ℓ in the coefficients of the two spin-dependent terms in Eq. (8) by S ℓ,Q 1 and S ℓ,Q 2 , respectively. The spin-splitting coefficient S ℓ,Q 1 Q 2 in the Hamiltonian can then be deduced from the Wigner-Eckardt theorem, as in Eq. (9). The spin-splitting coefficient for cc and bb baryons are the same as for c and b mesons, respectively: The spin-splitting coefficient for a bc baryon is Note that this prediction for S ℓ,bc is between S ℓ,b and S ℓ,c . In the Hamiltonian for doubly heavy baryons in Eq. (8), it is convenient to combine the diquark energy E Q 1 Q 2 with the energies E ℓ and K ℓ /2µ Q 1 Q 2 of the light QCD fields into a coefficient E ℓ,Q 1 Q 2 that depends on the heavy flavors: We allow for dependence of the spin-splitting coefficients on the heavy flavors Q 1 Q 2 by replacing S ℓ in Eq. (9) by S ℓ,Q 1 Q 2 . With these two changes, our Hamiltonian for doubly heavy baryons in Eq. (8) reduces to the same form as that for singly heavy hadrons in Eq. (5): The quantum number S for the total spin S of the two heavy quarks depends on the flavor symmetry of the diquark. If the L = 0 diquark has the symmetric heavy flavor {Q 1 , Q 2 }, it must be in a symmetric spin state with S = 1. The ground-state for the doubly heavy baryon with a given light flavor q is therefore a spin doublet consisting of two states with J P = 1 Heavy-quark-diquark symmetry gives simple predictions for the hyperfine splittings of QQ baryons. The predictions are encapsulated in the equality between the spin-splitting coefficients S ℓ,QQ and Sl ,Q in Eq. (10). They imply, for example, Our identity between S ℓ,QQ and Sl ,Q in Eq. (10) takes into account the 1/m Q corrections in the matching between pNRQCD for doubly heavy baryons and NRQCD. However there are other contributions to the hyperfine splittings. Mehen and Mohapatra have calculated perturbative and nonperturbative corrections to the hyperfine splittings for doubly heavy baryons with identical heavy flavors [32]. The corrections can be interpreted as contributions to the spin-splitting coefficient S ℓ,QQ in Eq (13). The perturbative correction arises from an effective five-point contact operator coupling four heavy quarks and a gluon. It gives a contribution to S ℓ,QQ of order α 2 s , where α s is the strong coupling constant evaluated at the scale m Q v of the relative momentum of the two heavy quarks in the doubly heavy hadron. The nonperturbative correction arises from an expansion to next-to-next-to-leading order in the inverse heavy quark mass 1/m Q . It gives a contribution to S ℓ,QQ of order Λ 2 QCD /m 2 Q . These corrections to the spin-splitting coefficient S ℓ are taken into account in S ℓ,QQ by allowing it to depend on the heavy flavor. Their dependence on the heavy quark mass differs from that of the 1/m Q correction from matching between pNRQCD and NRQCD at order 1/m 2 Q . The results of Ref. [32] suggest that the deviations of S ℓ,Q 1 Q 2 from the asymptotic spin-splitting coefficient Sl in a heavy meson may have no simple dependence on the heavy quark masses.

C. Doubly Heavy Tetraquark
The second simplest possibility for the flavor of the light QCD fields in the presence of a heavy diquark in the color state 3 * is the flavorqq ′ of two light antiquarks. The resulting hadron is a doubly heavy tetraquark with flavor Q 1 Q 2qq ′ . The lowest-energy states of the light QCD fields are discrete states ℓ that are the charge conjugates of the discrete statesl in a heavy baryon with flavor Qqq ′ . The states ℓ can be specified by the light flavor, which can be antisymmetric [qq ′ ] or symmetric {qq ′ }, the angular-momentum/parity quantum numbers j P , where j is an integer, and a principal quantum number. The light QCD states ℓ with the lowest principal quantum number in a doubly heavy tetraquark can be labeled by [qq ′ ], j P or {qq ′ }, j P .
The effective field theory pNRQCD for doubly heavy baryons developed by Brambilla, Vairo, and Rosch in Ref. [5] applies equally well to doubly heavy tetraquarks. The effective field theory can be used to organize the masses of doubly heavy tetraquarks into expansions in powers of inverse heavy quark masses. From the form of the pNRQCD Lagrangian for the triplet field T and the light QCD fields, we can infer the form of the Hamiltonian H Q 1 Q 2 ℓ for doubly heavy tetraquarks with light QCD fields in the state ℓ. Through first order in 1/m Q , the only operators in the Hamiltonian are the heavy-quark spins S 1 and S 2 , the relative orbital angular momentum L of the heavy quarks, and the total angular momentum j ℓ of the light QCD fields. The quantum number for j ℓ is the integer j. We consider only the L = 0 state of the heavy diquark, so the terms involving L can be omitted. In this case, the Hamiltonian H Q 1 Q 2 ℓ for doubly heavy tetraquarks reduces to Eq. (8). By heavy-quarkdiquark symmetry, the coefficients E ℓ , K ℓ , and S ℓ are identical to those in the Hamiltonian H Q ℓ in Eq. (1) for heavy baryons. We can allow for dependence of the coefficients of the spin-dependent terms in Eq. (8) on the heavy flavor Q 1 Q 2 by replacing S ℓ in Eq. (9) by S ℓ,Q 1 Q 2 . It is also convenient to combine the diquark energy E Q 1 Q 2 with the energies E ℓ and K ℓ /2µ Q 1 Q 2 of the light QCD fields into a coefficient E ℓ,Q 1 Q 2 that depends on the heavy flavors, as in Eq. (12). The resulting Hamiltonian for doubly heavy tetraquarks reduces to Eq. (13).
The quantum number S for the total spin S of the two heavy quarks depends on the flavor symmetry of the heavy diquark. If the L = 0 diquark has the antisymmetric heavy flavor [b c], it must be in a antisymmetric spin state with S = 0. The doubly heavy tetraquarks with light QCD states ℓ = [qq ′ ], 0 + and ℓ = {qq ′ }, 1 + are spin singlets with J P = 0 + and J P = 1 + , respectively. If the L = 0 diquark has the symmetric heavy flavor cc, {b c}, or bb, it must be in a symmetric spin state with S = 1. The doubly heavy tetraquark with light QCD state ℓ = [qq ′ ], 0 + is a spin singlet with J P = 1 + . The doubly heavy tetraquarks with light QCD states ℓ = {qq ′ }, 1 + are a spin triplet with J P = 0 + , 1 + , 2 + .
Using heavy-quark-diquark symmetry, the coefficients in the Hamiltonian in Eq. (13) for doubly heavy tetraquarks can be deduced from the corresponding coefficients for doubly heavy baryons together with coefficients in the Hamiltonian in Eq. (5) for heavy baryons and for heavy mesons. The spin-splitting coefficient S ℓ,Q 1 Q 2 for doubly heavy tetraquarks with light QCD fields ℓ can be determined to first order in 1/m Q from the spin-splitting coefficients S ℓ,c and S ℓ,b for heavy antibaryons. The spin-splitting coefficients for cc, bb, and bc tetraquarks are given in Eqs. (10) and (11). The energy E ℓ,Q 1 Q 2 for a Q 1 Q 2 tetraquark with light QCD fields ℓ can be determined from the energy E ℓ ′ ,Q 1 Q 2 for a Q 1 Q 2 baryon with light QCD fields ℓ ′ together with the energies E ℓ,c and E ℓ,b for heavy antibaryons and the energies E ℓ ′ ,c and E ℓ ′ ,b for heavy mesons. The energies E ℓ,Q 1 Q 2 for cc, bb, and bc tetraquarks through first order in 1/m Q are The diquark energy E Q 1 Q 2 cancels between E ℓ,Q 1 Q 2 and E ℓ ′ ,Q 1 Q 2 . The coefficients of E ℓ and E ℓ ′ on the left sides and right sides match. The coefficients of K ℓ and K ℓ ′ on the left sides and right sides match upon using the asymptotic expressions for these coefficients in Eq. (6). The strange-quark masses m s,ℓ,Q 1 Q 2 for Q 1 Q 2 tetraquarks are differences between energies in Eq. (12) obtained by replacing a u or d quark by a s quark. Since the heavy-diquark energy cancels in the difference, m s,ℓ,Q 1 Q 2 can be determined from strange-quark masses in heavy anti-baryons without the subtractions in Eq. (15). The strange-quark masses for cc, bb, and bc tetraquarks through first order in 1/m Q are

D. Coefficients for Doubly Heavy Baryons from Lattice QCD
The coefficients in the Hamiltonian in Eq. (13) for doubly heavy baryons can be determined from the masses of doubly heavy baryons calculated using lattice QCD. There have been several calculations of the masses of doubly heavy baryons using lattice QCD with all the major sources of systematic uncertainties under control. The systematic errors include those from the extrapolations to zero lattice spacing, to infinite volume, and to the physical light quark masses. In bc and bb baryons, there are additional systematic errors from using lattice NRQCD for the b quark. Briceno et al. [21], Alexandrou et al. [22], and Brown et al. [23] have calculated the masses of all the ground-state cc baryons with dynamical light quarks and with the systematic errors quantified. There are also calculations of these masses with only a single lattice spacing [25][26][27][28]. Mathur and Padmanath have calculated the masses of Ω cc and Ω * cc with the systematic errors quantified [33]. Brown et al. have calculated the masses of all the ground-state bb and bc baryons with dynamical light quarks and with the systematic errors quantified [23]. There is a previous calculation of the masses of ground-state bb baryons at a single lattice spacing [29]. There is also a recent calculation of the masses of the ground-state bc baryons with only statistical errors [30]. In all these lattice QCD calculations, the effects of electromagnetism and the u − d mass difference were ignored.
We proceed to apply the Hamiltonian in Eq. (13) to the ground-state doubly heavy baryons using masses calculated using lattice QCD. We include only those calculations in   which all the important systematic errors have been quantified. For cc baryons, we use the masses calculated in Refs. [21][22][23]33]. For bc and bb baryons, we use the masses calculated by Brown et al. [23]. Some mass differences, including all the hyperfine splittings within the two doublets of cc baryons, the two doublets of bc baryons, and the two doublets of bb baryons, were also calculated in Ref. [23]. Mathur and Padmanath calculated the hyperfine splitting between Ω * cc and Ω cc [33]. The calculations of the hyperfine splittings is important, because the errors in the mass differences are 5 to 10 times smaller than the errors that would be obtained by subtracting the masses.
The only term in the Hamiltonian for doubly heavy baryons in Eq. (13) that contributes to the hyperfine splittings within a spin multiplet is the spin-dependent term with coefficient S ℓ,Q 1 Q 2 . We determine S ℓ,Q 1 Q 2 for {Q 1 , Q 2 } by using hyperfine splittings in doubly heavybaryon doublets calculated using lattice QCD. The value of S cc is determined by minimizing the χ 2 for the 3 results for the hyperfine splittings in the two doublets (Ξ cc , Ξ * cc ) and (Ω cc , Ω * cc ) calculated in Refs. [23] and [33]. The value of S bb is determined by minimizing the χ 2 for the hyperfine splittings in the two doublets (Ξ bb , Ξ * bb ) and (Ω bb , Ω * bb ) calculated in Ref. [23]. The value of S {bc} is determined by minimizing the χ 2 for the hyperfine splittings in the two doublets (Ξ ′ bc , Ξ * bc ) and (Ω ′ bc , Ω * bc ) calculated in Ref. [23]. The results are given in Table IV. We determine the coefficients E u/d,Q 1 Q 2 and m s,Q 1 Q 2 in the Hamiltonian for doubly heavy baryons in Eq. (13) by minimizing the χ 2 for the difference between the predicted masses and the masses of doubly heavy baryons calculated using lattice QCD. The statistical and systematic errors in each mass calculated using lattice QCD are added linearly. In the Hamiltonian for heavy baryons with symmetric heavy flavor, we replace the coefficient S ℓ,Q 1 Q 2 by the central value of S cc , S {bc} , or S bb in Table IV. We ignore the errors in S ℓ,Q 1 Q 2 , because the errors in the hyperfine splittings used to determine these coefficients are small compared to the errors in the masses. In the Hamiltonian for cc baryons, we determine E u/d,cc and m s,cc by minimizing the χ 2 for the 14 values of the four cc baryon masses calculated in Refs. [21][22][23]33]. In the Hamiltonians for [bc], {b, c} and bb baryons, we determine the coefficients E u/d,Q 1 Q 2 and m s,Q 1 Q 2 by minimizing the χ 2 for the 2, 4, and 4 masses calculated in Ref. [23], respectively. The results are all given in Table IV. Our results for the energies of light QCD fields in the state ℓ = q, 1 2 + in Table IV for doubly heavy baryons can be compared to those of light QCD fields in the charge-conjugate statē ℓ =q, 1 2 − in Table I for heavy mesons. The energies E u/d,Q 1 Q 2 and E u/d,Q cannot be compared, because E u/d,Q 1 Q 2 includes the diquark energy E Q 1 Q 2 . The values of the strange-quark mass in Table IV are consistent with m s in Table I averaged Table IV is consistent within the errors with S b in Table I, in accord with Eq. (10). The value of S cc in Table IV is smaller than that of S c in Table I by about 4 error bars. This raises a question about the accuracy of the compact heavy-diquark approximation for cc baryons. The value of S {b c} in Table IV is only about one half of S cc , which is incompatible with the prediction in Eq. (11) that S {b c} is between S cc and S bb . The value of S {b c} is smaller than the prediction in Eq. (11) by 6.1 error bars. It is however consistent with S b for a heavy baryon in Table II. This suggests that it may be more appropriate to regard the charm quark in a bc baryon as a light quark instead of as a constituent of a bc diquark.
The first definitive discovery of a double-charm baryon was that of the Ξ ++ cc by the LHCb collaboration in 2017 [17]. Its mass is measured to be (3621.2 ± 0.7) MeV [17,34]. (We ignore the previous observations of double-charm baryons with smaller masses by the SELEX collaboration [35,36].) We can combine the 3 results on the Ξ ++ cc mass from lattice QCD in Refs. [21][22][23] to give a postdiction for the mass. If the statistical and systematic errors in each mass are added linearly, the predicted mass from averaging the 3 results is (3585 ± 25) MeV. This is 1.4 error bars below the mass measured by LHCb. Alternatively, the mass of the Ξ ++ cc can be predicted using the Hamiltonian in Eq. (13). The expression for the mass is Upon inserting the values of m c from Eq. (2a) and E u/d,cc and S cc from Table IV, we get the prediction (3585 ± 12) MeV, which is about 3 error bars below the measured mass. The diquark energy E Q 1 Q 2 in the doubly heavy baryon can be estimated using the energy E u/d,Q 1 Q 2 in Table IV and the energies E u/d,Q in heavy mesons from Table I. These energies are defined in Eqs. (12) and (4), respectively. The energy of a diquark with identical heavy flavors can be expressed up to errors of order 1/m 2 Q as The coefficients of E ℓ on the right side cancel. The coefficients of K ℓ on the right side cancel upon using the asymptotic expression for these coefficients in Eq. (6). The predicted energy of the cc diquark is E cc = (−5 ± 12) MeV. The predicted energy of the bb diquark is E bb = (−158 ± 25) MeV. The energy of a bc diquark up to errors of order 1/m 2 Q can be expressed as The coefficients of E ℓ and K ℓ on the right side cancel. The predicted energy of a bc diquark is E [bc] = (−41 ± 37) MeV in the case of antisymmetric heavy flavor and E {bc} = (−8 ± 27) MeV in the case of symmetric heavy flavor. A sufficiently negative value of E Q 1 Q 2 implies the binding of the two heavy quarks into a diquark. The bb diquark is the only one for which our estimate indicates a significant binding energy. It is interesting to compare lattice QCD results for hyperfine splittings with the predictions of heavy-quark-diquark symmetry. The prediction for the hyperfine splitting between Ξ * ++ cc and Ξ ++ cc from heavy-quark-diquark symmetry in Eq. (14) is 107 MeV. The result from lattice QCD in Ref. [23], with statistical and systematic errors added linearly, is 83±13 MeV, which is lower by about 2 error bars. The prediction for the hyperfine splitting between Ξ * ++ bb and Ξ ++ bb from heavy-quark-diquark symmetry analogous to Eq. (14) is 34 MeV. The result from lattice QCD in Ref. [23], with statistical and systematic errors added linearly, is 35 ± 10 MeV, which is consistent within the error. Our expression for the spin-splitting coefficient S ℓ,bc for bc baryons in Eq. (11) implies the inequalities S ℓ,b < S ℓ,bc < S ℓ,c . They imply corresponding inequalities between the hyperfine splittings of b mesons, bc baryons, and c mesons. They imply, for example, The hyperfine splitting between Ξ * + bc and Ξ ′+ bc is predicted to be between 34 MeV and 107 MeV. The hyperfine splitting calculated using lattice QCD in Ref. [23], with statistical and systematic errors added linearly, is 27 ± 12 MeV. The central value is below the lower bound in Eq. (20) but by less than an error bar. If the charm quark is treated as a light quark instead of as a constituent of the bc diquark, the prediction for the spin-splitting coefficient is S ℓ,bc = S ℓ,b , where S ℓ,b is the spin-splitting coefficient for the ground-state b baryons. The corresponding prediction for the hyperfine splitting between Ξ * + bc and Ξ ′+ bc is This hyperfine splitting is predicted to be 20 MeV, which is consistent with the lattice QCD prediction 27 ± 12 MeV.

E. Predictions for Doubly Heavy Tetraquarks
The masses of the ground-state doubly heavy tetraquarks can be predicted using the Hamiltonian in Eq. (13). Our choices for the heavy quark masses m c and m b are given in Eqs. (2). We determine the coefficients in the Hamiltonian for doubly heavy tetraquarks in Eq. (13) from the coefficients for doubly heavy baryons in Table IV, the coefficients for  heavy baryons in Table II, and the coefficients for heavy mesons in Table I. The energies E u/d,Q 1 Q 2 for doubly heavy tetraquarks in Table V are determined from E u/d,Q 1 Q 2 for doubly heavy baryons in Table IV, E u/d,c and E u/d,b for heavy baryons in Table II, and E u/d,c and E u/d,b for heavy mesons in Table I by using Eqs. (15). The strange-quark masses m s,Q 1 Q 2 for doubly heavy tetraquarks in Table V are determined from the strange-quark masses m s,c and m s,b for heavy baryons in Table II by using Eqs. (16). The spin-splitting coefficients S Q 1 Q 2 for doubly heavy tetraquarks in Table V are determined from the coefficients S b and S c for heavy baryons in Table II using Eqs. (10) and Eq. (11). The results are all given in Table V.
The masses of the ground-state doubly heavy tetraquarks can be predicted by inserting the coefficients in Table V into the Hamiltonian in Eq. (13). The resulting predictions for the masses of cc and bb tetraquarks are given in Table VI. The resulting predictions for the masses of bc tetraquarks are given in The predictions for the masses of doubly heavy tetraquarks by Eichten and Quigg (EQ) in Ref. [11] are also given in Tables VI and VII. They did not give any error bars on their   Tables I,  II, and IV using Eqs. (15) and (16). The spin-splitting coefficients S cc and S bb are equal to S c and S b in Table II   Only one member of any isospin multiplet is given. All masses are in MeV. The results labeled "Eichten-Quigg" are from Ref. [11]. The results labeled "this work" were obtained using the Hamiltonian in Eq. (13) with the coefficients in Table V. The last column is the strong-decay threshold. The bold-faced energies are below the strong-decay threshold.
predictions. The predictions of EQ for the masses of bb tetraquarks agree with our predictions to within our error bars.  was (3882 ± 12) MeV, which is 7 MeV above the strong-decay threshold. Our prediction in Table VII is higher by 81 MeV, which is about 6 of our error bars. There have been several calculations of the masses of doubly heavy tetraquarks using lattice QCD. Francis et al. presented strong evidence for the existence of deeply bound tetraquark states with flavors bbūd, bbsū, and bbsd [12]. Their simulations were extrapolated to the physical values of the light quark masses, but they were carried out at a single lattice spacing and volume. They used the correlators of two local operators to determine the masses. They also provided evidence for the existence of bound tetraquark states with flavor bcūd [13]. Junnarkar, Mathur and Padmanath verified the existence of deeply bound tetraquark states with flavors bbūd and bbsū [15]. They also found that the ccūd and ccsū tetraquark states were below but close to the relevant meson pair thresholds. Their simulations were carried out at three lattice spacings but only at a single volume, and they were extrapolated to the physical values of the light quark masses. They used the correlators of two local operators to determine the masses. Leskovec et al. calculated the mass of the ground-state bbūd tetraquark with quantum numbers 1 + and quantified all the major systematic errors [16]. Their simulations were carried out using five lattice gauge ensembles with various lattice spacings and various light quark masses, including the physical values, and they were extrapolated to infinite volume. They used the correlators of three local operators and two bilocal operators to determine the mass. Their result for the mass is 10476 ± 24 ± 10 MeV. This result is in excellent agreement with our prediction for the bb[ūd] tetraquark in Table VI, which has comparable errors.

IV. DISCUSSION
We have presented predictions for the masses of doubly heavy tetraquarks with error bars. We followed the general strategy that Eichten and Quigg used to provide convincing evidence that there are stable bb tetraquarks with masses below their strong-decay thresholds [11]. Our analysis was based on the Hamiltonian for doubly heavy hadrons in Eq. (8), which provides an expansion for their masses to first order in the inverse heavy quark masses 1/m Q . The analogous Hamiltonian for singly heavy hadrons is given in Eq. (1). The approximate heavy-quark-diquark symmetry of QCD relates the coefficients in the Hamiltonians for doubly heavy baryons and doubly heavy tetraquarks to those for heavy mesons and heavy baryons, respectively. We allowed for dependence of the coefficients of the spin-dependent terms in the Hamiltonians on the heavy flavors. The resulting Hamiltonians for heavy hadrons and for doubly heavy hadrons are given in Eqs. (5) and (13), respectively. The coefficients in the Hamiltonians for the ground-state heavy mesons and for the groundstate heavy baryons were determined from measured hadron masses and are given in Tables I  and II, respectively. The coefficients in the Hamiltonian for the ground-state doubly heavy baryons were determined from baryon masses calculated using lattice QCD and are given in Table IV. The coefficients in Tables I, II, and IV were used to determine the coefficients in the Hamiltonian for the ground-state doubly heavy tetraquarks, which are given with error bars in Table V. Those coefficients were then used to predict the masses of the ground-state cc and bb tetraquarks in Table VI and the masses of the ground-state bc tetraquarks in Table VII. Our Hamiltonian for doubly heavy hadrons in Eq. (8) was deduced from the Lagrangian for an effective field theory for doubly heavy hadrons developed by Brambilla, Vairo, and Rosch (BVR) [5] and by Fleming and Mehen [6]. The terms of order 1/m Q in the BVR Lagrangian that involve the triplet diquark field T and light QCD fields were used to deduce corresponding terms in the Hamiltonian in Eq. (8). We ignored the terms in the BVR Lagrangian that involve the sextet diquark field Σ. The coefficients in the Hamiltonian for doubly heavy hadrons in Eq. (8) have factors E ℓ , K ℓ , and S ℓ that depend on the discrete state ℓ of the light QCD fields. Those same factors appear in the Hamiltonian for singly heavy hadrons. The dependence of the coefficients on the heavy quark masses m c and m b is determined by the BVR Lagrangian.
Our analysis differs from that of Eichten and Quigg (EQ) in Ref. [11] in several important ways. Instead of setting the heavy quark masses m c and m b equal to half the masses of the quarkonium states J/ψ and Υ, respectively, we set them equal to the linear combinations of the heavy meson masses and the heavy baryon mass in Eqs. (2). We corrected errors in the coefficients of the kinetic and spin-dependent terms in the Hamiltonian for doubly heavy hadrons in Ref. [11]. EQ assumed the denominator of the coefficient of the kinetic term is twice the total mass m Q 1 + m Q 2 of the diquark instead of twice the reduced mass. The correct coefficent increases the size of the kinetic term in the Hamiltonian for doubly heavy hadrons by a factor of 4 for cc and bb hadrons and by a factor of 5.35 for bc hadrons. EQ assumed the denominator of the coefficient of the spin-dependent term in the Hamiltonian was also proportional to m Q 1 + m Q 2 . Their coefficient is fortuitously correct in the case of identical heavy quarks, but the correct coefficient is larger by a factor of 1.34 in the case of bc hadrons. We also improved upon the analysis of EQ by using expressions for the energies E ℓ,Q 1 Q 2 and the spin-splitting coefficients S ℓ,Q 1 Q 2 of doubly heavy tetraquarks that are accurate through first order in 1/m Q instead of only 0th order.
Our results are in surprisingly good agreement with those of Eichten and Quigg [11]. Most of their predicted masses agree with ours to within the errors, with the largest difference being 1.3 error bars. We agree with Eichten and Quigg that the only doubly heavy tetraquarks with masses below the strong-decay thresholds are the ground states with flavor bb[ūd], bb[sū], and bb[sd]. Our prediction for the mass of the bb[ūd] tetraquark is in very good agreement with the lattice QCD calculation in Ref. [16]. Our prediction for the mass is larger than that of Karliner and Rosner [10] by 3.5 of our error bars. Our predictions for the masses of the ground-state bc and cc tetraquarks are well above their strong-decay thresholds, in agreement with the results of Eichten and Quigg and contrary to the predictions of Karliner and Rosner [10].
Our predictions for the masses of doubly heavy tetraquarks could be made more precise by additional calculations of the masses of doubly heavy baryons using lattice QCD. Additional calculations of the masses of bc and bb baryons would be especially helpful. The only such calculations in which all the important systematic errors have been quantified are the pioneering calculations by Brown et al. [23]. In order for the lattice QCD calculations to provide the greatest possible insights, it is important that results are given not only for individual masses, but also for appropriate mass differences, as in Ref. [23]. The error bars in hyperfine splittings can be much smaller than those obtained by subtracting the masses.
Our analysis of doubly heavy baryons using masses calculated using lattice QCD raises questions about the accuracy of heavy-quark-diquark symmetry for doubly heavy hadrons containing charm quarks. The spin-splitting coefficient S bb for bb baryons in Table IV is compatible within errors with S b for bottom mesons in Table I. However S cc for cc baryons in Table IV is smaller than S c for charm mesons in Table I by several error bars. A more striking problem is that S {bc} for bc baryons in Table IV is significantly smaller than either S cc or S bb , instead of being intermediate between S cc and S bb as predicted by Eq. (11). The spin-splitting coefficient S {bc} is however compatible within errors with S b for bottom baryons in Table II. This suggests that it may be a better approximation to treat the charm quark in a bc hadron as a light quark instead of as a constituent of a bc diquark. In this case, bc baryons and bc tetraquarks would not be related by heavy-quark-diquark symmetry.
Our analysis was based on the assumption that the two heavy quarks in a doubly heavy hadron form a compact diquark in a 3 * color state. The diquark can also be in a 6 color state. In the effective field theory of Ref. [5], the color-sextet component of the heavy diquark is taken into account through terms in the BVR Lagrangian that involve the sextet field Σ. If the effects of Σ are significant, the Hamiltonian for doubly heavy hadrons would not have the simple form in Eq. (13), with coefficients related to those in the Hamiltonian for heavy hadrons by heavy-quark-diquark symmetry.
The two heavy quarks in a doubly heavy hadron can be treated as a diquark only if their separation is smaller than the length scale of the light QCD fields. It may be possible to take into account contributions from heavy quarks with larger separations using the Born-Oppenheimer approximation. The Born-Oppenheimer approximation for QCD was pioneered by Juge, Kuti, and Morningstar, who applied it to heavy quarkonium and to heavy quarkonium hybrids [37]. The application of the Born-Oppenheimer approximation to exotic heavy hadrons that contain a QQ pair and a light quark-antiquark pair was proposed in Ref. [38]. Brambilla et al. have developed an effective field theory framework based on the Born-Oppenheimer approximation [39]. Such a framework has been applied extensively to heavy quarkonium hybrids [40][41][42][43]. The Born-Oppenheimer approximation for QQ hadrons was pioneered by Bicudo et al., who applied it to doubly heavy tetraquarks [44]. Bicudo et al. used lattice QCD calculations of the static potentials for two heavy quarks together with the Born-Oppenheimer approximation to present evidence for the existence of stable bb tetraquarks [8,9].
The Born-Oppenheimer approximation suggests that deviations from the compact heavy diquark limit are very different for doubly heavy tetraquarks and doubly heavy baryons. The difference arises because the constituents of a doubly heavy tetraquark can be rearranged into two color-singlet clusters, while those of a doubly heavy baryon cannot. Forming two color-singlet clusters from the constituents of a doubly heavy baryon requires the creation of a quark-antiquark pair, and this is a process that is dynamically suppressed in low-energy QCD. The Born-Oppenheimer potentials between the two heavy quarks in a tetraquark and in a baryon therefore differ dramatically as their separation increases. For doubly heavy tetraquarks, the Born-Oppenheimer potentials display screening behavior, approaching constants near the thresholds for pairs of heavy mesons. For doubly heavy baryons, the Born-Oppenheimer potentials instead display string-breaking behavior, with avoided crossings near the thresholds for a heavy meson and a heavy baryon.
A new nonrelativistic effective field theory for doubly heavy hadrons beyond the compact heavy-diquark limit was recently developed by Soto and Tarrús Castellà [45] and applied to doubly heavy baryons [46]. At leading order in the 1/m Q expansion, the Lagrangian has heavy quark spin symmetry and corresponds to the Born-Oppenheimer approximation. At next-to-leading order in the 1/m Q expansion, the Lagrangian includes terms that depend on the total spin and the orbital angular momentum of the two heavy quarks. Soto and Tarrús Castellà calculated the spectra of cc and bb baryons at leading order in the 1/m Q expansion by solving a Schroedinger equation with three coupled channels [46] using Born-Oppenheimer potentials, extracted from lattice QCD calculations with unphysically large light-quark masses. They found that the compact heavy-diquark limit was not a good approximation for orbital-angular-momentum excitations of the heavy quark pair, even in the case of the b quark. Using the measured mass of the Ξ ++ cc as input [17,34], they predicted the hyperfine splitting between the Ξ * ++ cc and Ξ ++ cc to be 136 ± 44 MeV. This is consistent to within the errors with the prediction of 107 MeV from heavy-quark-diquark symmetry in Eq. (14).
While lattice QCD will provide definitive calculations of the masses of some doubly heavy tetraquarks, a more complete picture of their spectra can be obtained by using lattice QCD in conjunction with other theoretical methods. We have used lattice QCD calculations of masses of doubly heavy baryons in conjunction with constraints deduced from the effective field theory pNRQCD for doubly heavy hadrons to predict the masses of doubly heavy tetraquarks with error bars. A crucial assumption in our analysis is that the two heavy quarks form a compact diquark. That assumption is avoided in the effective field theory for doubly heavy hadrons developed in Refs. [45,46], which makes use of lattice QCD calculations of Born-Oppenheimer potentials. Such an effective field theory, along with the discovery of more doubly heavy hadrons in future experiments, should eventually provide a complete picture of the spectrum of doubly heavy tetraquarks.