Doubly heavy tetraquark bound and resonant states
Wei-Lin Wu Yao Ma Yan-Ke Chen Lu Meng Shi-Lin Zhu

	 
	wlwu@pku.edu.cn
	School of Physics, Peking University, Beijing 100871, China
	yaoma@pku.edu.cn
	School of Physics and Center of High Energy Physics,
		Peking University, Beijing 100871, China
	chenyanke@stu.pku.edu.cn
	School of Physics, Peking University, Beijing 100871, China
	lu.meng@rub.de
	Institut für Theoretische Physik II, Ruhr-Universität Bochum,  D-44780 Bochum, Germany 
	zhusl@pku.edu.cn
	School of Physics and Center of High Energy Physics,
		Peking University, Beijing 100871, China
	
	§ ABSTRACT 
		We calculate the energy spectrum of the S-wave doubly heavy tetraquark systems, including the QQ^(')q̅q̅,  QQ^(')s̅q̅, and QQ^(')s̅s̅ (Q^(')=b,c and q=u,d) systems within the constituent quark model. We use the complex scaling method to obtain bound states and resonant states simultaneously, and the Gaussian expansion method to solve the complex-scaled four-body Schrödinger equation. With a novel definition of the root-mean-square radii, we are able to distinguish between meson molecules and compact tetraquark states. The compact tetraquarks are further classified into three different types with distinct spatial configurations: compact even tetraquarks, compact diquark-antidiquark tetraquarks and compact diquark-centered tetraquarks. In the I(J^P)=0(1^+) QQq̅q̅ system, there exists the D^*D molecular bound state with a binding energy of -14 MeV, which is the candidate for T_cc(3875)^+. The shallow B̅^*B̅ molecular bound state is the bottom analogue of T_cc(3875)^+. Moreover, we identify two resonant states near the D^*D^* and B̅^*B̅^* thresholds. In the J^P=1^+ bbq̅q̅ (I=0) and bbs̅q̅ systems, we obtain deeply bound states with a compact diquark-centered tetraquark configuration and a dominant χ_3̅_c⊗ 3_c component, along with resonant states with similar configurations as their radial excitations. These states are the QCD analogue of the helium atom.  We also obtain some other bound states and resonant states with “QCD Hydrogen molecule" configurations. Moreover, we investigate the heavy quark mass dependence of the I(J^P)=0(1^+) QQq̅q̅ bound states. We strongly urge the experimental search for the predicted states.
		
		
	
	
	Doubly heavy tetraquark bound and resonant states
    Shi-Lin Zhu 0000-0002-4055-6906
    
=================================================



§ INTRODUCTION
 
	
	In 2021, the LHCb Collaboration discovered the first doubly charmed tetraquark state T_cc(3875)^+ in the D^0D^0π^+ invariant mass spectrum <cit.>. It is a narrow state with a mass extremely close to the  D^*+D^0 threshold, having a binding energy of only around -300 keV. The observation of T_cc(3875)^+ significantly advances the hadron spectroscopy and may open a new chapter for the discovery of other doubly heavy exotic states in the future. 
	
	Theoretical investigations on the possible existence of doubly heavy tetraquark bound states dated back to the 1980s <cit.>. Many studies aimed at predicting the masses of doubly charmed tetraquark states using various frameworks <cit.>, but their conclusions were wildly inconsistent, with the predicted masses ranging from -300 MeV to +300 MeV relative to the DD^* threshold. After the discovery of T_cc(3875)^+, its exotic properties have attracted much attention and reignited interest in doubly heavy tetraquark states <cit.>. 
	The proximity of T_cc(3875)^+ to the D^*+D^0 threshold favors its interpretation as a D^*D molecular state, while for other doubly heavy tetraquark systems, both the compact tetraquark picture and the hadronic molecular picture have been proposed. More discussions can be found in recent reviews <cit.>.
	
	The discovery of the doubly charmed tetraquark bound state implies the existence of other doubly heavy tetraquark states. While many studies focus on the existence and properties of bound states, some also explore possible resonant states. In Ref. <cit.>, the authors employed the complex scaling method to study the ccq̅q̅, bbq̅q̅, bcq̅q̅ (q=u,d) bound and resonant states in the chiral quark model. However, they predicted a deeply bound ccq̅q̅ state with a binding energy of around -150 MeV, which contradicts with the experimental results. In Ref. <cit.>, the author used the heavy quark spin symmetry to predict a resonant pole T_cc' below the D^*D^* channel as a partner of T_cc(3875). Similar results were obtained in the constituent quark model <cit.> and lattice QCD <cit.>.  In Ref. <cit.>, the authors adopted a constituent quark model including the one-pion exchange interaction to study the ccq̅q̅ and bbq̅q̅ systems using the real scaling method. They did not find any resonant states in the doubly charmed sector, but reported a bbq̅q̅ resonant state. 
	
	In this work, we conduct a comprehensive study on the S-wave doubly heavy tetraquark systems, including the QQ^(')q̅q̅,  QQ^(')s̅q̅, and QQ^(')s̅s̅ (Q^(')=b,c and q=u,d)  systems, within the constituent quark model. We utilize the complex scaling method <cit.> to obtain possible bound states and resonant states simultaneously. We employ the Gaussian expansion method <cit.> to solve the four-body Schrödinger equation, which has been successfully used in our previous work on tetraquark bound states <cit.> and resonant states <cit.>. Moreover, we calculate the root-mean-square radii of the tetraquark states to analyze their spatial structures and distinguish between meson molecular states and compact tetraquark states. 
    We further classify the compact tetraquark states into three different types with distinct spatial configurations, unraveling the rich internal structures and different forming mechanisms of the tetraquark states.
	
	This paper is organized as follows. In Sec. <ref>, we introduce the theoretical framework, including the constituent quark model, the complex scaling method and the wave function construction. In Sec. <ref>, we demonstrate different spatial structures of tetraquarks and how to distinguish them by calculating the root-mean-square radii. In Sec. <ref>, we present the numerical results and discuss the properties of doubly heavy tetraquark states. We summarize our findings in Sec. <ref>.
	
	
 
	

§ THEORETICAL FRAMEWORK

	

 §.§ Hamiltonian

	In a nonrelativistic quark potential model, the Hamiltonian of a tetraquark system in the center-of-mass frame reads
	
    H=∑_i=1^4 (m_i+p_i^2/2 m_i)+∑_i<j=1^4 V_ij
    		,

	where the last term represents the two-body interaction between the i-th and j-th (anti)quark. We adopt the AL1 potential <cit.>, which includes the one-gluon-exchange interaction and a linear quark confinement interaction,
	
    V_i j =-3/16λ_i ·λ_j(-κ/r_i j+λ r_i j-Λ. 
    .+8 πκ^'/3 m_i m_jexp(-r_i j^2 / r_0^2)/π^3 / 2 r_0^3S_i ·S_j),

	where  λ_i and S_i are the SU(3) color Gell-Mann matrix and the spin operator acting on quark i, respectively. The parameters of the AL1 model were determined by fitting the meson spectra across all flavor sectors  <cit.>. No additional free parameters are introduced. The theoretical masses as well as the root-mean-square (rms) radii of the corresponding mesons are listed in Table <ref>. It can be seen that the theoretical results for the 1S mesons agree with the experimental values within tens of MeV. We also list possible experimental candidates for the 2S excited mesons, whose theoretical masses deviate from the experimental values by up to 100 MeV. However, these candidates are not yet well established as 2S mesons. Moreover, relativistic effects and coupled-channel effects may play a crucial role in understanding these excited mesons, which is beyond the scope of this work. Nonetheless, since we focus only on the tetraquark states below the M(1S)M'(2S) dimeson thresholds in this work, we expect the uncertainties to be within tens of MeV, similar to those of the 1S mesons.
	


	

 §.§ Complex scaling method

	In contrast to bound states, the wave functions of resonant states are not square-integrable and cannot be obtained by solving the eigenequation of the Hermitian Hamiltonian directly. The complex scaling method (CSM) enables us to obtain possible bound states and resonant states simultaneously. In the CSM <cit.>, the coordinate r and its conjugate momentum p are transformed as
	
    U(θ) r=r e^i θ,    U(θ) p=p e^-i θ.

	Under such a transformation, the Hamiltonian is analytically continued to the complex plane and no longer Hermitian,
	
    H(θ)=∑_i=1^4 (m_i+p_i^2e^-2iθ/2 m_i)+∑_i<j=1^4 V_ij(r_ije^iθ).

	By solving the eigenvalue equation of the complex-scaled Hamiltonian, we can obtain the eigenenergies of bound states, scattering states and resonant states simultaneously. Bound states are located on the negative real axis in the energy plane and remain unchanged as θ varies. Scattering states align along rays starting from threshold energies with Arg(E)=-2θ.  If the complex scaling angle θ is chosen to be larger than the angle of the resonant state θ_r=1/2tan^-1(Γ_r/2M_r), where M_r and Γ_r represent its mass and width, the wave function of the resonant state becomes square integrable and can be obtained with eigenenergy E_r=M_r-iΓ_r/2. 
	
	
	

 §.§ Wave function

	The basis functions of the tetraquark wave function are expressed as
	
    ψ=𝒜(ϕ⊗χ),

	where 𝒜 is the antisymmetric operator of identical particles, ϕ and χ represent the spatial wave function and color-spin wave function, respectively. 
	
	For the spatial wave function, we employ the Gaussian expansion method (GEM) <cit.>. We consider three sets of spatial configurations (dimeson and diquark-antidiquark), which are denoted by (jac) = (a), (b), (c). In each configuration, there are three independent Jacobian coordinates r_jac, λ_jac, ρ_jac, as shown in Fig. <ref>. The S-wave spatial basis function is written as
	
    ϕ^( jac)_n_1,n_2,n_3=ϕ_n_1(r_jac)ϕ_n_2(λ_jac)ϕ_n_3(ρ_jac),

	where ϕ_n_i(r) takes the Gaussian form,
	
    ϕ_n_i(r)=N_n_ie^-ν_n_ir^2,
       ν_n_i=ν_1γ^n_i-1   (n_i=1∼ n_ max),

	N_n_i is the normalization factor. 
	
	For the color-spin wave function, we choose a complete set of bases written as
	
    χ^s_1,s_2,S_3̅_c⊗ 3_c=[(Q_1Q_2)_3̅_c^s_1(q̅_1q̅_2)_3_c^s_2]_1_c^S,
    χ^s_1,s_2,S_6_c⊗6̅_c=[(Q_1Q_2)_6_c^s_1(q̅_1q̅_2)_6̅_c^s_2]_1_c^S,

	where the subscripts and superscripts denote the color and spin representations, respectively.
	
	

 
	

§ SPATIAL STRUCTURES

	Tetraquark states are generally classified into meson molecules and compact tetraquarks. In the molecular scheme, the (anti)quarks cluster into two color-singlet mesons, and their relative distance is expected to be larger than the typical range of color confinement, Λ_QCD^-1∼ 1 fm, which is also the typical size of a meson. In contrast, in the compact tetraquark scheme, all four (anti)quarks are confined together, with their relative distances being on the order of Λ_QCD^-1. The compact tetraquark scheme can be further subdivided into several types, particularly for doubly heavy systems, as illustrated in Fig. <ref>. In compact even tetraquarks, the relative distance between four (anti)quarks is of similar size. In compact diquark-antidiquark tetraquarks, two quarks and two antiquarks form two clusters, respectively. The sizes of the clusters are smaller than the distance between them. In compact diquark-centered tetraquarks, the two heavy quarks form a very compact diquark cluster,  while the two light antiquarks orbit around the diquark, similar to the way two electrons orbit around the nucleus in a helium atom. The compact diquark-centered tetraquarks correspond to the “QCD Helium atom" in Ref. <cit.>. The compact even tetraquarks  roughly correspond the “QCD Hydrogen molecule" states as coined in Ref. <cit.>. The two light antiquarks are shared by two heavy quarks as in the Hydrogen molecule where two electrons are shared by two protons.  It should be noted that in this paper (anti)diquark only refers to two (anti)quarks that form a spatially compact cluster, rather than the color-spin-isospin configuration of two (anti)quarks as discussed in Ref. <cit.>.
	
	
	The root-mean-square (rms) radius is a commonly used criterion to analyze the spatial structures of tetraquark states and distinguish between different tetraquark configurations. However, we found that the conventional rms radii calculated using the complete wave function could be misleading and fail to identify molecular configuration due to the antisymmetrization of identical particles <cit.>. In order to eliminate the ambiguity, we proposed a new approach to calculate the rms radius. For systems without identical particle such as bcs̅q̅, such an ambiguity does not exist and we can calculate the rms radii using the complete wave function directly. For systems with one pair of identical particles ( bcq̅q̅ (qqb̅c̅) ,  bcs̅s̅ (ssb̅c̅) , ccs̅q̅, bbs̅q̅), we decompose the complete antisymmetric wave function as 
	
    Ψ(θ)=   [(q_1q̅'̅)_1_c(q_2q̅”̅)_1_c]_1_c⊗|ψ_1(θ)⟩
       +[(q_2q̅'̅)_1_c(q_1q̅”̅)_1_c]_1_c⊗|ψ_2(θ)⟩
    
    			=   𝒜 [(q_1q̅'̅)_1_c(q_2q̅”̅)_1_c]_1_c⊗|ψ_1(θ)⟩
    ≡   𝒜 Ψ_ nA(θ).

	For systems with two pairs of identical particles (ccq̅q̅, bbq̅q̅, ccs̅s̅, bbs̅s̅), we decompose the complete antisymmetric wave function as
	
    Ψ(θ)=   ∑_s_1≥ s_2([(Q_1q̅_1)^s_1_1_c(Q_2q̅_2)^s_2_1_c]^S_1_c⊗|ψ_1^s_1s_2(θ)⟩.
       +[(Q_1q̅_1)^s_2_1_c(Q_2q̅_2)^s_1_1_c]^S_1_c⊗|ψ_2^s_1s_2(θ)⟩
       +[(Q_1q̅_2)^s_1_1_c(Q_2q̅_1)^s_2_1_c]^S_1_c⊗|ψ_3^s_1s_2(θ)⟩
       +.[(Q_1q̅_2)^s_2_1_c(Q_2q̅_1)^s_1_1_c]^S_1_c⊗|ψ_4^s_1s_2(θ)⟩)
    
    			=   𝒜∑_s_1≥ s_2[(Q_1q̅_1)^s_1_1_c(Q_2q̅_2)^s_2_1_c]^S_1_c⊗|ψ_1^s_1s_2(θ)⟩
    ≡   𝒜 Ψ_ nA(θ),

	where s_1, s_2 sum over spin configurations with total spin S. Instead of using the complete wave function Ψ(θ), we use the decomposed non-antisymmetric wave function Ψ_nA(θ) to calculate the rms radius:
	
    r^rms_ij≡Re[√(⟨Ψ_nA(θ) | r_ij^2 e^2iθ|Ψ_nA(θ)⟩/⟨Ψ_nA(θ) | Ψ_nA(θ)⟩)].

 For a meson molecule, r^rms_13 and r^rms_24 are expected to be the sizes of the constituent mesons, and much smaller than the other rms radii. The novel definition of the rms radius is useful for distinguishing the molecular configuration from the others.  For a compact tetraquark state, it may seem more reasonable to calculate the rms radii using the complete wave function. However, we find that the results from the conventional definition and the novel definition are qualitatively the same for compact tetraquarks, as shown in Appendix <ref>. Therefore, we will only use the novel definition to calculate the rms radii. 
	
The inner products in the CSM are defined using the c-product <cit.>, 
	
    ⟨ϕ_n |ϕ_m⟩≡∫ϕ_n(r)ϕ_m(r)d^3r,

where the complex conjugate of the “bra" state is not taken. The c-product definition ensures that the expectation values of physical quantities are independent of the complex scaling angle  θ. The rms radius calculated by the c-product is generally not real, but its real part can still reflect the internal quark clustering behavior if the width of the resonant state is not too large, as discussed in Refs. <cit.>.




§ RESULTS AND DISCUSSIONS



 §.§ QQ^(')q̅q̅

With the CSM, we calculate the complex eigenenergies of the S-wave QQ^(')q̅q̅ systems. The results for ccq̅q̅, bbq̅q̅ and bcq̅q̅ are shown in Figs. <ref>, <ref> and <ref>, respectively. We only focus on the energy spectra below the M(1S)M'(2S) dimeson thresholds. We identify resonant states by varying complex scaling angles θ. We obtain a series of QQ^(')q̅q̅ bound and resonant states, which are labeled as T_QQ^('),I(J)(M) in the following discussions, where M represents the mass of the state. The complex energies, proportions of different color configurations and rms radii of these states are summarized in Table <ref>.
	
	

  §.§.§ ccq̅q̅

	
	
	For the ccq̅q̅ system, we obtain one bound state, three resonant states with quantum number I(J^P)=0(1^+) and a resonant state with I(J^P)=1(2^+). The isoscalar bound state T_cc,0(1)(3864) is located 14 MeV below the D^*D threshold. It has a molecular configuration, with r_c_1q̅_1^ rms and r_c_2q̅_2^ rms close to the sizes of D^* and D mesons, respectively, indicating that it is a D^*D molecular state. Comparing with the experimentally observed T_cc(3875)^+, whose binding energy is around -300 keV and characteristic size is around 7 fm <cit.>, T_cc,0(1)(3864) is lower in energy and smaller in size. However, considering that the uncertainties of quark model are up to tens of MeV, it may still serve as a candidate for T_cc(3875)^+. The discrepancies suggest that improvements of the constituent quark model are needed to better describe the long-range interactions between hadrons. Additionally, in Section <ref>, we calculate the rms radii of T_QQ,0(1) with different binding energies by varying the heavy quark mass. When the binding energy is adjusted to match that of the experimental T_cc(3875)^+, the distance between the (anti)quarks in D and D^* is approximately 6 fm. This result is in excellent agreement with the characteristic size estimated in Ref. <cit.>, further supporting the interpretation of the T_cc(3875)^+ as a very loosely bound molecular state.
	
	We obtain an isoscalar resonant state T_cc,0(1)(4031) with a width of Γ=54 MeV, located 1 MeV below the D^*D^* threshold. This state is in agreement with the results from previous works using heavy quark spin symmetry <cit.>, constituent quark model <cit.> and lattice QCD <cit.>. It can be searched for in the D^*D and DDπ channel in future experiments. One may expect that T_cc,0(1)(4031) is the partner of the D^*D molecular bound state T_cc,0(1)(3864) and have a D^*D^* molecular configuration. However, the rms radii for T_cc,0(1)(4031) is highly numerically unstable and change dramatically as θ varies. The numerical instability of the rms radii may result from the state being sandwiched between the D^*D and D^*D^* thresholds, which makes it strongly coupled to the scattering states. Calculations with higher numerical precision are needed to identify its spatial structure.
	
	We also obtain three higher resonant states T_cc,0(1)(4466), T_cc,0(1)(4542), and T_cc,1(2)(4673), which have compact tetraquark configurations.  The state T_cc,1(2)(4673) has 10% isovector (q̅q̅)_6_c^1 component, where the spatial wave function between two light antiquarks must be antisymmetric. Similar configurations also appear in other 2^+ QQ^(')q̅q̅ states. Antisymmetric spatial wave function is possible in our S-wave calculations since we only restrict the total orbital angular momentum to be S-wave. The states T_cc,0(1)(4466) and T_cc,0(1)(4542) can decay into the D^*D and D^*D^* channels while the isovector state can only decay into the D^*D^* channel. These three states were not found in Ref. <cit.>, where the authors used the same constituent quark model as we do. The reason for the discrepancies may be that we use a larger set of wave function bases compared to the previous work, as discussed in Ref. <cit.>. Future experimental explorations of these states may help resolve the discrepancies between different methods.
	
	

  §.§.§ bbq̅q̅

	We first focus on the I(J^P)=0(1^+) bbq̅q̅ system, which is the bottom analogue of the T_cc state. We obtain a deeply bound state T_bb,0(1)(10491) with a binding energy of -153 MeV, and a shallow bound state T_bb,0(1)(10642) with a binding energy of -1 MeV. The lower state is a compact diquark-centered tetraquark, which is similar to the helium atom. The two bottom quarks form a very compact diquark cluster, with the rms radius between them (r_b_1b_2^ rms=0.33 fm) being approximately the size of the bottomonia, as listed in Table <ref>. The two light antiquarks orbit around the diquark, similar to the way two electrons orbit around the nucleus in a helium atom. The rms radii between the bottom quarks and light antiquarks r_bq̅^ rms are around 0.7 fm. The dominant color configuration of the state is χ_3̅_c⊗ 3_c, where the color electric interactions are attractive between two (anti)quarks. Note that the χ_3̅_c⊗ 3_c component reaches 97%! The strong attraction between two bottom quarks contribute to the deep binding energy of the state.  A deeply bound state in the I(J^P)=0(1^+) bbq̅q̅ system has been anticipated since the 1980s <cit.> and was recently predicted in lattice QCD studies <cit.>. Our findings agree with the previous works. The bound state T_bb,0(1)(10491) is below the B̅B̅ threshold and can only decay weakly. 
	
	The higher bound state T_bb,0(1)(10642) has a molecular configuration, with r_b_1q̅_1^ rms and r_b_2q̅_2^ rms being close to the sizes of B̅^* and B̅ mesons, respectively, indicating that it is a B̅^*B̅ molecular state. It is the bottom analogue of the D^*D molecular bound state T_cc,0(1)(3864). The bound state T_bb,0(1)(10642) can decay radiatively to B̅B̅γ.
	
	

	In addition to bound states, we also obtain six isoscalar bbq̅q̅ resonant states with J^P=1^+. The state T_bb,0(1)(10700) is located near the B̅^*B̅^* threshold and identified as a B̅^*B̅^* molecular state. It is dominated by the B̅^*B̅^* component, resulting in a small decay width of only 2 MeV to the B̅^*B̅ channel. In contrast, its charmed partner T_cc,0(1)(4031) is a relatively broad state with a width of 54 MeV. This discrepancy may indicate that the nature of the two states is different, namely T_cc,0(1)(4031) may have sizable contributions from both the D^*D and D^*D^* channels. The state T_bb,0(1)(10700) can be searched for in the B̅^*B̅ channel.
	
	The state T_bb,0(1)(11025) is a compact diquark-antidiquark tetraquark state. Comparing its internal structures with those of the bound state T_bb,0(1)(10491), we find an interesting resemblance. Both states are dominated by the χ_3̅_c⊗3_c color configuration. The rms radii between two bottom quarks r^ rms_b_1b_2=0.33 fm are the same for these two states, while the rms radii between the bottom quarks and light antiquarks r^ rms_bq̅ are around 1 fm for T_bb,0(1)(11025), larger than those for T_bb,0(1)(10491). This clearly suggests that T_bb,0(1)(11025) is the radial excitation in the light degree of freedom of T_bb,0(1)(10491). The state T_bb,0(1)(11025) can decay into the B̅^*B̅ and B̅^*B̅^* channels.
	
	Two isoscalar bound states T_bb,0(0)(11195) and T_bb,0(2)(11370) are found in the 0^+ and 2^+ systems, respectively. The T_bb,0(0)(11195) does not decay into the S-wave BB pair and T_bb,0(2)(11370) does not decay into the S-wave B^*B^* pair, since the total wave functions of two identical bosons must be interchange symmetric. Hence, their lowest S-wave dimeson thresholds are M(1S) M'(2S). However, it should be noted that P-wave dimeson channels exist below these states. The coupling between the P-wave channels and these states may alter their positions and nature.
	
	In contrast to the isoscalar systems, no bound state is obtained in the isovector bbq̅q̅ systems. The absence of isovector bound state suggests that the configuration (q̅q̅)^0_3_c with isospin 0, which is referred to as “good antidiquark" in Ref. <cit.>, is essential for forming the doubly bottomed tetraquark bound states. We obtain two resonant states with J^P=0^+, four resonant states with J^P=1^+, and three resonant states with J^P=2^+. They are identified as either compact diquark-centered tetraquarks or compact even tetraquarks. The two lowest states T_bb,1(1)(10685) and T_bb,1(2)(10715) may be more likely to be observed in experiments. It is worth noting that in Ref. <cit.>, the authors used a different constituent quark model including the one-pion exchange interaction, and identified an isovector resonant state with J^P=1^+, whose mass and width are in agreement with the state T_bb,1(1)(10685) obtained in our calculations. We recommend experimental exploration of this state in the B̅^*B̅ channel.

	

	

  §.§.§ bcq̅q̅

	We identify three isoscalar bound states T_bc,0(0)(7129), T_bc,0(1)(7185), and T_bc,0(2)(7363) in the bcq̅q̅ systems. T_bc,0(0)(7129) and T_bc,0(1)(7185) are compact even tetraquark states, where the rms radii between four (anti)quarks are of similar size. Unlike the χ_3̅_c⊗3_c dominated compact diquark-centered tetraquark T_bb,0(1)(10491), these two compact even tetraquarks have sizable components of both χ_3̅_c⊗3_c and χ_6_c⊗6̅_c configurations. This suggests that both the interactions between two (anti)quarks and the interactions between quarks and antiquarks are important for these two states. The shallow bound state T_bc,0(2)(7363) with a binding energy of -3 MeV is identified as a B̅^*D^* molecular state. The scalar bound state T_bc,0(0)(7129) can only decay weakly, while T_bc,0(1)(7185) can decay radiatively to B̅ Dγ, and T_bc,0(2)(7363) can decay strongly to B̅^*Dπ. In contrast to the isoscalar case, no isovector bcq̅q̅ bound state is found. This again shows that the (q̅q̅)^0_3_c configuration plays a vital role in forming doubly heavy tetraquark bound states, as in the bbq̅q̅ systems. In Refs. <cit.>, the authors investigated possible bcq̅q̅ bound states from the molecular picture. It was shown in  Ref. <cit.> that for the isovector system the interactions between charmed and bottomed mesons are repulsive and no bound state exists, which is consistent with our findings. Moreover, both studies <cit.> suggested that shallow bound states with binding energies of several MeV may exist in the I(J^P)=0(0^+),0(1^+), and 0(2^+) bcq̅q̅ systems. For comparison, we also obtain a loosely bound molecular state in the isoscalar 2^+ system. But for the isoscalar 0^+ and 1^+ systems, our quark model calculations lead to bcq̅q̅  compact tetraquark states with binding energies around -30 MeV. 
	
	We also obtain a series of resonant states with different quantum numbers. The two lowest states T_bc,0(0)(7301) and T_bc,1(2)(7430) may be more likely to be observed in experiments. They can be searched for in the B̅ D and B̅^*D^* channels, respectively. Moreover, in the I(J^P)=0(1^+) bcq̅q̅ system, we find that the scattering states of B̅D^* and B̅^*D^* deviate from the continuum lines, indicating possible existence of a resonant state, which may be the bcq̅q̅ analogue of T_cc,0(1)(4031) and T_bb,0(1)(10700). However, we cannot determine this state accurately due to the limitations of numerical precision in the present calculations. 
	
	

	
	

 §.§ QQ^(')s̅q̅

	

	

	The complex eigenenergies of the S-wave ccs̅q̅, bbs̅q̅, and bcs̅q̅ systems are shown in Figs. <ref>, <ref> and <ref>, respectively. We obtain a series of QQ^(')s̅q̅ bound and resonant states, which are labeled as T_QQ^(')s̅,J(M) in the following discussions, where M represents the mass of the state. The complex energies, proportions of different color configurations and rms radii of these states are summarized in Table <ref>.
	
	We obtain two resonant states in the ccs̅q̅ system. They are located near the M(1S)M'(2S) dimeson threshold and have compact even tetraquark configurations. T_ccs̅,1(4578) can decay into the D_s^*D, D_sD^*, and D_s^*D^* channel while T_ccs̅,2(4723) can only decay into the D_s^*D^* channel.
	
	For the bbs̅q̅ system, we identify a bound state T_bbs̅,1(10647) with a binding energy of -64 MeV. It is a deeply bound state with compact diquark-centered tetraquark configuration, and can be considered as the strange partner of T_bb,0(1)(10491). It is below the B̅_sB̅ threshold and can only decay weakly. Moreover, we also obtain four resonant states with  J^P=1^+. The state T_bbs̅,1(10766) is located 3 MeV below the B̅_s^*B̅^* threshold and have a compact diquark-centered tetraquark configuration. It is the radial excitation in the light degree of freedom of T_bbs̅,1(10647), similar to T_bb,0(1)(11025) being the radial excitation of T_bb,0(1)(10491). It can decay into the B̅_s^*B̅ and B̅_sB̅^* channel, while the three higher resonant states can also decay into the B̅_s^*B̅^* channel. In addition, we obtain two resonant states with J^P=0^+ and three resonant states with J^P=2^+. The 0^+ states can decay into the B̅_sB̅ and B̅_s^*B̅^* channels while the 2^+ states can only decay into the B̅_s^*B̅^* channel.
	
	For the bcs̅q̅ system, we obtain two resonant states each in the 0^+, 1^+, and 2^+ systems. Similar to the resonant state T_cc,0(1)(4031), T_bcs̅,1(7437) is sandwiched between two thresholds and its rms radii results are numerically unstable. The lowest 0^+ resonant state can decay into the B̅_sD and B̅D_s channels. The lowest 1^+ resonant state can decay into the B̅_s^*D, B̅^*D_s, B̅_sD^*, and B̅D_s^* channels. The lowest 2^+ resonant state can decay into the B̅_s^*D^* and B̅^*D_s^* channels.
	

	
	
	

	
	
	
	

 §.§ QQ^(')s̅s̅

	

	The complex eigenenergies of the S-wave ccs̅s̅, bbs̅s̅, and bcs̅s̅ systems are shown in Figs. <ref>, <ref> and <ref>, respectively. We obtain a series of QQ^(')s̅s̅ resonant states, which are labeled as T_QQ^(')s̅s̅,J(M) in the following discussions, where M represents the mass of the state. The complex energies, proportions of different color configurations and rms radii of these states are summarized in Table <ref>.
	
	The QQ^(')s̅s̅ systems and the isovector QQ^(')q̅q̅ systems share the same internal symmetries, and their energy spectra bear a resemblance. In the ccs̅s̅ system, we obtain a resonant state T_ccs̅s̅,2(4808) below the D_s^*(2S)D_s^* threshold. It can be considered as the strange partner of T_cc,1(2)(4673). In the bbs̅s̅ system, the two lowest resonant states T_bbs̅s̅,1(10846) and T_bbs̅s̅,2(10877) are the strange partners of T_bb,1(1)(10685) and T_bb,1(2)(10715), respectively. In the bcs̅s̅ system, two resonant states with J^P=2^+ are found. The lowest state T_bcs̅s̅,2(7605) is the strange partner of T_bc,1(2)(7430). The mass differences between the isovector QQ^(')q̅q̅ resonant states and their strange partners QQ^(')s̅s̅ resonant states are around 150 MeV. For higher resonant states, the correspondence between QQ^(')s̅s̅ states and isovector QQ^(')q̅q̅ states is less clear. Compared to the isovector QQ^(')q̅q̅ systems, the QQ^(')s̅s̅ systems have fewer resonant states near the M(1S)M'(2S) dimeson thresholds. All of the QQ^(')s̅s̅ states are identified as compact tetraquarks and can be searched for in corresponding dimeson decay channels.
	
	The state T_bcs̅s̅,2(7605) has a negative proportion of color configuration χ_6_c⊗6̅_c, which is rather confusing. However, it should be emphasized that the inner products calculated by the c-product in Eq. (<ref>) are not positive-definite and generally not real. Moreover, the physical quantities calculated by the c-product may only have probabilistic interpretation for resonant states that are not too broad <cit.>. Therefore, the negative proportions of the broad resonance T_bcs̅s̅,2(7605) should not be taken too seriously.
	
	
	
	

	

 §.§ Heavy quark mass dependence of T_QQ,0(1) bound states
 
	The  I(J^P)=0(1^+) QQq̅q̅ systems, in which both loosely bound molecular states and deeply bound compact states are observed in our calculations, have attracted great attention. We shall investigate very carefully how the properties of the bound states change with the heavy quark mass m_Q. The binding energies of the ground state and the first excited state with varying m_Q are shown in Fig. <ref>. We see that the ground state is bound for m_Q larger than around 1200 MeV, and the binding energy increases as m_Q increases. When m_Q reaches around 4600 MeV, the first excited state also becomes a bound state. The rms radii of the two bound states with varying m_Q are shown in Fig. <ref>. As m_Q increases, r^ rms_Q_1q̅_1 and r^ rms_Q_2q̅_2 remain stable and match the sizes of the ground state Qq̅ mesons with spin 0 and 1 (denoted as M and M^*) respectively, while the other four rms radii decreases significantly. For the ground state, the rms radii between the clusters Q_1q̅_1 and Q_2q̅_2 are larger than 1 fm when m_Q is less than around 2 GeV, indicating a M^*M molecular configuration. In this scenario, the ground state is a loosely bound state with a binding energy |Δ E |< 20 MeV. For a larger m_Q, the two heavy quarks form a compact diquark, and the ground state transforms into a deeply bound compact diquark-centered tetraquark state. For the first excited state, the rms radii between the clusters Q_1q̅_1 and Q_2q̅_2 decreases more slowly compared to the ground state. The first excited state remains a molecular state for m_Q < 6 GeV.   
	
	The physical interpretation for the above findings is straightforward. As m_Q increases, the kinetic energies of the heavy quarks are suppressed. The distance between two heavy quarks gets smaller and the attractive color electric interaction gets stronger, leading to a deeper bound state and a potential second bound state. Interestingly, a similar conclusion on the existence of two types of tetraquark bound states in the limit of large heavy quark mass was reached in Ref. <cit.> using the Born-Oppenheimer approximation within the framework of large N QCD. Future experimental explorations of the two T_bb,0(1) bound states may help test theoretical predictions.
	
	

	



	

§ SUMMARY

	In summary, we calculate the energy spectra of the S-wave doubly heavy tetraquark systems QQ^(')q̅q̅, QQ^(')s̅q̅, and QQ^(')s̅s̅ (Q^(')=b,c) using the AL1 quark potential model.  We apply the complex scaling method to study possible bound states and resonant states simultaneously, and the Gaussian expansion method to solve the four-body Schrödinger equation. We focus on the low-lying states below the lowest M(1S)M'(2S) dimeson threshold. The uncertainties of these tetraquark states are expected to be of the same order as those of the 1S mesons, which are around tens of MeV. 	
	
	We obtain bound states in the ccq̅q̅, bbq̅q̅, bcq̅q̅, and bbs̅q̅ systems. The shallow bound state T_cc,0(1)(3864) serves as a candidate for the experimental T_cc(3875)^+ state. The bound state T_bc,0(2)(7363) can decay strongly to B̅^*Dπ. The bound states T_bb,0(1)(10642) and T_bc,0(1)(7185) can decay radiatively. The bound states T_bb,0(1)(10491), T_bc,0(0)(7129), and T_bbs̅,1(10647) can only decay weakly. In addition, a series of doubly heavy resonant states are found. We urge future experimental explorations of these predicted states.
	
    We use the rms radii to distinguish between meson molecular states and compact tetraquark states. The compact tetraquark states are further classified into three different configurations: compact even tetraquark, compact diquark-antidiquark tetraquark and compact diquark-centered tetraquark. The shallow bound states T_cc,0(1)(3864), T_bb,0(1)(10642), and T_bc,0(2)(7363) have molecular configurations, which are QCD molecules. 
    The deeply bound states T_bb,0(1)(10491) and T_bbs̅,1(10647) are compact diquark-centered tetraquarks, which are coined as the “QCD Helium atom" in Ref. <cit.>. The T_bc,0(0)(7129) and  T_bc,0(1)(7185) are compact even tetraquarks, which are ideal candidates of the “QCD Hydrogen molecule" as noted in Ref. <cit.>. 
 
 Most of the resonant states are compact tetraquark states, except that T_bb,0(1)(10700) is a B̅^*B̅^* molecular state. The resonant states T_bb,0(1)(11025) and T_bbs̅,1(10766) are considered as the radial excitations in the light degree of freedom of the bound states T_bb,0(1)(10491) and T_bbs̅,1(10647), respectively. It is worth noting that all of the compact diquark-centered tetraquarks and compact diquark-antidiquark tetraquarks identified in our calculations are dominated by the χ_3̅_c⊗3_c color configuration, except for the broad resonant state T_bc,0(0)(7301). In these states, the attractive color electric interactions between two heavy quarks play an important role. On the other hand, mixing effect between χ_3̅_c⊗3_c and χ_6_c⊗6̅_c configurations is important in compact even tetraquarks and meson molecules. Similar classifications of tetraquark configurations were made in Ref. <cit.>. The compact diquark-centered tetraquarks and compact diquark-antidiquark tetraquarks resemble the type-1 tetraquarks in Ref. <cit.>, and the compact even tetraquarks and meson molecules resemble the type-2 tetraquarks. The classifications of tetraquarks based on their color-spatial configurations help unravel the rich internal structures and various forming mechanisms of tetraquark states. 
	
	We also explore the heavy quark mass dependence of the T_QQ,0(1) bound states. As the heavy quark mass increases from 1.2 GeV to 6 GeV, the ground state transforms from a loosely bound molecular state to a deeply bound compact diquark-centered tetraquark state, with the emergence of a second loosely bound molecular state.  The future experimental explorations of the two T_bb,0(1) bound states may help test theoretical predictions and deepen our understanding of quantum chromodynamics.
	
	
	
	
	

§ ACKNOWLEDGMENTS

	
	We thank Zi-Yang Lin for the helpful discussions. This project was supported by the National
	Natural Science Foundation of China (No. 12475137). This
	project was also funded by the Deutsche Forschungsgemeinschaft (DFG,
	German Research Foundation, Project ID 196253076-TRR 110). The computational resources were supported by High-performance Computing Platform of Peking University.
	
	

§ TWO DEFINITIONS OF ROOT-MEAN-SQUARE RADIUS 

	In our calculations, we use the decomposed non-antisymmetric wave function to calculate the rms radii. It seems more reasonable to calculate the rms radii of compact tetraquark states using the complete wave function, which considers contributions from both direct terms and exchange terms. However, our primary interest lies in the general clustering behavior of tetraquark states, rather than in specific numerical results of the rms radii, which are not experimentally observables at present. The rms radii calculated using the decomposed non-antisymmetric wave function are already capable of distinguishing between different tetraquark configurations.  To illustrate this, we compare the results of rms radii calculated using the complete wave function Ψ and the decomposed non-antisymmetric term Ψ^J_ nA in Table <ref>. We take four states with different configurations as examples. For the compact tetraquark states T_bc,0(1)(7185), T_cc,0(1)(4466), and T_bb,0(1)(10491), the results from Ψ and Ψ_ nA are similar. We can draw the same conclusion on their spatial structures from both results. However, for the molecular state T_bb,0(1)(10642), the results from Ψ_ nA can clearly demonstrate the clustering behaviour of a molecular state, while the results from  Ψ are more ambiguous due to the antisymmetrization. 
	
	
	In conclusion, the novel definition of rms radii, which are calculated using only the decomposed non-antisymmetric wave function, can reflect the internal spatial structure of tetraquark states more transparently. 
	
	




	
	
