Axial Vector $cc$ and $bb$ Diquark Masses from QCD Laplace Sum-Rules

TL;DR

This paper predicts the masses of axial vector $cc$ and $bb$ diquarks using QCD Laplace sum-rules and employs these results to estimate tetraquark masses within a specific model.

Contribution

It provides the first NLO QCD sum-rule calculations of axial vector diquark masses including higher-dimensional condensate effects.

Findings

$cc$ diquark mass is 3.51 GeV

$bb$ diquark mass is 8.67 GeV

Tetraquark masses estimated using diquark inputs

Abstract

Constituent mass predictions for axial vector (i.e., $J^P=1^+$) $cc$ and $bb$ colour antitriplet diquarks are generated using QCD Laplace sum-rules. We calculate the diquark correlator within the operator product expansion to NLO, including terms proportional to the four- and six-dimensional gluon and six-dimensional quark condensates. The sum-rules analyses stabilize, and we find that the mass of the $cc$ diquark is 3.51~GeV and the mass of the $bb$ diquark is 8.67~GeV. Using these diquark masses as inputs, we calculate several tetraquark masses within the Type-II diquark-antidiquark tetraquark model.

Tables

Table 1. Table 1: constituent mass predictions and sum rule parameters for axial vector c ​ c 𝑐 𝑐 cc and b ​ b 𝑏 𝑏 bb diquarks. The theoretical uncertainties are obtained by varying the the QCD input parameters in Eqs. ( 29 )–( 35 ).

$QQ$	$M_P$ (GeV)	$s_0$ (GeV2)	${\tau }_{\text{min}}$ (GeV-2)	${\tau }_{\text{max}}$ (GeV-2)
$cc$	$3.51\pm 0.35$	$17.5\pm 3.4$	$0.10\pm 0.02$	$0.71\pm 0.07$
$bb$	$8.67\pm 0.69$	$80.0\pm 9.2$	$0.02\pm 0.01$	$0.21\pm 0.02$

Equations

$$j_{\mu ,\alpha }={ϵ}_{\alpha \beta \gamma }{Q}_{\beta }^{T}C{\gamma }_{\mu }{Q}_{\gamma }$$

$$j_{\mu ,\alpha }^{†}=-{ϵ}_{\alpha \beta \gamma }{\overline{Q}}_{\beta }{\gamma }_{\mu }C{\overline{Q}}_{\gamma }^T$$

$$\Pi (q^2)=\frac{i}{D-1}\left(\frac{q_{\mu }{q}_{\nu }}{q^2}-g_{\mu \nu }\right)\int d^{D}x{e}^{iq\cdot x}⟨\Omega |\tau [j_{\mu ,\alpha }(x)S_{\alpha \omega }(x,0)j_{\nu ,\omega }^{†}(0)]|\Omega ⟩$$

$$S_{\alpha\omega}(x,0)=\hat{\mathcal{P}}\;\text{exp}\bigg{[}ig_{s}\frac{\lambda^{a}_{\alpha\omega}}{2}\int_{0}^{x}dz^{\mu}A_{\mu}^{a}(z)\bigg{]}$$

$$\big{\langle}\alpha G^{2}\big{\rangle}=\alpha_{s}\big{\langle}\!\colon\!G_{\omega\phi}^{a}G_{\omega\phi}^{a}\!\colon\!\big{\rangle}$$

$$\displaystyle\big{\langle}g^{3}G^{3}\big{\rangle}=g_{s}^{3}f^{abc}\big{\langle}\!\colon\!G^{a}_{\omega\zeta}\;G^{b}_{\zeta\rho}G^{c}_{\rho\omega}\!\colon\!\big{\rangle}$$

$$\big{\langle}J^{2}\big{\rangle}=\frac{D}{6}\,\kappa\,g_{s}^{4}\big{\langle}\overline{q}q\big{\rangle}^{2}$$

$${\Pi }^{\text{OPE}}(q^2)=\sum _{\text{i}=\text{I}}^{\text{XI}}{\Pi }^{\text{(i)}}(q^2)$$

$${\Pi }^{(\mathrm{I})}(z)=\frac{4m^2}{3{\pi }^2}z(2z+1){\text{H}}_1(z)$$

$$z=\frac{q^2}{4m^2}.$$

$${\text{H}}_1(z)={}_2{F}_1\left(1,1;\frac{5}{2};z\right),$$

$${\Pi }^{(\mathrm{II})}$$

$$m^{2ϵ}\left(4z^2ϵ(ϵ+1)+z\left(8{ϵ}^3+18{ϵ}^2+13ϵ+2\right)+1\right)\Gamma (-ϵ-1)+$$

$$4z\left(3(4\pi {)}^{ϵ}(2ϵ+1)(2zϵ+1)-m^{2ϵ}(2z(ϵ+1)(2z+2ϵ-1)+2ϵ+1)\Gamma (1-ϵ)\right){\text{H}}_2(z;ϵ)-$$

$$m^{2ϵ}\left(-4z^2(ϵ+1)(3ϵ+2)+zϵ(6ϵ+5)+2(ϵ+1)\right)\Gamma (-ϵ-1){\text{H}}_3(z;ϵ)+$$

$$\displaystyle m^{2\epsilon}\left(8z^{3}\epsilon(\epsilon+1)-8z^{2}(\epsilon+1)(3\epsilon+1)+2z(\epsilon-1)(2\epsilon+1)+2\epsilon+1\right)\Gamma(-\epsilon-1)\text{H}_{4}(z;\epsilon)\bigg{]}$$

$${\Pi }^{(\mathrm{III})}$$

$$(4ϵ(ϵ+2)+3)(-z(-8z(ϵ+1)+ϵ(4ϵ(ϵ+2)+7)+2)-2ϵ-3)-$$

$$8z(ϵ+1)(2ϵ+1)(2ϵ+3)\left(4z^2(ϵ+1)+zϵ(2ϵ+3)+2{ϵ}^2+ϵ-1\right){\text{H}}_2(z;ϵ)-$$

$$4z(ϵ+1)(2ϵ+1)(4ϵ(ϵ+2)+3)\left(1-2z(ϵ+1)\left(2(z-1){ϵ}^2+zϵ+2z+ϵ\right)\right){\text{H}}_2(z;ϵ{)}^2+$$

$$(2ϵ+3)(4ϵ(ϵ+2)+3)\left(-8z^2(ϵ+1{)}^2+zϵ(2ϵ+1)+2(ϵ+1)\right){\text{H}}_3(z;ϵ)-$$

$$\displaystyle(4\epsilon(\epsilon+2)+3)\Big{(}-16z^{3}\begin{aligned} &(\epsilon+1)-8z^{2}(\epsilon+1)(\epsilon(2\epsilon+7)+2)-\\ &2z(\epsilon(4\epsilon(\epsilon+2)+5)-1)+4\epsilon(\epsilon+2)+3\Big{)}\text{H}_{4}(z;\epsilon)\bigg{]},\end{aligned}$$

$${\text{H}}_2(z;ϵ)={}_2{F}_1\left(1,-ϵ;\frac{3}{2};z\right)$$

$${\text{H}}_3(z;ϵ)={}_3{F}_2\left(1,-2ϵ-1,-ϵ;\frac{1}{2}-ϵ,ϵ+2;z\right)$$

$${\text{H}}_4(z;ϵ)={}_3{F}_2\left(1,-2ϵ,-ϵ;\frac{1}{2}-ϵ,ϵ+2;z\right).$$

$$\displaystyle\Pi^{(\mathrm{IV})}(z)=\frac{-3\left(8z^{2}-17z+6\right)+\left(2z^{2}-11z+6\right)\text{H}_{1}(z)}{288\pi m^{2}(z-1)^{3}}\big{\langle}\alpha G^{2}\big{\rangle}$$

$$\displaystyle\Pi^{(\mathrm{V})}(z)=\frac{12z-15-(2z-3)\text{H}_{1}(z)}{576\pi m^{2}(z-1)^{2}}\big{\langle}\alpha G^{2}\big{\rangle}$$

$$\displaystyle\begin{aligned} \Pi^{(\mathrm{VI})}(z)=\frac{\big{\langle}g^{3}G^{3}\big{\rangle}}{92160\pi^{2}m^{4}(z-1)^{5}z}\Big{(}416z^{5}&-1888z^{4}+3078z^{3}-1836z^{2}+90z+35+\\ &5\left(8z^{5}-36z^{4}+42z^{3}-20z^{2}+20z-7\right)\text{H}_{1}(z)\Big{)}\end{aligned}$$

$$\displaystyle\Pi^{(\mathrm{VII})}(z)=\frac{32z^{3}-89z^{2}+19z+8-\left(12z^{4}-66z^{3}+73z^{2}-37z+8\right)\text{H}_{1}(z)}{55296\pi^{2}m^{4}(z-1)^{4}z}\big{\langle}g^{3}G^{3}\big{\rangle}$$

$$\displaystyle\begin{aligned} \Pi^{(\mathrm{VIII})}(z)=\frac{\alpha_{s}^{2}\big{\langle}\overline{q}q\big{\rangle}^{2}}{4860m^{4}(z-1)^{5}z}\Big{(}3&\left(576z^{5}-2608z^{4}+4458z^{3}-3316z^{2}+765z+20\right)-\\ &5\left(8z^{5}-36z^{4}+66z^{3}-74z^{2}+3z+12\right)\text{H}_{1}(z)\Big{)}\end{aligned}$$

$$\displaystyle\Pi^{(\mathrm{IX})}(z)=\frac{160z^{3}-478z^{2}+386z+7-\left(56z^{4}-196z^{3}+226z^{2}-68z+7\right)\text{H}_{1}(z)}{1944m^{4}(z-1)^{4}z}\alpha_{s}^{2}\big{\langle}\overline{q}q\big{\rangle}^{2}$$

$$\displaystyle\Pi^{(\mathrm{X})}(z)=\frac{-8z^{3}+19z^{2}+7z-3+\left(4z^{4}-22z^{3}+23z^{2}-13z+3\right)\text{H}_{1}(z)}{55296\pi^{2}m^{4}(z-1)^{4}z}\big{\langle}g^{3}G^{3}\big{\rangle}$$

$$\displaystyle\Pi^{(\mathrm{XI})}(z)=\frac{-3\left(16z^{4}-56z^{3}+57z^{2}-z-1\right)-\left(4z^{4}-14z^{3}-5z^{2}-3z+3\right)\text{H}_{1}(z)}{27648\pi^{2}m^{4}(z-1)^{4}z}\big{\langle}g^{3}G^{3}\big{\rangle}.$$

$$\displaystyle\alpha_{s}\rightarrow\alpha_{s}(\mu)=\frac{\alpha_{s}(M_{\tau})}{1+\frac{25\alpha_{s}(M_{\tau})}{12\pi}\log\!{\Big{(}\frac{\mu^{2}}{M_{\tau}^{2}}\Big{)}}}$$

$$\displaystyle m\rightarrow m_{c}(\mu)=\overline{m}_{c}\bigg{(}\frac{\alpha_{s}(\mu)}{\alpha_{s}(\overline{m}_{c})}\bigg{)}^{12/25}$$

$$\displaystyle\alpha_{s}\rightarrow\alpha_{s}(\mu)=\frac{\alpha_{s}(M_{Z})}{1+\frac{23\alpha_{s}(M_{Z})}{12\pi}\log\!{\Big{(}\frac{\mu^{2}}{M_{Z}^{2}}\Big{)}}}$$

$$\displaystyle m\rightarrow m_{b}(\mu)=\overline{m}_{b}\bigg{(}\frac{\alpha_{s}(\mu)}{\alpha_{s}(\overline{m}_{b})}\bigg{)}^{12/23},$$

$${\alpha }_s(M_{\tau })=0.330\pm 0.014$$

$${\alpha }_s(M_Z)=0.1185\pm 0.0006$$

$${\overline{m}}_c=(1.275\pm 0.025)\text{GeV}$$

$${\overline{m}}_b=(4.18\pm 0.03)\text{GeV}.$$

$$\displaystyle\big{\langle}\alpha G^{2}\big{\rangle}=(0.075\pm 0.02)\ \text{GeV}^{4}$$

Full text

1 Introduction

Outside-the-quark-model hadrons consisting of four (or more) valence quarks have been theorized for decades. For example, the concept of tetraquarks, hadrons composed of four quarks $(qq\overline{q}\overline{q})$, was introduced in [1, 2] in 1977. Jump forward to 2003 and the discovery of the X(3872) by the Belle collaboration [3] and its subsequent confirmation by several other experimental collaborations [4, 5, 6, 7] places us in a new era of hadron spectroscopy. Since then more and more of these hadrons have been discovered in the heavy quarkonium spectra. These hadrons, now collectively referred to as the XYZ resonances, are difficult to explain within the quark model [8]. These XYZ resonances have served as a strong motivator for research into beyond-the-quark-model hadrons. See [9, 10] for a review of experimental findings and [11, 12] for a review of several multiquark systems.

Looking at four-quark states in particular, there are several interpretations of what their internal quark structure might resemble. One possibility is that there are no particularly strong correlations between any of the quarks. However, another possible interpretation is that these states could be meson-meson molecule states in which two colour-singlet mesons form a weakly bound conglomerate state. See [13, 14, 15, 16, 17, 18, 19, 20] for discussions about the X(3872) in this configuration. Yet another possible interpretation is that four-quark states are diquark-antidiquark states. Diquarks are strongly correlated, colour antitriplet pairs of quarks within a hadron. (As such, their colour configurations are identical to those of antiquarks.) See [21] for applications of diquarks and [22] for a discussion of possible diquark configurations. In a diquark-antidiquark configuration, the diquark constituents are strongly bound together in a four-quark configuration. See [23, 24, 25, 26, 27] for discussions about the X(3872) in the diquark-antidiquark configuration. Also, see [28] for additional discussions on the differences between the molecular and tetraquark models in the context of a QCD sum-rules analysis.

QCD sum-rules analyses of diquarks in several channels have been presented in [29, 30, 31, 32, 33, 34]. Lattice QCD analyses of light diquarks have also been performed [35, 36, 37]. In this paper, we use QCD Laplace sum-rules (LSRs) to calculate the constituent masses of axial vector (i.e., $J^{P}=1^{+}$) $cc$ and $bb$ diquarks. The axial vector is the only quantum number that can be realized for colour antitriplet diquarks of identical flavours in an S-wave configuration. We use the operator product expansion (OPE) [38] to compute the correlation function between a pair of diquark currents (1)–(2). In this calculation, in addition to leading-order (LO) pertubative contributions, we also include next-to-leading-order (NLO) perturbative contributions and non-perturbative corrections proportional to the four-dimensional (4d) and 6d gluon condensates as well as the 6d quark condensate. The results of these calculations are summarized in Table 1. In particular, we find that the constituent mass of the $cc$ diquark is $(3.51\pm 0.35)$ GeV and the constituent mass of the $bb$ diquark is $(8.67\pm 0.69$) GeV. Substituting these diquark constituent masses into the Type-II diquark-antidiquark tetraquark model of Ref. [39], we calculate masses of several $[cc][\bar{c}\bar{c}]$, $[cc][\bar{b}\bar{b}]$, and $[bb][\bar{b}\bar{b}]$ tetraquarks.

2 The Correlator

The axial vector, colour antitriplet diquark current is given by [30, 31]

[TABLE]

with adjoint

[TABLE]

where $C$ denotes the charge conjugation operator, $\epsilon_{\alpha\beta\gamma}$ is a Levi-Civita symbol in quark colour space, and $Q$ is a heavy (charm or bottom) quark field.

Using (2), we consider the diquark correlator

[TABLE]

where $D$ is the spacetime dimension. In (3), $S_{\alpha\omega}(x,0)$ is a path-ordered exponential, or Schwinger string, given by

[TABLE]

where $\hat{\mathcal{P}}$ is the path-ordering operator. The Schwinger string allows gauge-invariant information to be extracted from the gauge-dependent current (1) [30, 31]. The explicit cancellation of the gauge parameter has been shown for perturbative contributions up to NLO [40], and in Landau gauge the NLO contributions from the Schwinger string are zero [30, 31]; hence $S_{\alpha\omega}(x,0)\to\delta^{\alpha\omega}$. For non-perturbative contributions of QCD condensates, gauge-invariance of the correlator (3) implies that fixed-point gauge methods used to obtain OPE coefficients are equivalent to other methods [41]. As observed in Refs. [30, 31], the Schwinger string will not contribute to the QCD condensate contributions in the fixed-point gauge, and hence $S_{\alpha\omega}(x,0)\to\delta^{\alpha\omega}$. Thus, using Landau gauge for pertubative contributions and fixed-point gauge methods for QCD condensate contributions, we can simplify (3) by setting $S_{\alpha\omega}(x,0)\to\delta^{\alpha\omega}$ (as in [28]). Lattice QCD analyses of constituent light diquark masses are also based on correlation functions of (coloured) diquark operators [35, 36, 37]. Instead of the Schwinger string, gauge dependence of the correlation function is addressed in lattice analyses either through gauge fixing or coupling to a heavy colour source.

We evaluate the correlator (3) within the OPE to NLO in perturbation theory and include non-perturbative corrections proportional to the 4d and 6d gluon condensates and the 6d quark condensate. Each non-perturbative correction is the product of a LO perturbatively computed Wilson coefficient and a QCD condensate. The 4d and 6d gluon and 6d quark condensates are defined respectively by

[TABLE]

[TABLE]

[TABLE]

where $\kappa$ in (7) quantifies deviation from vacuum saturation. As in [42, 43], we set $\kappa=2$ for the remainder of this calculation, e.g., see [44] and references contained therein.

The diagrams computed in the simplification of (3) are given in Figure 1. Each diagram has a (base) multiplicity of two associated with interchanging the quark fields contracted on the top and bottom quark lines. Diagrams II, IV, VI, VIII, X, and XI receive an additional factor of two to account for vertical reflections. As noted earlier, Wilson coefficients are calculated in the Landau gauge. Divergent integrals are handled using dimensional regularization in $D=4+2\epsilon$ dimensions at $\overline{\text{MS}}$ renormalization scale $\mu$. We use a dimensionally regularized $\gamma^{5}$ satisfying $(\gamma^{5})^{2}=1$ and $\{\gamma^{\mu},\,\gamma^{5}\}=0$ [46]. The recurrence relations of Refs. [47, 48] are implemented via the Mathematica package TARCER [49] resulting in expressions phrased in terms of master integrals with known solutions including those of [50, 51].

The OPE computation of $\Pi$, denoted $\Pi^{\text{OPE}}$, can be written as

[TABLE]

where the superscript in (8) corresponds to the labels of the diagrams in Figure 1. Evaluating the first term in this sum, $\Pi^{(\mathrm{I})}$, expanding the result in $\epsilon$, and dropping a polynomial in $q^{2}$ (which does not contribute to the sum rules—see Section 3), we find

[TABLE]

where $m$ is a heavy quark mass and

[TABLE]

Also,

[TABLE]

where functions of the form ${}_{p}F_{q}\left(\cdots;\cdots;z\right)$ are generalized hypergeometric functions, e.g., [52]. Note that hypergeometric functions of the form ${}_{p}F_{p-1}(\cdots;\cdots;z)$ have a branch point at $z=1$ and a branch cut extending along the positive real semi-axis. In evaluating $\Pi^{(\mathrm{II})}(z)$, we find a nonlocal divergence which is eliminated through the inclusion of the counterterm diagram, Diagram C1, of Figure 1. From this point forward, we refer to the renormalized contribution arising from the sum of Diagrams II and C1 as $\Pi^{(\mathrm{II})}(z)$. Note that, in Landau gauge, Diagram III does not have a nonlocal divergence corresponding to the fact that the (multiplicative) vector diquark renormalization constant is trivial [40]. The Mathematica package HypExp [53] is used to generate the $\epsilon$-expansions of $\Pi^{(\mathrm{II})}(z)$ and $\Pi^{(\mathrm{III})}(z)$. These expansions are lengthy, and so we omit them for the sake of brevity; instead, we present the exact ($\epsilon$-dependent) results,

[TABLE]

[TABLE]

where

[TABLE]

The $\epsilon$-expanded results for the remaining terms in (8) can be written more concisely and are given by

[TABLE]

Finally, substituting (9), (12), (13) and (17)–(24) into (8) gives us $\Pi^{\text{OPE}}$.

Renormalization-group improvement requires that the strong coupling and quark mass be replaced by their corresponding running quantities evaluated at renormalization scale $\mu$ [54]. At one-loop in the $\overline{\text{MS}}$ renormalization scheme, for $cc$ diquarks, we have

[TABLE]

and for $bb$ diquarks,

[TABLE]

where [55]

[TABLE]

For $cc$ diquarks, $\mu\to\overline{m}_{c}$ and for $bb$ diquarks, $\mu\to\overline{m}_{b}$. Finally, the following values are used for the gluon and quark condensates [56, 57, 58]:

[TABLE]

3 QCD Laplace Sum-Rules, Analysis, and Results

We now proceed with the QCD LSRs analysis of axial vector $cc$ and $bb$ diquarks. Laplace sum-rules analysis techniques were originally introduced in [59, 60]. Subsequently, the LSRs methodology was reviewed in [61, 62].

The function $\Pi(q^{2})$ of (3) satisfies a dispersion relation

[TABLE]

for $q^{2}<0$. In (36), $t_{0}$ is an effective threshold and $\cdots$ represents a polynomial in $q^{2}$. On the left-hand side of (36), $\Pi$ is identified with $\Pi^{\text{OPE}}$ computed in Section 2. On the right-hand side of (36), we express $\mathrm{Im}\Pi(t)$, i.e., the spectral function, using a single narrow resonance plus continuum model,

[TABLE]

where $M$ is the diquark constituent mass and $h_{+}$ is the diquark coupling defined by

[TABLE]

which aligns with the notation of Ref. [31]. Also, $\theta(t)$ is a Heaviside step function and $s_{0}$ is the continuum threshold parameter. Constituent diquark masses are key input parameters of Type I & II diquark-antidiquark models of tetraquarks [23, 39] (see Section 4). However, as the couplings are not parameters of Type I & II tetraquark models, we eliminate them by working with ratios of LSRs (e.g., see (45)). Though not relevant for our purposes here, we note that knowledge of the coupling $h_{+}$ for light diquarks allows estimation of baryon matrix elements of the effective weak Hamiltonian [30, 31].

As discussed in Ref. [30], the duality relation (37) for diquarks is more subtle than for hadrons because diquarks are constituent degrees of freedom rather than hadron states. Ref. [30] argues that, similar to constituent quarks, the diquark mass and coupling should be regarded as effective quantities which describe the correlator at intermediate scales. Above the threshold $s_{0}$, the diquark loses its meaning as a constituent degree of freedom, and the correlator is dominated by the parton-level quark description (see Diagram I in Fig. 1). In the context of lattice QCD, the coupling $h_{+}$ is proportional to the signal strength, and Ref. [35] finds a remarkably clean exponential decay indicative of a single narrow resonance below the lattice cutoff $1/a^{2}$. In (37), $s_{0}$ is analogous to the lattice cutoff $1/a^{2}$. Thus, in the light quark sector studied in [35], there exists direct lattice QCD evidence supporting the spectral decomposition (37).

Laplace sum-rules are obtained by Borel transforming (36) weighted by powers of $Q^{2}$ (see [59, 60] as well as, e.g., [63, 44]). For a function such as $\Pi^{\text{OPE}}$ computed in Section 2, details on how to evaluate the Borel transform can be found in [42, 43] for instance. We find

[TABLE]

where $\mathcal{R}_{k}(\tau)$ are unsubtracted LSRs of (usually non-negative) integer order $k$ evaluated at Borel scale $\tau$ and where $\Gamma$ is the integration contour depicted in Figure 2. Subtracting the continuum contribution,

[TABLE]

from the right-hand sides of (39) and (40), we find

[TABLE]

where $\mathcal{R}_{k}(\tau,\,s_{0})$ are (continuum-)subtracted LSRs.

In (42), explicitly parametrizing each $\Gamma_{i}$ of $\Gamma$, we have

[TABLE]

which is then calculated numerically. In (44), $R$ is set to $2m^{2}$. Also, it is intended that $\delta\rightarrow 0^{+}$. In practice, this can be achieved by setting $\delta=10^{-12}\ \text{GeV}^{2}$. Finally, using (43), we find

[TABLE]

To use (45) to predict diquark constituent masses, we must first select an acceptable range of $\tau$ values, i.e., a Borel window $(\tau_{\text{min}},\,\tau_{\text{max}})$. To determine the Borel window, we follow the methodology outlined in [28]. To generate $\tau_{\text{max}}$, we require OPE convergence of the $k=0$ LSRs as $s_{0}\to\infty$. By OPE convergence, we mean that the total perturbative contribution to the LSRs (pert), the total 4d contribution to the LSRs (4d), and the total 6d contribution to the LSRs (6d) must obey the inequality

[TABLE]

The lowest value of $\tau$ for which (46) is violated as $s_{0}\to\infty$ becomes $\tau_{\text{max}}$. Additionally, $\tau_{\text{max}}$ is constrained by the requirement

[TABLE]

where this inequality results from requiring that individually both $\mathcal{R}_{1}(\tau,s_{0})$ and $\mathcal{R}_{0}(\tau,s_{0})$ satisfy the Hölder inequalities [64, 65] as per [28]. For the specific LSRs being studied here, it turns out that the condition (46) is more restrictive than the condition (47). For both diquark channels under consideration, the values of $\tau_{\text{max}}$ obtained are given in the last column of Table 1. To select $\tau_{\text{min}}$, in addition to the Hölder inequality constraint (47), we require that

[TABLE]

i.e., that the resonance contribution to $\mathcal{R}_{1}/\mathcal{R}_{0}$ must be at least 50%. The highest value of $\tau$ which does not violate (47)–(48) becomes $\tau_{\text{min}}$. For both diquark channels under consideration, the values of $\tau_{\text{min}}$ obtained are given in the second-to-last column of Table 1.

The procedure described above for choosing a Borel window is $s_{0}$-dependent. However, $s_{0}$ is a parameter that is predicted using the optimization procedure described below. As such, choosing a Borel window and predicting $s_{0}$ are actually handled iteratively. Typically, the Borel window widens as $s_{0}$ increases. As such, we begin by selecting the minimum value of $s_{0}$ for which a Borel window exists. The corresponding Borel window is then used to predict a new, updated $s_{0}$. This new $s_{0}$ is then used to update the Borel window which, in turn, is used to update $s_{0}$ and so on until both the Borel window and $s_{0}$ settle. This iterative process has been taken into account in reporting diquark constituent masses, continuum thresholds, and Borel windows in Table 1.

To predict $s_{0}$ and $M$, we optimize the agreement between left- and right-hand sides of (45) by minimizing

[TABLE]

where we have partitioned the Borel window into 20 equal length subintervals with $\{\tau_{j}\}_{j=0}^{20}$. For both diquark channels under consideration, the optimized values of $s_{0}$ obtained are given in the third column of Table 1. As a consistency check on our methodology, we require that the optimized mass $M$ from (49) actually yields a good fit to (45) and that the left-hand side of (45) exhibits $\tau$ stability [28], that is

[TABLE]

within the Borel window. And so, in Figures 3 and 4, we plot the left-hand side of (45) at the appropriate optimized $s_{0}$ versus $\tau$ over the appropriate Borel window for both diquark channels under consideration. For the $bb$ diquark, the optimized $M$ from (49) does indeed yield a good fit to (45). Specifically, $M=8.67$ GeV in agreement with Figure 4. Regarding condition (50), over the Borel window,

[TABLE]

implying that the plot in Figure 4 can be considered flat to an excellent approximation. For the $cc$ diquarks, it is clear from Figure 3 that the fitted value of $M$ will be biased by the rapid increase at large $\tau$ values. We thus use the critical point $\frac{\mathrm{d}}{\mathrm{d}\tau}\sqrt{\mathcal{R}_{1}/\mathcal{R}_{0}}=0$ for our $cc$ diquark mass prediction, i.e., $M=3.51$ GeV. For both diquark channels under consideration, predicted diquark constituent masses $M$ are given in the second column of Table 1. The theoretical uncertainties associated with the mass predictions take into account the uncertainties arising from the strong coupling and mass parameters (29)–(32) as well as those associated with the QCD condensate values (33)–(35). The dominant theoretical uncertainty is associated with the quark masses.

In the $s_{0}\to\infty$ limit, the left-hand side of (45) corresponds to an upper bound on $M$ for a wide variety of resonance shapes [66], allowing the sensitivity to the threshold $s_{0}$ and resonance model to be explored. As shown in Figs. 5–6, within the Borel window $\tau<\tau_{\text{max}}$, we find $M\lesssim 3.6\ \text{GeV}$ for the $cc$ case and $M\lesssim 8.8\ \text{GeV}$ for the $bb$ case, remarkably close to the Table 1 predictions.

4 Discussion

Compared with potential model approaches [67, 68, 69] (and others cited therein) our $cc$ central value diquark constituent mass prediction is slightly larger and $bb$ is slightly smaller. For Bethe-Salpeter approaches [70], there is closer alignment in the $cc$ constituent mass prediction, but the $bb$ constituent mass prediction is still slightly smaller. However, taking into account theoretical uncertainties, we find good agreement between our QCD LSRs mass predictions and those of Refs. [67, 68, 69, 70], providing QCD evidence to support the study of diquark-antidiquark tetraquarks and doubly-heavy baryons with diquark cluster models.

Constituent diquark masses are key inputs into chromomagnetic interaction (CMI) models of diquark-antidiquark tetraquarks. For example, consider the Type-II model of Ref. [39] in which colour-spin interactions are ignored except between the quarks (antiquarks) within the diquark (antidiquark). This simplification assumes that the diquark and antidiquark within the tetraquark are point-like and well-separated. Focusing on S-wave combinations of doubly-heavy, equal mass diquarks and antidiquarks, the Type-II CMI Hamiltonian reduces to [39]

[TABLE]

where $m_{[Q_{1}Q_{1}]}$ and $m_{[\bar{Q}_{2}\bar{Q}_{2}]}$ are constituent diquark and antidiquark masses respectively and where $\kappa_{Q_{1}Q_{1}}$ and $\kappa_{\bar{Q}_{2}\bar{Q}_{2}}$ are colour-spin interaction coefficients. (Note that $\kappa_{\bar{Q}\bar{Q}}$ and $\kappa_{QQ}$ are equal as are $m_{[QQ]}$ and $m_{[\bar{Q}\bar{Q}]}$.) As the (anti-)diquarks have $J=1$, they must have $S=1$ for $L=0$ (where $J,\,L,\,S$ are the usual angular momentum quantum numbers). Hence, the Hamiltonian (52) simplifies to

[TABLE]

Our predictions for $m_{[cc]}$ and $m_{[bb]}$ are in Table 1; however, the coefficients $\kappa_{cc}$ and $\kappa_{bb}$ are not known. In [39], the $X(3870)$, $Z(3900)$, and $Z(4020)$ resonances were interpreted as Type-II diquark-antidiquark tetraquarks and were used to predict $\kappa_{cq}=67\ \text{MeV}$ where $q$ is a light quark. As the $\kappa$ coefficients are expected to decrease with increasing quark masses [23], we assume here that

[TABLE]

The absolute uncertainties in our diquark constituent mass predictions in Table 1 are significantly larger than 67 MeV, and so, as a first approximation, we simply ignore the $\kappa$ contributions to (53). Therefore, within the Type-II diquark-antidiquark model, we predict $J^{P}\in\{0^{+},\,1^{+},\,2^{+}\}$ tetraquark masses of 7.0 GeV for $[cc][\bar{c}\bar{c}]$, 12.2 GeV for $[cc][\bar{b}\bar{b}]$, and 17.3 GeV for $[bb][\bar{b}\bar{b}]$. The relative uncertainty in these mass predictions is roughly 10%. Furthermore, note that the $[cc][\bar{c}\bar{c}]$ and $[bb][\bar{b}\bar{b}]$ tetraquarks are charge conjugation eigenstates where $C=+$ for $J=0,\,2$ and $C=-$ for $J=1$ [71, 72]. The $[cc][\bar{b}\bar{b}]$ tetraquarks are not charge conjugation eigenstates.

Regarding $[cc][\bar{c}\bar{c}]$ tetraquarks, taking into account 10% theoretical uncertainty, our Type-II model mass predictions are in reasonable agreement with those of [71, 72, 73], although our central values are higher. However, our results are much higher than those of [74]. Furthermore, our tetraquark mass predictions are above both the $\eta_{c}(1S)$-$\eta_{c}(1S)$ and $J/\psi$-$J/\psi$ thresholds indicating that the corresponding decay modes should be accessible as fall-apart decays.

Regarding $[cc][\bar{b}\bar{b}]$ tetraquarks, again factoring in 10% uncertainty, our Type-II model mass predictions are in reasonable agreement with those of [71, 72], although our central values are lower. With an electric charge of $+2$, two charm quarks, and two bottom antiquarks, such a state would be easy to identify through its decay products, and could not be misinterpreted as a conventional meson. Unfortunately, within theoretical uncertainty, we are unable to say whether or not our tetraquark mass predictions lie above or below the $B_{c}^{+}$-$B_{c}^{+}$ threshold.

Regarding $[bb][\bar{b}\bar{b}]$ tetraquarks, taking into account theoretical uncertainty, our Type-II model mass predictions are in reasonable agreement with those of [73] although our central values are lower. Our results are about 10% lower than those of [74, 75], and are much lower than those of [71, 72]. Tetraquarks with $bb\bar{b}\bar{b}$ quark composition (so-called “beauty-full” tetraquarks) have attracted considerable attention recently due to the possibility that some might have masses below the $\Upsilon(1S)$-$\Upsilon(1S)$ threshold and perhaps even the $\eta_{b}(1S)$-$\eta_{b}(1S)$ threshold. For $bb\bar{b}\bar{b}$ tetraquarks with masses below the $\eta_{b}(1S)$-$\eta_{b}(1S)$ threshold, fall-apart modes would be inaccessible and decays would instead proceed through OZI-suppressed processes. Central values of our Type-II diquark-antidiquark mass estimates put the $0^{++}$, $1^{+-}$, and $2^{++}$ states about 9% below the $\Upsilon(1S)$-$\Upsilon(1S)$ threshold and about 7% below the $\eta_{b}(1S)$-$\eta_{b}(1S)$ threshold.

In summary, we have used QCD LSRs to predict the axial vector doubly-heavy $cc$ and $bb$ diquark constituent masses. Our results are summarized in Table 1. These results were obtained from a calculation of the diquark correlation function at NLO in perturbation theory and to LO in the 4d and 6d gluon condensates as well as the 6d quark condensate. That the LSRs analyses stabilized in both the double charm and double bottom diquark channels provides QCD-based evidence for the existence of these structures. Within the Type-II diquark-antidiquark tetraquark model of Ref. [39], we predict, with an uncertainty of roughly 10%, $0^{++}$, $1^{+-}$, and $2^{++}$ $[cc][\bar{c}\bar{c}]$ tetraquarks of mass 7.0 GeV; $0^{+}$, $1^{+}$, and $2^{+}$ $[cc][\bar{b}\bar{b}]$ tetraquarks of mass 12.2 GeV; and $0^{++}$, $1^{+-}$, and $2^{++}$ $[bb][\bar{b}\bar{b}]$ tetraquarks of mass 17.3 GeV. Central values of our $[bb][\bar{b}\bar{b}]$ tetraquark mass predictions are well below the $\Upsilon(1S)$-$\Upsilon(1S)$ and $\eta_{b}(1S)$-$\eta_{b}(1S)$ thresholds, providing support for the possibility that fall-apart decay modes are inaccessible to some $bb\bar{b}\bar{b}$ tetraquarks.

Acknowledgements

We are grateful for financial support from the National Sciences and Engineering Research Council of Canada (NSERC).