All-heavy tetraquarks

TL;DR

This paper calculates the mass spectra of fully-heavy tetraquarks using a potential model, finding that no bound states exist below meson pair thresholds, implying such states are unlikely to be observed experimentally.

Contribution

It provides the first detailed potential model analysis of all-heavy tetraquark systems, emphasizing the importance of the linear confining potential in mass predictions.

Findings

No bound states below meson pair thresholds.

Linear confining potential significantly affects mass spectra.

Fully-heavy tetraquarks unlikely to be observed as narrow states.

Abstract

In this work we study the mass spectra of the fully-heavy tetraquark systems, i.e. $cc\bar{c}\bar{c}$, $bb\bar{b}\bar{b}$, $bb\bar{c}\bar{c}/cc\bar{b}\bar{b}$, $bc\bar{c}\bar{c}/cc\bar{b}\bar{c}$, $bc\bar{b}\bar{b}/bb\bar{b}\bar{c}$, and $bc\bar{b}\bar{c}$, within a potential model by including the linear confining potential, Coulomb potential, and spin-spin interactions. It shows that the linear confining potential has important contributions to the masses and is crucial for our understanding of the mass spectra of the fully-heavy tetraquark systems. For the fully-heavy tetraquarks $Q_1Q_2\bar{Q}_3\bar{Q}_4$ our explicit calculations suggest that no bound states can be formed below the thresholds of any meson pairs $(Q_1\bar{Q}_3)$-$(Q_2\bar{Q}_4)$ or $(Q_1\bar{Q}_4)$-$(Q_2\bar{Q}_3)$. Thus, we do not expect narrow fully-heavy tetraquark states to be existing in experiments.

Tables

Table 1. Table 1: Configurations of all-heavy tetraquarks.

System	$J^{P\left(C\right)}$	Configuration
$cc\overline{c}\overline{c}$	$0^{++}$	${|{\left\{cc\right\}}_0^6{\left\{\overline{c}\overline{c}\right\}}_0^{\overline{6}}⟩}_0^0$	${|{\left\{cc\right\}}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_0^0$	$\cdots$
	$1^{+-}$	${|{\left\{cc\right\}}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_1^0$	$\cdots$	$\cdots$
	$2^{++}$	${|{\left\{cc\right\}}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_2^0$	$\cdots$	$\cdots$
$bb\overline{b}\overline{b}$	$0^{++}$	${|{\left\{bb\right\}}_0^6{\left\{\overline{b}\overline{b}\right\}}_0^{\overline{6}}⟩}_0^0$	${|{\left\{bb\right\}}_1^{\overline{3}}{\left\{\overline{b}\overline{b}\right\}}_1^3⟩}_0^0$	$\cdots$
	$1^{+-}$	${|{\left\{bb\right\}}_1^{\overline{3}}{\left\{\overline{b}\overline{b}\right\}}_1^3⟩}_1^0$	$\cdots$	$\cdots$
	$2^{++}$	${|{\left\{bb\right\}}_1^{\overline{3}}{\left\{\overline{b}\overline{b}\right\}}_1^3⟩}_2^0$	$\cdots$	$\cdots$
$bb\overline{c}\overline{c}$	$0^{+}$	${|{\left\{bb\right\}}_0^6{\left\{\overline{c}\overline{c}\right\}}_0^{\overline{6}}⟩}_0^0$	${|{\left\{bb\right\}}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_0^0$	$\cdots$
	$1^{+}$	${|{\left\{bb\right\}}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_1^0$	$\cdots$	$\cdots$
	$2^{+}$	${|{\left\{bb\right\}}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_2^0$	$\cdots$	$\cdots$
$bc\overline{c}\overline{c}$	$0^{+}$	${|{\left(bc\right)}_0^6{\left\{\overline{c}\overline{c}\right\}}_0^{\overline{6}}⟩}_0^0$	${|{\left(bc\right)}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_0^0$	$\cdots$
	$1^{+}$	${|{\left(bc\right)}_1^6{\left\{\overline{c}\overline{c}\right\}}_0^{\overline{6}}⟩}_1^0$	${|{\left(bc\right)}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_1^0$	${|{\left(bc\right)}_0^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_1^0$
	$2^{+}$	${|{\left(bc\right)}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_2^0$	$\cdots$	$\cdots$
$bc\overline{b}\overline{b}$	$0^{+}$	${|{\left(bc\right)}_0^6{\left\{\overline{b}\overline{b}\right\}}_0^{\overline{6}}⟩}_0^0$	${|{\left(bc\right)}_1^{\overline{3}}{\left\{\overline{b}\overline{b}\right\}}_1^3⟩}_0^0$	$\cdots$
	$1^{+}$	${|{\left(bc\right)}_1^6{\left\{\overline{b}\overline{b}\right\}}_0^{\overline{6}}⟩}_1^0$	${|{\left(bc\right)}_1^{\overline{3}}{\left\{\overline{b}\overline{b}\right\}}_1^3⟩}_1^0$	${|{\left(bc\right)}_0^{\overline{3}}{\left\{\overline{b}\overline{b}\right\}}_1^3⟩}_1^0$
	$2^{+}$	${|{\left(bc\right)}_1^{\overline{3}}{\left\{\overline{b}\overline{b}\right\}}_1^3⟩}_2^0$	$\cdots$	$\cdots$
$bc\overline{b}\overline{c}$	$0^{++}$	${|{\left(bc\right)}_1^6{\left(\overline{b}\overline{c}\right)}_1^{\overline{6}}⟩}_0^0$	${|{\left(bc\right)}_0^6{\left(\overline{b}\overline{c}\right)}_0^{\overline{6}}⟩}_0^0$	$\cdots$
		${|{\left(bc\right)}_1^{\overline{3}}{\left(\overline{b}\overline{c}\right)}_1^3⟩}_0^0$	${|{\left(bc\right)}_0^{\overline{3}}{\left(\overline{b}\overline{c}\right)}_0^3⟩}_0^0$	$\cdots$
	$1^{+-}$	${|{\left(bc\right)}_1^6{\left(\overline{b}\overline{c}\right)}_1^{\overline{6}}⟩}_1^0$	$\frac{1}{\sqrt{2}}{|{\left(bc\right)}_1^6{\left(\overline{b}\overline{c}\right)}_0^{\overline{6}}⟩}_1^0-{|{\left(bc\right)}_0^6{\left(\overline{b}\overline{c}\right)}_1^{\overline{6}}⟩}_1^0$	$\cdots$
		${|{\left(bc\right)}_1^{\overline{3}}{\left(\overline{b}\overline{c}\right)}_1^3⟩}_1^0$	$\frac{1}{\sqrt{2}}{|{\left(bc\right)}_1^{\overline{3}}{\left(\overline{b}\overline{c}\right)}_0^3⟩}_1^0-{|{\left(bc\right)}_0^{\overline{3}}{\left(\overline{b}\overline{c}\right)}_1^3⟩}_1^0$	$\cdots$
	$1^{++}$	$\frac{1}{\sqrt{2}}{|{\left(bc\right)}_1^6{\left(\overline{b}\overline{c}\right)}_0^{\overline{6}}⟩}_1^0+{|{\left(bc\right)}_0^6{\left(\overline{b}\overline{c}\right)}_1^{\overline{6}}⟩}_1^0$	$\frac{1}{\sqrt{2}}{|{\left(bc\right)}_1^{\overline{3}}{\left(\overline{b}\overline{c}\right)}_0^3⟩}_1^0+{|{\left(bc\right)}_0^{\overline{3}}{\left(\overline{b}\overline{c}\right)}_1^3⟩}_1^0$	$\cdots$
	$2^{++}$	${|{\left(bc\right)}_1^6{\left(\overline{b}\overline{c}\right)}_1^{\overline{6}}⟩}_2^0$	${|{\left(bc\right)}_1^{\overline{3}}{\left(\overline{b}\overline{c}\right)}_1^3⟩}_2^0$	$\cdots$

Table 2. Table 2: Quark model parameters used in this work.

$m_c$ (GeV)	1.483
$m_b$ (GeV)	4.852
${\alpha }_{cc}$	0.5461
${\alpha }_{bb}$	0.4311
${\alpha }_{bc}$	0.5021
${\sigma }_{cc}$ (GeV)	1.1384
${\sigma }_{bb}$ (GeV)	2.3200
${\sigma }_{bc}$ (GeV)	1.3000
$b$  (GeV 2)	0.1425

Table 3. Table 3: The masses (MeV) of bottomonium mesons. Experimental date are taken from PDG  Patrignani:2016xqp .

Meson	${\rm Y}$	${\eta }_b$	${\rm Y}\left(2S\right)$	${\eta }_b\left(2S\right)$	$h_b\left(1P\right)$	${\chi }_{b0}\left(1P\right)$	${\chi }_{b1}\left(1P\right)$	${\chi }_{b2}\left(1P\right)$
Ours	9460	9390	10024	10005	9941	9859	9933	9957
Expt.	9460	9399	10023	9999	9899	9859	9893	9912

Table 4. Table 4: Color matrix elements.

		$⟨{𝝀}_1\cdot {𝝀}_2⟩$	$⟨{𝝀}_3\cdot {𝝀}_4⟩$	$⟨{𝝀}_1\cdot {𝝀}_3⟩$	$⟨{𝝀}_2\cdot {𝝀}_4⟩$	$⟨{𝝀}_1\cdot {𝝀}_4⟩$	$⟨{𝝀}_2\cdot {𝝀}_3⟩$
	$⟨{\zeta }_1\left|\hat{O}\right|{\zeta }_1⟩$	4/3	4/3	$-10/3$	$-10/3$	$-10/3$	$-10/3$
	$⟨{\zeta }_2\left|\hat{O}\right|{\zeta }_2⟩$	$-8/3$	$-8/3$	$-4/3$	$-4/3$	$-4/3$	$-4/3$
	$⟨{\zeta }_1\left|\hat{O}\right|{\zeta }_2⟩$	0	0	$-2\sqrt{2}$	$-2\sqrt{2}$	$2\sqrt{2}$	$2\sqrt{2}$

Table 5. Table 5: Spin matrix elements.

		$⟨{𝝈}_1\cdot {𝝈}_2⟩$	$⟨{𝝈}_3\cdot {𝝈}_4⟩$	$⟨{𝝈}_1\cdot {𝝈}_3⟩$	$⟨{𝝈}_2\cdot {𝝈}_4⟩$	$⟨{𝝈}_1\cdot {𝝈}_4⟩$	$⟨{𝝈}_2\cdot {𝝈}_3⟩$
	$⟨{\chi }_0^{00}\left|\hat{O}\right|{\chi }_0^{00}⟩$	-3	-3	0	0	0	0
	$⟨{\chi }_0^{11}\left|\hat{O}\right|{\chi }_0^{11}⟩$	1	1	-2	-2	-2	-2
	$⟨{\chi }_0^{00}\left|\hat{O}\right|{\chi }_0^{11}⟩$	0	0	$-\sqrt{3}$	$-\sqrt{3}$	$\sqrt{3}$	$\sqrt{3}$
	$⟨{\chi }_1^{01}\left|\hat{O}\right|{\chi }_1^{01}⟩$	-3	1	0	0	0	0
	$⟨{\chi }_1^{10}\left|\hat{O}\right|{\chi }_1^{10}⟩$	1	-3	0	0	0	0
	$⟨{\chi }_1^{11}\left|\hat{O}\right|{\chi }_1^{11}⟩$	1	1	-1	-1	-1	-1
	$⟨{\chi }_1^{01}\left|\hat{O}\right|{\chi }_1^{10}⟩$	0	0	1	1	-1	-1
	$⟨{\chi }_1^{01}\left|\hat{O}\right|{\chi }_1^{11}⟩$	0	0	$-\sqrt{2}$	$\sqrt{2}$	$-\sqrt{2}$	$\sqrt{2}$
	$⟨{\chi }_1^{10}\left|\hat{O}\right|{\chi }_1^{11}⟩$	0	0	$\sqrt{2}$	$-\sqrt{2}$	$-\sqrt{2}$	$\sqrt{2}$
	$⟨{\chi }_2^{11}\left|\hat{O}\right|{\chi }_2^{11}⟩$	1	1	1	1	1	1

Table 6. Table 6: Predicted mass spectra for the c ​ c ​ c ¯ ​ c ¯ 𝑐 𝑐 ¯ 𝑐 ¯ 𝑐 cc\bar{c}\bar{c} , b ​ b ​ b ¯ ​ b ¯ 𝑏 𝑏 ¯ 𝑏 ¯ 𝑏 bb\bar{b}\bar{b} and b ​ b ​ c ¯ ​ c ¯ 𝑏 𝑏 ¯ 𝑐 ¯ 𝑐 bb\bar{c}\bar{c} systems.

$J^{P\left(C\right)}$	Configuration	$⟨H⟩$ (MeV)	Mass (MeV)	Eigenvector
$0^{++}$	${|{\left\{cc\right\}}_0^6{\left\{\overline{c}\overline{c}\right\}}_0^{\overline{6}}⟩}_0^0$		$\left[\begin{array}{c}6550\\ 6455\end{array}\right]$	$\left[\begin{array}{c}(0.81,-0.58)\\ (0.58,0.81)\end{array}\right]$
	${|{\left\{cc\right\}}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_0^0$
$1^{+-}$	${|{\left\{cc\right\}}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_1^0$	$\left(\begin{array}{c}6500\end{array}\right)$	6500	1
$2^{++}$	${|{\left\{cc\right\}}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_2^0$	$\left(\begin{array}{c}6524\end{array}\right)$	6524	1
$0^{++}$	${|{\left\{bb\right\}}_0^6{\left\{\overline{b}\overline{b}\right\}}_0^{\overline{6}}⟩}_0^0$		$\left[\begin{array}{c}19355\\ 19306\end{array}\right]$	$\left[\begin{array}{c}(0.81,-0.58)\\ (0.58,0.81)\end{array}\right]$
	${|{\left\{bb\right\}}_1^{\overline{3}}{\left\{\overline{b}\overline{b}\right\}}_1^3⟩}_0^0$
$1^{+-}$	${|{\left\{bb\right\}}_1^{\overline{3}}{\left\{\overline{b}\overline{b}\right\}}_1^3⟩}_1^0$	$\left(\begin{array}{c}19329\end{array}\right)$	19329	1
$2^{++}$	${|{\left\{bb\right\}}_1^{\overline{3}}{\left\{\overline{b}\overline{b}\right\}}_1^3⟩}_2^0$	$\left(\begin{array}{c}19341\end{array}\right)$	19341	1
$0^{+}$	${|{\left\{bb\right\}}_0^6{\left\{\overline{c}\overline{c}\right\}}_0^{\overline{6}}⟩}_0^0$		$\left[\begin{array}{c}13039\\ 12947\end{array}\right]$	$\left[\begin{array}{c}(0.97,-0.26)\\ (0.26,0.97)\end{array}\right]$
	${|{\left\{bb\right\}}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_0^0$
$1^{+}$	${|{\left\{bb\right\}}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_1^0$	$\left(\begin{array}{c}12960\end{array}\right)$	12960	1
$2^{+}$	${|{\left\{bb\right\}}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_2^0$	$\left(\begin{array}{c}12972\end{array}\right)$	12972	1

Table 7. Table 7: Our predicted masses (MeV) for the c ​ c ​ c ¯ ​ c ¯ 𝑐 𝑐 ¯ 𝑐 ¯ 𝑐 cc\bar{c}\bar{c} system compared with others.

State	Ours	Ref. Wu:2016vtq	Ref. Lloyd:2003yc	Ref. Chen:2016jxd	Ref. Ader:1981db	Ref. Iwasaki:1975pv	Ref. Karliner:2016zzc	Ref. Barnea:2006sd	Refs. Wang:2017jtz ; Wang:2018poa	Ref. Debastiani:2017msn	Ref. Berezhnoy:2011xn	Ref. Anwar:2017toa
$0^{++}$	6487	6797	6477	6460-6470	6437	6200	6192	6038-6115	5990	5969	5966
$0^{++}$	6518	7016	6695	6440-6820	6383	…	…	…	…	…	…
$1^{+-}$	6500	6899	6528	6370-6510	6437	…	…	6101-6176	6050	6021	6051
$2^{++}$	6524	6956	6573	6370-6510	6437	…	…	6172-6216	6090	6115	6223

Table 8. Table 8: Our predicted masses (MeV) for the b ​ b ​ b ¯ ​ b ¯ 𝑏 𝑏 ¯ 𝑏 ¯ 𝑏 bb\bar{b}\bar{b} system compared with others.

State	Ours	Ref. Wu:2016vtq	Ref. Wang:2017jtz ; Wang:2018poa	Ref. Karliner:2016zzc	Ref. Berezhnoy:2011xn	Ref. Anwar:2017toa	Ref. Bai:2016int	Ref. Chen:2016jxd	Ref. Hughes:2017xie	Ref. Anwar:2017toa
$0^{++}$	19322	20155	18840	18826	18754	18720	18690	18460-18490	18798
$0^{++}$	19338	20275	…	…	…	…	…	18450-19640	…
$1^{+-}$	19329	20212	18840	…	18808	…	…	18320-18540	…
$2^{++}$	19341	20243	18850	…	18916	…	…	18320-18530	…

Table 9. Table 9: The contributions from each part of the Hamiltonian of the c ​ c ​ c ¯ ​ c ¯ 𝑐 𝑐 ¯ 𝑐 ¯ 𝑐 cc\bar{c}\bar{c} and b ​ b ​ b ¯ ​ b ¯ 𝑏 𝑏 ¯ 𝑏 ¯ 𝑏 bb\bar{b}\bar{b} systems in units of MeV.

$J^{PC}$	Configuration	$M$	$⟨T⟩$	$⟨V^{Conf}⟩$	$⟨V_{coul}^{OGE}⟩$	$⟨V_{CM}^{OGE}⟩$
$0^{++}$	${|{\left\{cc\right\}}_0^6{\left\{\overline{c}\overline{c}\right\}}_0^{\overline{6}}⟩}_0^0$	6518	715	664	$-811$	18
	${|{\left\{cc\right\}}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_0^0$	6487	756	646	$-834$	$-13$
$1^{+-}$	${|{\left\{cc\right\}}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_1^0$	6500	739	653	$-825$	0
$2^{++}$	${|{\left\{cc\right\}}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_2^0$	6524	708	667	$-806$	23
$0^{++}$	${|{\left\{bb\right\}}_0^6{\left\{\overline{b}\overline{b}\right\}}_0^{\overline{6}}⟩}_0^0$	19338	768	356	$-1203$	9
	${|{\left\{bb\right\}}_1^{\overline{3}}{\left\{\overline{b}\overline{b}\right\}}_1^3⟩}_0^0$	19322	796	350	$-1225$	$-6$
$1^{+-}$	${|{\left\{bb\right\}}_1^{\overline{3}}{\left\{\overline{b}\overline{b}\right\}}_1^3⟩}_1^0$	19329	785	353	$-1216$	0
$2^{++}$	${|{\left\{bb\right\}}_1^{\overline{3}}{\left\{\overline{b}\overline{b}\right\}}_1^3⟩}_2^0$	19341	763	357	$-1199$	12

Table 10. Table 10: Predicted mass spectra for the b ​ c ​ c ¯ ​ c ¯ 𝑏 𝑐 ¯ 𝑐 ¯ 𝑐 bc\bar{c}\bar{c} and b ​ c ​ b ¯ ​ b ¯ 𝑏 𝑐 ¯ 𝑏 ¯ 𝑏 bc\bar{b}\bar{b} systems.

$J^P$	Configuration	$⟨H⟩$ (MeV)	Mass (MeV)	Eigenvector
$0^{+}$	${|{\left(bc\right)}_0^6{\left\{\overline{c}\overline{c}\right\}}_0^{\overline{6}}⟩}_0^0$		$\left[\begin{array}{c}9787\\ 9715\end{array}\right]$	$\left[\begin{array}{c}(0.81,-0.58)\\ (0.58,0.81)\end{array}\right]$
	${|{\left(bc\right)}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_0^0$
$1^{+}$	${|{\left(bc\right)}_1^6{\left\{\overline{c}\overline{c}\right\}}_0^{\overline{6}}⟩}_1^0$		$\left[\begin{array}{c}9773\\ 9752\\ 9727\end{array}\right]$	$\left[\begin{array}{c}(0.80,-0.21,0.56)\\ (-0.07,0.90,0.43)\\ (-0.59,-0.38,0.71)\end{array}\right]$
	${|{\left(bc\right)}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_1^0$
	${|{\left(bc\right)}_0^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_1^0$
$2^{+}$	${|{\left(bc\right)}_1^{\overline{3}}{\left\{\overline{c}\overline{c}\right\}}_1^3⟩}_2^0$	$\left(\begin{array}{c}9768\end{array}\right)$	9768	1
$0^{+}$	${|{\left(bc\right)}_0^6{\left\{\overline{b}\overline{b}\right\}}_0^{\overline{6}}⟩}_0^0$		$\left[\begin{array}{c}16190\\ 16141\end{array}\right]$	$\left[\begin{array}{c}(0.81,-0.58)\\ (0.58,0.81)\end{array}\right]$
	${|{\left(bc\right)}_1^{\overline{3}}{\left\{\overline{b}\overline{b}\right\}}_1^3⟩}_0^0$
$1^{+}$	${|{\left(bc\right)}_1^6{\left\{\overline{b}\overline{b}\right\}}_0^{\overline{6}}⟩}_1^0$		$\left[\begin{array}{c}16176\\ 16164\\ 16148\end{array}\right]$	$\left[\begin{array}{c}(0.82,0.04,0.57)\\ (0,1,-0.06)\\ (-0.57,0.05,0.82)\end{array}\right]$
	${|{\left(bc\right)}_1^{\overline{3}}{\left\{\overline{b}\overline{b}\right\}}_1^3⟩}_1^0$
	${|{\left(bc\right)}_0^{\overline{3}}{\left\{\overline{b}\overline{b}\right\}}_1^3⟩}_1^0$
$2^{+}$	${|{\left(bc\right)}_1^{\overline{3}}{\left\{\overline{b}\overline{b}\right\}}_1^3⟩}_2^0$	$\left(\begin{array}{c}16176\end{array}\right)$	16176	1

Table 11. Table 11: Predicted mass spectra for the b ​ c ​ b ¯ ​ c ¯ 𝑏 𝑐 ¯ 𝑏 ¯ 𝑐 bc\bar{b}\bar{c} system.

$J^{PC}$	Configuration	$⟨H⟩$ (MeV)	Mass (MeV)	Eigenvector
$0^{++}$	${|{\left(bc\right)}_1^6{\left(\overline{b}\overline{c}\right)}_1^{\overline{6}}⟩}_0^0$		$\left[\begin{array}{c}12854\\ 12931\\ 12975\\ 13024\end{array}\right]$	$\left[\begin{array}{c}(0.75,0.33,0.43,0.38)\\ (-0.49,0.68,-0.12,0.53)\\ (-0.27,0.37,0.69,-0.57)\\ (-0.34,-0.54,0.58,0.51)\end{array}\right]$
	${|{\left(bc\right)}_0^6{\left(\overline{b}\overline{c}\right)}_0^{\overline{6}}⟩}_0^0$
	${|{\left(bc\right)}_1^{\overline{3}}{\left(\overline{b}\overline{c}\right)}_1^3⟩}_0^0$
	${|{\left(bc\right)}_0^{\overline{3}}{\left(\overline{b}\overline{c}\right)}_0^3⟩}_0^0$
$1^{+-}$	${|{\left(bc\right)}_1^6{\left(\overline{b}\overline{c}\right)}_1^{\overline{6}}⟩}_1^0$		$\left[\begin{array}{c}12881\\ 12909\\ 13004\\ 13020\end{array}\right]$	$\left[\begin{array}{c}(0.62,0.54,0.37,0.44)\\ (0.61,-0.55,0.35,-0.46)\\ (-0.50,0.06,0.86,-0.09)\\ (-0.05,-0.63,0.10,0.77)\end{array}\right]$
	$\frac{1}{\sqrt{2}}{|{\left(bc\right)}_1^6{\left(\overline{b}\overline{c}\right)}_0^{\overline{6}}⟩}_1^0-{|{\left(bc\right)}_0^6{\left(\overline{b}\overline{c}\right)}_1^{\overline{6}}⟩}_1^0$
	${|{\left(bc\right)}_1^{\overline{3}}{\left(\overline{b}\overline{c}\right)}_1^3⟩}_1^0$
	$\frac{1}{\sqrt{2}}{|{\left(bc\right)}_1^{\overline{3}}{\left(\overline{b}\overline{c}\right)}_0^3⟩}_1^0-{|{\left(bc\right)}_0^{\overline{3}}{\left(\overline{b}\overline{c}\right)}_1^3⟩}_1^0$
$1^{++}$	$\frac{1}{\sqrt{2}}{|{\left(bc\right)}_1^6{\left(\overline{b}\overline{c}\right)}_0^{\overline{6}}⟩}_1^0+{|{\left(bc\right)}_0^6{\left(\overline{b}\overline{c}\right)}_1^{\overline{6}}⟩}_1^0$		$\left[\begin{array}{c}12933\\ 12992\end{array}\right]$	$\left[\begin{array}{c}(0.82,0.58)\\ (-0.58,0.82)\end{array}\right]$
	$\frac{1}{\sqrt{2}}{|{\left(bc\right)}_1^{\overline{3}}{\left(\overline{b}\overline{c}\right)}_0^3⟩}_1^0+{|{\left(bc\right)}_0^{\overline{3}}{\left(\overline{b}\overline{c}\right)}_1^3⟩}_1^0$
$2^{++}$	${|{\left(bc\right)}_1^6{\left(\overline{b}\overline{c}\right)}_1^{\overline{6}}⟩}_2^0$		$\left[\begin{array}{c}12933\\ 13021\end{array}\right]$	$\left[\begin{array}{c}(0.82,0.58)\\ (-0.58,0.82)\end{array}\right]$
	${|{\left(bc\right)}_1^{\overline{3}}{\left(\overline{b}\overline{c}\right)}_1^3⟩}_2^0$

Equations

$$|1⟩=|[Q_1{Q}_2{]}_1^6[\overline{Q_3}\overline{Q_4}{]}_1^{\bar{6}}{⟩}_0^0,|2⟩=|\{Q_1{Q}_2{\}}_0^6\{\overline{Q_3}\overline{Q_4}{\}}_0^{\bar{6}}{⟩}_0^0,$$

$$|3⟩=|\{Q_1{Q}_2{\}}_1^{\bar{3}}\{\overline{Q_3}\overline{Q_4}{\}}_1^3{⟩}_0^0,|4⟩=|[Q_1{Q}_2{]}_0^{\bar{3}}[\overline{Q_3}\overline{Q_4}{]}_0^3{⟩}_0^0,$$

$$|5⟩=|[Q_1{Q}_2{]}_1^6[\overline{Q_3}\overline{Q_4}{]}_1^{\bar{6}}{⟩}_1^0,|6⟩=|[Q_1{Q}_2{]}_1^6\{\overline{Q_3}\overline{Q_4}{\}}_0^{\bar{6}}{⟩}_1^0,$$

$$|7⟩=|\{Q_1{Q}_2{\}}_0^6[\overline{Q_3}\overline{Q_4}{]}_1^{\bar{6}}{⟩}_1^0,|8⟩=|\{Q_1{Q}_2{\}}_1^{\bar{3}}\{\overline{Q_3}\overline{Q_4}{\}}_1^3{⟩}_1^0,$$

$$|9⟩=|\{Q_1{Q}_2{\}}_1^{\bar{3}}[\overline{Q_3}\overline{Q_4}{]}_0^3{⟩}_1^0,|10⟩=|[Q_1{Q}_2{]}_0^{\bar{3}}\{\overline{Q_3}\overline{Q_4}{\}}_1^3{⟩}_1^0,$$

$$|11⟩=|[Q_1{Q}_2{]}_1^6[\overline{Q_3}\overline{Q_4}{]}_1^{\bar{6}}{⟩}_2^0,|12⟩=|\{Q_1{Q}_2{\}}_1^{\bar{3}}\{\overline{Q_3}\overline{Q_4}{\}}_1^3{⟩}_2^0,$$

$$|6^{\prime }⟩=\frac{1}{\sqrt{2}}(|(bc{)}_1^6(\bar{b}\bar{c}{)}_0^{\bar{6}}{⟩}_1^0-|(bc{)}_0^6(\bar{b}\bar{c}{)}_1^{\bar{6}}{⟩}_1^0),$$

$$|7^{\prime }⟩=\frac{1}{\sqrt{2}}(|(bc{)}_1^6(\bar{b}\bar{c}{)}_0^{\bar{6}}{⟩}_1^0+|(bc{)}_0^6(\bar{b}\bar{c}{)}_1^{\bar{6}}{⟩}_1^0),$$

$$|9^{\prime }⟩=\frac{1}{\sqrt{2}}(|(bc{)}_1^{\bar{3}}(\bar{b}\bar{c}{)}_0^3{⟩}_1^0-|(bc{)}_0^{\bar{3}}(\bar{b}\bar{c}{)}_1^3{⟩}_1^0),$$

$$|10^{\prime }⟩=\frac{1}{\sqrt{2}}(|(bc{)}_1^{\bar{3}}(\bar{b}\bar{c}{)}_0^3{⟩}_1^0+|(bc{)}_0^{\bar{3}}(\bar{b}\bar{c}{)}_1^3{⟩}_1^0),$$

$$H=(\sum _{i=1}^4{m}_i+T_i)-T_G+\sum _{i<j}{V}_{ij}(r_{ij}),$$

$$V_{ij}(r_{ij})=V_{ij}^{OGE}(r_{ij})+V_{ij}^{Conf}(r_{ij}),$$

$$V^{OGE}_{ij}=\frac{\alpha_{ij}}{4}({\mbox{\boldmath$\lambda$\unboldmath}}_{i}\cdot{\mbox{\boldmath$\lambda$\unboldmath}}_{j})\left\{\frac{1}{r_{ij}}-\frac{\pi}{2}\cdot\frac{\sigma^{3}_{ij}e^{-\sigma^{2}_{ij}r_{ij}^{2}}}{\pi^{3/2}}\cdot\frac{4}{3m_{i}m_{j}}({\mbox{\boldmath$\sigma$\unboldmath}}_{i}\cdot{\mbox{\boldmath$\sigma$\unboldmath}}_{j})\right\},$$

$$V^{Conf}_{ij}(r_{ij})=-\frac{3}{16}({\mbox{\boldmath$\lambda$\unboldmath}}_{i}\cdot{\mbox{\boldmath$\lambda$\unboldmath}}_{j})\cdot br_{ij},$$

$${\zeta }_1=|6\bar{6}⟩=|(Q_1{Q}_2{)}^6({\bar{Q}}_3{\bar{Q}}_4{)}^{\bar{6}}{⟩}^0,$$

$${\zeta }_2=|\bar{3}3⟩=|(Q_1{Q}_2{)}^{\bar{3}}({\bar{Q}}_3{\bar{Q}}_4{)}^3{⟩}^0.$$

$$\begin{split}\zeta_{1}=&\frac{1}{2\sqrt{6}}\bigg{[}(rb+br)(\bar{b\mathstrut}\bar{r\mathstrut}+\bar{r\mathstrut}\bar{b\mathstrut})+(gr+rg)(\bar{g\mathstrut}\bar{r\mathstrut}+\bar{r\mathstrut}\bar{g\mathstrut})\\ &+(gb+bg)(\bar{b\mathstrut}\bar{g\mathstrut}+\bar{g\mathstrut}\bar{b\mathstrut})\\ &+2(rr)(\bar{r}\bar{r})+2(gg)(\bar{g}\bar{g})+2(bb)(\bar{b}\bar{b})\bigg{]},\end{split}$$

$$\begin{split}\zeta_{2}=&\frac{1}{2\sqrt{3}}\bigg{[}(br-rb)(\bar{b\mathstrut}\bar{r\mathstrut}-\bar{r\mathstrut}\bar{b\mathstrut})-(rg-gr)(\bar{g\mathstrut}\bar{r\mathstrut}-\bar{r\mathstrut}\bar{g\mathstrut})\\ &+(bg-gb)(\bar{b\mathstrut}\bar{g\mathstrut}-\bar{g\mathstrut}\bar{b\mathstrut})\bigg{]},\end{split}$$

$${\chi }_0^{00}=|(Q_1{Q}_2{)}_0({\bar{Q}}_3{\bar{Q}}_4{)}_0{⟩}_0,$$

$${\chi }_0^{11}=|(Q_1{Q}_2{)}_1({\bar{Q}}_3{\bar{Q}}_4{)}_1{⟩}_0,$$

$${\chi }_1^{01}=|(Q_1{Q}_2{)}_0({\bar{Q}}_3{\bar{Q}}_4{)}_1{⟩}_1,$$

$${\chi }_1^{10}=|(Q_1{Q}_2{)}_1({\bar{Q}}_3{\bar{Q}}_4{)}_0{⟩}_1,$$

$${\chi }_1^{11}=|(Q_1{Q}_2{)}_1({\bar{Q}}_3{\bar{Q}}_4{)}_1{⟩}_1,$$

$${\chi }_2^{11}=|(Q_1{Q}_2{)}_1({\bar{Q}}_3{\bar{Q}}_4{)}_1{⟩}_2.$$

$${\chi }_0^{00}$$

$${\chi }_0^{11}$$

$${\chi }_1^{01}$$

$${\chi }_1^{10}$$

$${\chi }_1^{11}$$

$${\chi }_2^{11}$$

$$\psi ({\mathbf{r}}_1,{\mathbf{r}}_2,{\mathbf{r}}_3,{\mathbf{r}}_4)=\prod _{i=1}^4\sum _{\ell =1}^n{\mathcal{C}}_{i\ell }{\left(\frac{1}{\pi b_{i\ell }^2}\right)}^{3/4}\exp \left[-\frac{1}{2b_{i\ell }^2}{r}_i^2\right],$$

$$\psi ({\mathbf{r}}_1,{\mathbf{r}}_2,{\mathbf{r}}_3,{\mathbf{r}}_4)$$

$$⟨\psi ({\omega }_{\ell },{\mathbf{r}}_i,{\mathbf{r}}_j)\left|\frac{1}{r_{ij}}\right|\psi ({\omega }_{{\ell }^{\prime }},{\mathbf{r}}_i,{\mathbf{r}}_j)⟩=2\sqrt{\frac{m_{ij}}{\pi }}\frac{({\omega }_{\ell }{\omega }_{{\ell }^{\prime }}{)}^{3/2}}{(\frac{{\omega }_{\ell }+{\omega }_{{\ell }^{\prime }}}{2}{)}^{5/2}},$$

$$⟨\psi ({\omega }_{\ell },{\mathbf{r}}_i,{\mathbf{r}}_j)\left|e^{-{\sigma }_{ij}^2{r}_{ij}^2}\right|\psi ({\omega }_{{\ell }^{\prime }},{\mathbf{r}}_i,{\mathbf{r}}_j)⟩={\left(\frac{m_{ij}(\frac{2{\omega }_{\ell }{\omega }_{{\ell }^{\prime }}}{{\omega }_{\ell }+{\omega }_{{\ell }^{\prime }}})}{m_{ij}\frac{{\omega }_{\ell }+{\omega }_{{\ell }^{\prime }}}{2}+{\sigma }_{ij}^2}\right)}^{\frac{3}{2}},$$

$$⟨\psi ({\omega }_{\ell },{\mathbf{r}}_i,{\mathbf{r}}_j)\left|r_{ij}\right|\psi ({\omega }_{{\ell }^{\prime }},{\mathbf{r}}_i,{\mathbf{r}}_j)⟩=2\sqrt{\frac{1}{\pi m_{ij}}}\frac{({\omega }_{\ell }{\omega }_{{\ell }^{\prime }}{)}^{3/2}}{(\frac{{\omega }_{\ell }+{\omega }_{{\ell }^{\prime }}}{2}{)}^{7/2}},$$

$$\displaystyle\mbox{\boldmath$\xi$\unboldmath}_{1}$$

$$\displaystyle\mbox{\boldmath$\xi$\unboldmath}_{2}$$

$$\displaystyle\mbox{\boldmath$\xi$\unboldmath}_{3}$$

$$\displaystyle\mbox{\boldmath$\xi$\unboldmath}_{4}$$

$$\psi(\mbox{\boldmath$\xi$\unboldmath}_{1},\mbox{\boldmath$\xi$\unboldmath}_{2},\mbox{\boldmath$\xi$\unboldmath}_{3},\mbox{\boldmath$\xi$\unboldmath}_{4})=\sum^{n}_{\ell}\mathcal{C}_{\ell}\prod_{i=1}^{4}\left(\frac{\mu_{i}\omega_{\ell}}{\pi}\right)^{3/4}\exp\left[-\frac{\mu_{i}\omega_{\ell}}{2}\xi^{2}_{i}\right],$$

$$⟨\sum _{i=1}^4{T}_i-T_G⟩=\frac{9}{4}\sum _{\ell }^n\sum _{{\ell }^{\prime }}^n{\mathcal{C}}_{\ell }{\mathcal{C}}_{{\ell }^{\prime }}\frac{({\omega }_{\ell }{\omega }_{{\ell }^{\prime }}{)}^4}{(\frac{{\omega }_{\ell }+{\omega }_{{\ell }^{\prime }}}{2}{)}^7}.$$

Full text

I Introduction

Experimental searches for and theoretical studies of exotic hadrons beyond the conventional quark model are an important test of nonperturbative properties of the strong interaction theory QCD. Since the discovery of quark model GellMann:1964nj and QCD, the progresses on the experimental tools have brought to us a lot of novel phenomena in hadron physics. In particular, during the past 15 years, there have been a sizeable number of candidates for QCD exotics Patrignani:2016xqp ; Olsen:2017bmm ; Lebed:2016hpi ; Chen:2016qju ; Ali:2017jda ; Esposito:2016noz ; Guo:2017jvc . Interestingly, but also puzzlingly, it shows that the number of exotic candidates is far less than what we have expected for the hadron spectroscopy, where the internal effective degrees of freedom of a hadron may contain quarks and gluons beyond the conventional quark model prescription. Strong evidences for such exotic hadrons include some of those recently observed $XYZ$ states, e.g, $X(3872)$, $Z_{c}(3900)$, $Z_{c}(4020)$, $Z_{b}(10610)$, and $Z_{b}(10650)$ Patrignani:2016xqp . In particular, these charged quarkoniumlike states, $Z_{c}$ and $Z_{b}$, contain not only the hidden heavy flavor $c\bar{c}$ or $b\bar{b}$, but also charged light flavors of $u\bar{d}$ or $d\bar{u}$. Since at least four constituent quarks are confined inside these $Z_{c}$ or $Z_{b}$ states, it makes them the best candidates for QCD exotic hadrons.

Recently, the tetraquarks of all-heavy systems, such as $cc\bar{c}\bar{c}$ and $bb\bar{b}\bar{b}$, have received considerable attention with the development of experiments. If there are stable tetraquark $cc\bar{c}\bar{c}$ and/or $bb\bar{b}\bar{b}$ states, they are most likely to be observed at LHC Eichten:2017ual . In fact, a search for the tetraquark $bb\bar{b}\bar{b}$ states is being carried out by the LHCb Collaboration although no confirmed information has been observed Aaij:2018zrb . Other study interests for physicists arise from the special aspects of the all-heavy tetraquark systems Chen:2016jxd . They may favor to form genuine tetraquark configurations rather than loosely bound hadronic molecules, since the light mesons cannot be exchanged between two heavy mesons. Furthermore, it will be very easy to distinguish the all-heavy tetraquark states from the states which have been observed because their masses should be far away from the mass regions of the observed states. Thus, besides some previous works on the all-heavy tetraquark states Ader:1981db ; Iwasaki:1975pv ; Zouzou:1986qh ; Heller:1985cb ; Lloyd:2003yc ; Barnea:2006sd , many new studies have been carried out in recent years Wang:2017jtz ; Karliner:2016zzc ; Berezhnoy:2011xn ; Bai:2016int ; Anwar:2017toa ; Esposito:2018cwh ; Chen:2016jxd ; Wu:2016vtq ; Hughes:2017xie ; Richard:2018yrm ; Debastiani:2017msn ; Wang:2018poa ; Richard:2017vry ; Vijande:2009kj , although some of the conclusions are quite different from each other. In some works, it is predicted that there exist stable bound tetraquark $cc\bar{c}\bar{c}$ states and/or bound tetraquark $bb\bar{b}\bar{b}$ states with relatively smaller masses below the thresholds of heavy charmonium pairs Wang:2017jtz ; Karliner:2016zzc ; Berezhnoy:2011xn ; Bai:2016int ; Anwar:2017toa ; Esposito:2018cwh ; Chen:2016jxd ; Debastiani:2017msn ; Wang:2018poa . Thus, their decays into heavy quarkonium pairs through quark rearrangements will be hindered. In contrast, in some other works it is predicted that there should be no stable bound tetraquark $cc\bar{c}\bar{c}$ and $bb\bar{b}\bar{b}$ states Wu:2016vtq ; Lloyd:2003yc ; Ader:1981db ; Hughes:2017xie ; Richard:2018yrm because the predicted masses are large enough for them to decay into heavy quarkonium pairs. To some extent, a better understanding of the possible mass locations is not only crucial for understanding their underlying dynamics, but also useful for experimental searches for their existence.

In this work, we systemically study the mass spectra of the all-heavy tetraquark $Q_{1}Q_{2}\bar{Q}_{3}\bar{Q}_{4}$ systems with a potential model widely used in the literature Godfrey:2004ya ; Eichten:1994gt ; Godfrey:1985xj ; Swanson:2005 ; Godfrey:2015dia ; Eichten:1978tg ; Gupta:1984jb ; Kwong:1988ae ; Chao:2009 ; Li:2009zu ; Segovia:2016xqb ; Wei-Zhao:2013sta ; Lakhina:2006fy ; Lu:2016bbk ; Li:2010vx ; Deng:2016stx ; Deng:2016ktl ; Song:2015nia ; Song:2015fha . Our purpose is to understand two key issues based on the knowledge collected in the study of heavy quarkonium spectrum. The first one is what a quark potential model can tell about the all-heavy tetraquark system. The second one is where the masses of the ground states could be located if the all-heavy tetraquark states do exist.

At this moment, we do not consider any orbital or radial excitations of the all-heavy tetraquarks. Instead, we would like to address where and how the all-heavy tetraquarks would manifest themselves in their lowest states. For a spectrum of multiquark states, a correct identification of the ground state should be the first step towards a better understanding of the multiquark dynamics in the nonperturbative regime.

The potentials between the quarks, such as the linear confining potential, color Coulomb potential, and spin-spin interactions, are adopted the standard forms of the potential models. The model parameters are determined by fitting the mass spectra of charmonium, bottomonium, and $B_{c}$ meson. In our calculations, we find both the confining potential and color Coulomb potential are very crucial for understanding the masses of the all-heavy tetraquarks. The linear confining potential as well as the kinetic energy contributes a quite large positive mass term to the all-heavy tetraquarks $Q_{1}Q_{2}\bar{Q}_{3}\bar{Q}_{4}$, which leads to a large mass far above the threshold of the meson pair $Q_{1}\bar{Q}_{3}$-$Q_{2}\bar{Q}_{4}$ or $Q_{1}\bar{Q}_{4}$-$Q_{2}\bar{Q}_{3}$, although the color Coulomb potential contributes a very large negative mass term. As a consequence, we find no bound all-heavy tetraquarks $Q_{1}Q_{2}\bar{Q}_{3}\bar{Q}_{4}$ below the threshold of any meson pairs $Q_{1}\bar{Q}_{3}$-$Q_{2}\bar{Q}_{4}$ or $Q_{1}\bar{Q}_{4}$-$Q_{2}\bar{Q}_{3}$.

The paper is organized as follows: a brief introduction to the framework is given in Sec. II. In Sec. III, the numerical results and discussions are presented. A short summary is given in Sec. IV.

II Framework

II.1 Quark model classification

In the charm and bottom quark sector, there are nine different fully-heavy $Q_{1}Q_{2}\bar{Q}_{3}\bar{Q}_{4}$ systems: $cc\bar{c}\bar{c}$, $bb\bar{b}\bar{b}$, $bc\bar{c}\bar{c}$, $bc\bar{b}\bar{b}$, $bb\bar{c}\bar{c}$, $bc\bar{b}\bar{c}$, $cc\bar{b}\bar{c}$, $bb\bar{b}\bar{c}$ and $cc\bar{b}\bar{b}$. Note that $cc\bar{b}\bar{c}$, $bb\bar{b}\bar{c}$, and $cc\bar{b}\bar{b}$ are the antiparticles of $bc\bar{c}\bar{c}$, $bc\bar{b}\bar{b}$, and $bb\bar{c}\bar{c}$, respectively. Thus, we need only consider six systems, $cc\bar{c}\bar{c}$, $bb\bar{b}\bar{b}$, $bc\bar{c}\bar{c}$, $bc\bar{b}\bar{b}$, $bb\bar{c}\bar{c}$, and $bc\bar{b}\bar{c}$, in our calculations.

To calculate the spectroscopy of a $Q_{1}Q_{2}\bar{Q}_{3}\bar{Q}_{4}$ system, first we construct the configurations in the space of flavor $\otimes$color$\otimes$spin. Considering the Pauli principle and color confinement for the four-quark system $Q_{1}Q_{2}\bar{Q}_{3}\bar{Q}_{4}$, we have 12 configurations as follows:

[TABLE]

where $\{~{}\}$ and $[~{}]$ denote the symmetric and antisymmetric flavor wave functions of the two quarks (antiquarks) subsystems, respectively. The subscripts and superscripts are the spin quantum numbers and representations of the color SU(3) group, respectively. A symmetric spatial wave function is implied for the ground states under investigation.

It should be emphasized that for the $bc\bar{b}\bar{c}$ systems the $J=1$ states can have both $C=\pm 1$, which can be constructed by the linear combinations of $|6\rangle$, $|7\rangle$, $|9\rangle$ and $|10\rangle$,

[TABLE]

where configurations $|6^{\prime}\rangle$ and $|9^{\prime}\rangle$ have $C=-1$, and $|7^{\prime}\rangle$ and $|10^{\prime}\rangle$ have $C=+1$. Since the permutation symmetries are lost for the bc and $\bar{b}\bar{c}$ subsystems, in this work, we use $(~{})$ denote no permutation symmetries for these quark pair subsystems.

In Table 1, all possible configurations and corresponding quantum numbers for the $cc\bar{c}\bar{c}$, $bb\bar{b}\bar{b}$, $bb\bar{c}\bar{c}$, $bc\bar{c}\bar{c}$, $bc\bar{b}\bar{b}$ and $bc\bar{b}\bar{c}$ systems are listed.

II.2 Hamiltonian for the multiquark system

The following nonrelativistic Hamiltonian is adopted for the calculation of the masses of the all-heavy $Q_{1}Q_{2}\bar{Q}_{3}\bar{Q}_{4}$ system:

[TABLE]

where $m_{i}$ and $T_{i}$ stand for the constituent quark mass and kinetic energy of the $i$th quark, respectively; $T_{G}$ stands for the center-of-mass (c.m.) kinetic energy of the $Q_{1}Q_{2}\bar{Q}_{3}\bar{Q}_{4}$ system; $r_{ij}\equiv|\mathbf{r}_{i}-\mathbf{r}_{j}|$ is the distance between the $i$th quark and $j$th quark; and $V_{ij}(r_{ij})$ stands for the effective potential between the $i$th and $j$th quark. In this work, we adopt a widely used potential form for $V_{ij}(r_{ij})$ Eichten:1978tg ; Godfrey:1985xj ; Swanson:2005 ; Godfrey:2015dia ; Godfrey:2004ya ; Lakhina:2006fy ; Lu:2016bbk ; Li:2010vx ; Deng:2016stx ; Deng:2016ktl , i.e,

[TABLE]

where $V^{OGE}_{ij}$ stands for the one-gluon-exchange (OGE) potential which describes the short-range quark-quark interactions, while $V^{Conf}_{ij}(r_{ij})$ stands for the confinement potential which describes the long-range interaction behaviors. The form of $V^{OGE}_{ij}$ is given by

[TABLE]

where \mbox{\boldmath\sigma\unboldmath}_{i} are the Pauli matrices, and $\alpha_{ij}$ stands for the strong coupling strength between two quarks. If the interaction occurs between two quarks or antiquarks, the \mbox{\boldmath\lambda\unboldmath}_{i}\cdot\mbox{\boldmath\lambda\unboldmath}_{j} operator appearing in Eq. (7) is defined as \mbox{\boldmath\lambda\unboldmath}_{i}\cdot\mbox{\boldmath\lambda\unboldmath}_{j}\equiv\sum_{a=1}^{8}\lambda_{i}^{a}\lambda_{j}^{a}, while if the interaction occurs between a quark and antiquark, the \mbox{\boldmath\lambda\unboldmath}_{i}\cdot\mbox{\boldmath\lambda\unboldmath}_{j} operator is defined as \mbox{\boldmath\lambda\unboldmath}_{i}\cdot\mbox{\boldmath\lambda\unboldmath}_{j}\equiv\sum_{a=1}^{8}-\lambda_{i}^{a}\lambda_{j}^{a*}, where $\lambda^{a*}$ is the complex conjugate of the Gell-Mann matrix $\lambda^{a}$. The OGE potential $V^{OGE}_{ij}$ is composed of the Coulomb type potential V^{OGE}_{coul}\propto(\mbox{\boldmath\lambda\unboldmath}_{i}\cdot\mbox{\boldmath\lambda\unboldmath}_{j})(1/r_{ij}) which provides the short-range interaction, and the color-magnetic interaction V^{OGE}_{CM}\propto(\mbox{\boldmath\lambda\unboldmath}_{i}\cdot\mbox{\boldmath\lambda\unboldmath}_{j})({\mbox{\boldmath\sigma\unboldmath}}_{i}\cdot{\mbox{\boldmath\sigma\unboldmath}}_{j}) which provides mass splittings. The form of $V^{Conf}_{ij}(r_{ij})$ is given by

[TABLE]

where the parameter $b$ denotes the strength of the confinement potential.

There are nine parameters $m_{c}$, $m_{b}$, $\alpha_{cc}$, $\alpha_{bb}$, $\alpha_{bc}$, $\sigma_{cc}$, $\sigma_{bb}$, $\sigma_{bc}$, and $b$ to be determined in the calculations. In Refs. Deng:2016stx ; Li:2019tbn , the masses of $c\bar{c}$ and $b\bar{c}$ spectrum are calculated by using the three-point difference central method cai:2003ktl from the center ($r=0$) towards outside ($r\to\infty$) point by point. The parameters $m_{c}$, $\alpha_{cc}$, $\sigma_{cc}$, $b$, $m_{b}$, $\alpha_{bc}$, $\sigma_{bc}$ have been determined. In this work, we use the same method to determine the parameters $\alpha_{bb}$, $\sigma_{bb}$, by fitting the masses of $b\bar{b}$ spectrum. The parameter set is listed in Table 2. The corresponding theoretical results for the masses of heavy quarkonia $b\bar{b}$ are shown in Table 3.

II.3 Matrix elements in color and spin spaces

In order to obtain the mass of a tetraquark state from the nonrelativistic Hamiltonian defined in Eq. (5), first one needs to calculate the matrix elements of \langle\mbox{\boldmath\lambda\unboldmath}_{i}\cdot\mbox{\boldmath\lambda\unboldmath}_{j}\rangle and \langle\mbox{\boldmath\sigma\unboldmath}_{i}\cdot\mbox{\boldmath\sigma\unboldmath}_{j}\rangle in the color and spin spaces, respectively.

In the color space, one has two kinds of a color-singlet state,

[TABLE]

According to the SU(3) Clebsch-Gordan coefficients, one easily obtains the expressions of the color wave functions as follows Liu:2004kd ; Liu:2016ogz ; deSwart:1963pdg ; Kaeding:1995vq :

[TABLE]

[TABLE]

with these color wave functions one can work out the matrix elements \langle\mbox{\boldmath\lambda\unboldmath}_{i}\cdot\mbox{\boldmath\lambda\unboldmath}_{j}\rangle Vijande:2009ac , which have been summarized in Table 4.

In the spin space, one has six spin wave functions,

[TABLE]

According to the SU(2) Clebsch-Gordan coefficients, we easily obtain the expressions of the spin wave functions as follows:

[TABLE]

with these spin wave functions one can work out the matrix elements of \langle\mbox{\boldmath\sigma\unboldmath}_{i}\cdot\mbox{\boldmath\sigma\unboldmath}_{j}\rangle Vijande:2009ac , which have been listed in Table 5.

II.4 Matrix elements in the coordinate space

The trail wave function of the four-quark states without any spatial excitations in the coordinate space is expanded by a series of Gaussian functions,

[TABLE]

where the parameters $b_{i\ell}$ are related to the harmonic oscillator frequencies $\omega_{\ell}$ with $1/b^{2}_{i\ell}=m_{i}\omega_{\ell}$. With the ansatz of the mass independent frequency $\omega_{\ell}$ for a quark of mass $m_{i}$ Straub:1988gj , i.e., $1/b^{2}_{i\ell}=1/b^{2}_{\ell}=m_{u}\omega_{\ell}$ ($m_{u}$=313 MeV), the trail wave function of the four quark states can be simplified to be

[TABLE]

which is often adopted in the calculations of the multiquark systems Zhang:2007mu ; Zhang:2005jz .

In the coordinate space, we need work out the matrix elements of $\langle 1/r_{ij}\rangle$, $\langle e^{-\sigma^{2}_{ij}r_{ij}^{2}}\rangle$, and $\langle r_{ij}\rangle$. Combing the trail wave functions defined in Eq. (26), we obtain

[TABLE]

[TABLE]

[TABLE]

where $\psi(\omega_{i},\mathbf{r}_{i},\mathbf{r}_{j})\equiv\phi(\omega_{\ell},\mathbf{r}_{i})\phi(\omega_{\ell},\mathbf{r}_{j})$, $m_{ij}=m_{i}m_{j}/(m_{i}+m_{j})$.

To separate out the center-of-mass kinetic energy $T_{G}$ and finally work out the kinetic energy matrix element $\langle\sum_{i=1}^{4}T_{i}-T_{G}\rangle$ , we need to redefine the coordinates by the following Jacobi coordinates,

[TABLE]

with these one can rewritten the Eq. (26) as

[TABLE]

where $\mu_{1}\equiv m_{1}m_{2}/(m_{1}+m_{2})$, $\mu_{2}\equiv m_{3}m_{4}/(m_{3}+m_{4})$, $\mu_{3}\equiv(m_{1}+m_{2})(m_{3}+m_{4})/M$, $\mu_{4}=M\equiv m_{1}+m_{2}+m_{3}+m_{4}$. With the trail wave function defined in Eq. (34), the kinetic energy matrix element is worked out to be

[TABLE]

III results and discussions

In this work, we adopt the variation principle to solve the Schrödinger equation. Following the method used in Ref. Hiyama:2003cu , the oscillator length $b_{\ell}$ are set to be

[TABLE]

where $n$ is the number of Gaussian functions, and $a$ is the ratio coefficient. There are three parameters $\{b_{1},b_{n},n\}$ to be determined through variation method. It is found that when we take $b_{1}=0.1$ fm, $b_{n}=4$ fm, $n=15$, we will obtain stable solutions for the four-quark systems.

When all the matrix elements have been worked out, we can solve the generalized matrix eigenvalue problem,

[TABLE]

where

[TABLE]

[TABLE]

with $\psi(\omega_{\ell})=\prod_{i=1}^{4}\left(\frac{\mu_{i}\omega_{\ell}}{\pi}\right)^{3/4}\exp\left[-\frac{\mu_{i}\omega_{\ell}}{2}\xi^{2}_{i}\right]$. $\varphi$, $\zeta$, and $\chi$ stand for the flavor, color, and spin wave functions, respectively. The physical state corresponds to the solution with a minimum energy $E_{m}$. By solving this generalized matrix eigenvalue problem, the mass of the tetraquark configuration and its spacial wave function can be determined.

III.1 The $cc\bar{c}\bar{c}$ and $bb\bar{b}\bar{b}$ systems

The predicted mass spectrum for the $cc\bar{c}\bar{c}$ system has been given in Table 6 and also shown in Fig. 1 (a). From Table 6, it is found that in the two $I^{G}(J^{PC})=0^{+}(0^{++})$ states there is a sizable configuration mixing between $|\{cc\}^{6}_{0}\{\bar{c}\bar{c}\}^{\bar{6}}_{0}\rangle^{0}_{0}$ and $|\{cc\}^{\bar{3}}_{1}\{\bar{c}\bar{c}\}^{3}_{1}\rangle^{0}_{0}$. The $J^{PC}=0^{++}$ state with a higher mass, 6550 MeV, is dominant by the $|\{cc\}^{6}_{0}\{\bar{c}\bar{c}\}^{\bar{6}}_{0}\rangle^{0}_{0}$ configuration, while another $J^{PC}=0^{++}$ state with a lower mass, 6455 MeV, is dominant by the $|\{cc\}^{\bar{3}}_{1}\{\bar{c}\bar{c}\}^{3}_{1}\rangle^{0}_{0}$ configuration. The mass splitting between these two $J^{PC}=0^{++}$ states is about $95$ MeV. The other two states $J^{PC}=1^{+-}$ and $J^{PC}=2^{++}$ are also located in a similar mass region, i.e. $\sim 6.5$ GeV, and the mass splitting between them is about 20 MeV. As shown in Fig. 1 (a), the two $J^{PC}=0^{++}$ states are above the thresholds of the charmonium pairs for about $260\sim 580$ MeV. It suggests that the $J^{PC}=0^{++}$ states are unstable, and they can easily decay into the $\eta_{c}\eta_{c}$ and $J/\psi J/\psi$ final states through quark rearrangements. The $J^{PC}=1^{+-}$ state lies about 420 MeV above the mass threshold of $\eta_{c}J/\psi$, while $J^{PC}=2^{++}$ is about 330 MeV above the mass threshold of $J/\psi J/\psi$, they can also easily decay into $\eta_{c}J/\psi$ and $J/\psi J/\psi$; respectively, through the quark rearrangements.

As a comparison, our predicted masses and some other typical results from other works are collected in Table 7. It shows that our predicted masses for the $cc\bar{c}\bar{c}$ system are roughly compatible with the nonrelativistic quark model predictions of Refs. Lloyd:2003yc ; Ader:1981db , where both confining and Coulomb potentials are considered. It is also interesting to find that similar results are given by the QCD sum rules Chen:2016jxd . In contrast, the masses predicted by us are much larger than those predicted in Refs. Berezhnoy:2011xn ; Debastiani:2017msn ; Wang:2017jtz ; Barnea:2006sd ; Karliner:2016zzc ; Anwar:2017toa . These methods which obtained small masses have some common features: either no confining potentials were explicitly included  Berezhnoy:2011xn ; Barnea:2006sd ; Karliner:2016zzc ; Anwar:2017toa or a diquark picture was adopted in the calculations Wang:2017jtz ; Debastiani:2017msn . Recently, Wu et al. also obtained a large mass $\sim 6.8-7.0$ GeV for the $cc\bar{c}\bar{c}$ system with the heavier constituent $c$-quark mass $1.72$ GeV adopted Wu:2016vtq .

We further analyze the contributions from each part of the Hamiltonian for the $cc\bar{c}\bar{c}$ system. The results are listed in Table 9. It shows that the averaged kinetic energy $\langle T\rangle$, the confining potential $\langle V^{Conf}\rangle$, and the Coulomb potential $\langle V^{OGE}_{coul}\rangle$ have the same order of magnitude. In particular, the contributions from the confining potential are sizeable and apparently cannot be neglected. Note that the confining potential contributes a positive energy to the system. Thus, neglecting this contribution will lead to much lower masses for the all-heavy system. In Refs. Berezhnoy:2011xn ; Barnea:2006sd ; Karliner:2016zzc ; Anwar:2017toa , the confining potential was explicitly neglected. Although part of the confining potential effects can be taken into account by the effective constituent quark masses in the ground states, our calculation shows that the impact from the inclusion of the confining potential seems not to be on the constituent quark masses in the heavy quark sector, but rather on the relative strengths of the averaged matrix elements among the terms of the nonrelativistic Hamiltonian.

In order to examine the role played by the confining potential in the spectrum of heavy quark system, we compare the contributions from the OGE and confining potential for the $\eta_{c}$ meson, i.e, $\langle V^{OGE}_{coul}\rangle\simeq-637$ MeV and $\langle V^{Conf}\rangle\simeq 233$ MeV, which are consistent with our previous study in Ref. Deng:2016stx . The ratio between the confining potential $\langle V^{Conf}\rangle$ and color Coulomb potential $\langle V^{OGE}_{coul}\rangle$ can reach up to

[TABLE]

This explicit result suggests that the neglect of confining potential cannot be justified for the $c\bar{c}$ system.

As a general conclusion, we find that the confining potential has significant contributions to the masses of the $cc\bar{c}\bar{c}$ system, and are the same order of magnitude as the color Coulomb potential. This will enhance the masses of the $cc\bar{c}\bar{c}$ system and does not support the existence of a bound tetraquark of $cc\bar{c}\bar{c}$ with narrow widths.

The predicted mass spectrum for the $bb\bar{b}\bar{b}$ system is very similar to that for the $cc\bar{c}\bar{c}$ one. The results are given in Table 6 and shown in Fig. 1 (b). The configuration mixing effects between $|\{bb\}^{6}_{0}\{\bar{b}\bar{b}\}^{\bar{6}}_{0}\rangle^{0}_{0}$ and $|\{bb\}^{\bar{3}}_{1}\{\bar{b}\bar{b}\}^{3}_{1}\rangle^{0}_{0}$ should be obvious in the two $J^{PC}=0^{++}$ states, the higher mass state with mass 19355 MeV is dominant by the $|\{bb\}^{6}_{0}\{\bar{b}\bar{b}\}^{\bar{6}}_{0}\rangle^{0}_{0}$ configuration, and the lower one of 19306 MeV is dominant by the $|\{bb\}^{\bar{3}}_{1}\{\bar{b}\bar{b}\}^{3}_{1}\rangle^{0}_{0}$ configuration. Due to the heavier mass of the $b$ quark, relatively smaller mass splittings among these states are found. The pattern is also similar to that of the $cc\bar{c}\bar{c}$ system. Note that the predicted masses are above the thresholds of the bottomonium pairs for about $380\sim 560$ MeV. It suggests that bound states of the $bb\bar{b}\bar{b}$ system with narrow widths are not favored.

In Table 8 we compare our results with other model calculations. It shows that our predicted masses are higher than most of the other predictions which are either calculated without including the confining potential explicitly  Berezhnoy:2011xn ; Barnea:2006sd ; Karliner:2016zzc ; Anwar:2017toa , or based on the diquark picture Wang:2017jtz . Similarly, based on the diquark picture, the lightest mass of $bb\bar{b}\bar{b}$ is estimated at $18.8$ GeV by Ref. Esposito:2018cwh . In these calculations the tetraquark states of $J^{PC}=0^{++}$, $1^{+-}$, or $2^{++}$ are either below or slightly above the thresholds of $\eta_{b}\eta_{b}$, $\eta_{b}\Upsilon(1S)$ or $\Upsilon(1S)\Upsilon(1S)$, respectively. Thus, they can become stable with narrow decay widths. In contrast, our calculations with the inclusion of the confining potential result in higher masses for the $bb\bar{b}\bar{b}$ system and do not favor the existence of such narrow tetraquark states. We note that a rather large mass $\sim 20.2$ GeV for the $bb\bar{b}\bar{b}$ system is estimated by Ref. Wu:2016vtq , where a heavier constituent $b$-quark mass $5.05$ GeV is adopted.

In Table 9, the contributions from each part of the Hamiltonian for the $bb\bar{b}\bar{b}$ system are listed. It shows that the kinetic energy $\langle T\rangle\simeq 800$ MeV, the confining potential $\langle V^{Conf}_{ij}(r_{ij})\rangle\simeq 400$ MeV, and the coulomb potential $\langle V^{OGE}_{coul}\rangle\simeq-1200$ MeV, have the same order of magnitude. As shown in Fig. 1 (b), the mass splittings among these $J^{PC}=0^{++}$, $1^{+-}$, and $2^{++}$ states follow a similar pattern as in the $cc\bar{c}\bar{c}$ system. Also similar to that for the $cc\bar{c}\bar{c}$ system, the neglect of the confining potential will lead to much lower masses for the $bb\bar{b}\bar{b}$ system, and this may explain the low masses obtained in Refs. Berezhnoy:2011xn ; Barnea:2006sd ; Karliner:2016zzc ; Anwar:2017toa . Although it is often argued that the confining potential contributions are perturbative for the bottomonium system, explicit calculations seem not to support this phenomenon. In Ref. Deng:2016stx , we have studied the $b\bar{b}$ spectrum and find that $\langle V^{Conf}(r)\rangle\simeq 122$ MeV and $\langle V^{OGE}_{coul}\rangle\simeq-970$ MeV for the $\eta_{b}$ meson. The ratio between $\langle V^{Conf}\rangle$ and $\langle V^{OGE}_{coul}\rangle$ can reach up to

[TABLE]

For the four heavy quark system of $bb\bar{b}\bar{b}$, the increase of the displacements between the two quarks (antiquarks) or quark-antiquark will experience larger confining forces. Thus, the confining potential contributions cannot be neglected in the calculations. As a consequence, our study does not support the existence of the tetraquark $bb\bar{b}\bar{b}$ bound states with narrow widths.

Finally, it should be mentioned that for a simplicity, in our calculation, the variational wave functions of the coordinate space are only adopted an $s$-wave form. Thus, the color wave functions for the $J^{PC}=1^{+-}$ and $2^{++}$ states is color $\bar{3}3$. However, the $\bar{3}3$ color wave functions for the $J^{PC}=1^{+-}$ and $2^{++}$ states might slightly mix with the color $6\bar{6}$ when one considers the orbital excitations in the coordinate space Richard:2017vry ; Vijande:2009kj ; Richard:2018yrm . With a color mixing effect, the mass of the $J^{PC}=1^{+-}$ and $2^{++}$ states might become slightly lower Richard:2017vry ; Vijande:2009kj ; Richard:2018yrm , which does not affect our conclusions.

III.2 The $bb\bar{c}\bar{c}$ system

The $bb\bar{c}\bar{c}$ system is similar to the $cc\bar{c}\bar{c}$ and $bb\bar{b}\bar{b}$ ones except that it does not have determined $C$ parity, and there is no contributions from the annihilation potential. The predicted mass spectrum for the $bb\bar{c}\bar{c}$ system is also listed in Table 6 and shown in Fig. 1 (c). From Table 6, a small configuration mixing between $|\{bb\}^{6}_{0}\{\bar{c}\bar{c}\}^{\bar{6}}_{0}\rangle^{0}_{0}$ and $|\{bb\}^{\bar{3}}_{1}\{\bar{c}\bar{c}\}^{3}_{1}\rangle^{0}_{0}$ can be identified. The higher mass state (13039 MeV) of $J^{P}=0^{+}$ state is dominant by the $|\{cc\}^{6}_{0}\{\bar{c}\bar{c}\}^{\bar{6}}_{0}\rangle^{0}_{0}$ configuration, while the lower mass one (12947 MeV) of the same quantum numbers is dominant by the $|\{cc\}^{\bar{3}}_{1}\{\bar{c}\bar{c}\}^{3}_{1}\rangle^{0}_{0}$ configuration. The mass splitting between these two $J^{P}=0^{+}$ states is about $92$ MeV. The other two states with $J^{P}=1^{+}$ and $2^{+}$ have a small mass splitting of about 10 MeV and are located around $12.96$ GeV. The masses predicted by us are about 600 MeV systematically smaller than those predicted in the recent work Wu:2016vtq , where relatively large constituent quark masses for the $b$ quark $5.05$ GeV and $c$ quark $1.72$ GeV are adopted.

As shown in Fig. 1 (c), all these states are above their lowest open flavor decay channels for about 300 MeV. Therefore, they can decay into the $B_{c}B_{c}$, $B_{c}^{*}B_{c}^{*}$, or $B_{c}B_{c}^{*}$ final states via the quark rearrangement quite easily.

III.3 The $bc\bar{c}\bar{c}$ and $bc\bar{b}\bar{b}$ systems

The states of both $bc\bar{c}\bar{c}$ and $bc\bar{b}\bar{b}$ systems do not have determined $C$ parity and they share some common features in terms of heavy quark symmetry. The predicted mass spectra for these two configurations are listed in Table 10 and shown in Fig. 1 (d) and Fig. 1 (e), respectively. It shows that both $bc\bar{c}\bar{c}$ and $bc\bar{b}\bar{b}$ systems have sizeable configuration mixings between the color $6\otimes\bar{6}$ and $3\otimes\bar{3}$ configurations. For the $bc\bar{c}\bar{c}$ system the mixing occurs between the $|\{bc\}^{6}_{0}\{\bar{c}\bar{c}\}^{\bar{6}}_{0}\rangle^{0}_{0}$ and $|\{bc\}^{\bar{3}}_{1}\{\bar{c}\bar{c}\}^{3}_{1}\rangle^{0}_{0}$ configurations. It shows that the higher and lower mass states of $J^{P}=0^{+}$ are dominated by the $|\{bc\}^{6}_{0}\{\bar{c}\bar{c}\}^{\bar{6}}_{0}\rangle^{0}_{0}$ and $|\{bc\}^{\bar{3}}_{1}\{\bar{c}\bar{c}\}^{3}_{1}\rangle^{0}_{0}$ configuration, respectively.

The configuration mixing effects among these three $J^{P}=1^{+}$ states are also sizeable which are shown in Table 10. The typical mass splitting is about 20 MeV and the predicted masses are about 500 MeV systematically smaller than those predicted in the recent work Wu:2016vtq . Again, we note that rather large constituent quark masses for the $c$ and $b$ quarks are adopted in Ref. Wu:2016vtq .

As a consequence of the high masses predicted by our model, namely, the states of $bc\bar{c}\bar{c}$ system are about $290-350$ MeV above the mass threshold of $B_{c}^{*}J/\psi$, we find that these states can easily decay into the $B_{c}\eta_{c}$, $B_{c}J/\psi$ or $B_{c}^{*}J/\psi$ final states via the quark rearrangements. Thus, we do not expect narrow states of $bc\bar{c}\bar{c}$ to be observed in experiment.

For the $bc\bar{b}\bar{b}$ system its main properties is very similar to that of the $bc\bar{c}\bar{c}$ system as shown in Table 10 and Fig. 1 (e). Instead of repeating the features seen in the $bc\bar{c}\bar{c}$ system, we only note the main features arising from the heavy constituent quark masses. Namely, the mass splittings among the multiplets with the same quantum numbers are expected to be smaller than that for the $bc\bar{c}\bar{c}$ system. For instance, the mass splitting among the $J^{P}=1^{+}$ states is about 10 MeV.

As shown in Fig. 1 (e), our results show that the states of the $bc\bar{b}\bar{b}$ system are about $350-390$ MeV above the mass threshold of $B_{c}^{*}\Upsilon$. Thus, these states with different quantum numbers can also easily decay into the $B_{c}\eta_{b}$, $B_{c}\Upsilon$ or $B_{c}^{*}\Upsilon$ final states via the quark rearrangement. Narrow states made of the $bc\bar{b}\bar{b}$ are not favored in our model.

III.4 The $bc\bar{b}\bar{c}$ system

The $bc\bar{b}\bar{c}$ system has no constraints from the Pauli principle, and there are 12 different configurations allowed by this system, namely, four $J^{PC}=0^{++}$ states, four $J^{PC}=1^{+-}$ states, two $J^{PC}=1^{++}$ states, and two $J^{PC}=2^{++}$ states. The predicted mass spectrum is listed in Table 11 and shown in Fig. 1 (f).

A main feature of the $bc\bar{b}\bar{c}$ system is that the configuration mixing appears to play an important role. For example, the highest mass $J^{PC}=0^{++}$ state is a mixed state containing comparable components from three configurations $|(bc)^{6}_{0}(\bar{b}\bar{c})^{\bar{6}}_{0}\rangle^{0}_{0}$, $|(bc)^{\bar{3}}_{1}(\bar{b}\bar{c})^{3}_{1}\rangle^{0}_{0}$ and $|(bc)^{\bar{3}}_{0}(\bar{b}\bar{c})^{3}_{0}\rangle^{0}_{0}$. As a consequence, the predicted masses for these tetraquark states are in the range of $12939\pm 85$ MeV. We note that our predicted masses are about 600 MeV systematically smaller than those predicted in Ref. Wu:2016vtq and about $350-600$ MeV systematically larger than those predicted with diquark picture in Ref. Berezhnoy:2011xn . Also, these states of $bc\bar{b}\bar{c}$ are about $200-300$ MeV above the mass threshold of $B_{c}^{*}B_{c}^{*}$. It suggests that these tetraquark states may easily decay into the $B_{c}B_{c}^{*}$, $B_{c}^{*}B_{c}^{*}$, $\eta_{b}J/\psi$, $\eta_{b}\eta_{c}$, $\Upsilon\eta_{c}$, or $\Upsilon J/\psi$ channels via quark rearrangements. Thus, they are expected to be broad in width.

IV Summary

In this work, we study the mass spectra of the all-heavy $cc\bar{c}\bar{c}$, $bb\bar{b}\bar{b}$, $bb\bar{c}\bar{c}/cc\bar{b}\bar{b}$, $bc\bar{c}\bar{c}/cc\bar{b}\bar{c}$, $bc\bar{b}\bar{b}/bb\bar{b}\bar{c}$, and $bc\bar{b}\bar{c}$ systems in the potential quark model with the linear confining potential, Coulomb potential, and spin-spin interactions included. We find that the linear confining potential contributes large positive energies to the eigenvalues of the ground states of these tetraquark systems. This is different from some existing calculations in the literature in which the neglect of the confining potential contributions leads to relatively low masses for the all-heavy systems and some of those can be lower than the two-body decay thresholds. In our case, all these states are found to have masses above the corresponding two meson decay thresholds via the quark rearrangement. This implies that narrow all-heavy tetraquark states may not exist in reality. Nevertheless, our explicit calculations suggest that the confining potential still plays an important role in the heavy flavor multiquark system, and it is crucial to include it in dynamical calculations in order to gain a better understanding of the multiquark dynamics. More experimental information from the Belle-II and LHCb analyses would be able to clarify these issues in the near future.

Acknowledgement

This work is supported by the National Natural Science Foundation of China (Grants No. 11775078, No. U1832173, No. 11705056, No. 11425525 No. 11521505). Q.Z. is also supported in part, by the DFG and NSFC funds to the Sino-German CRC 110 “Symmetries and the Emergence of Structure in QCD” (NSFC Grant No. 11261130311), National Key Basic Research Program of China under Contract No. 2015CB856700.