Benchmark calculations of fully heavy compact and molecular tetraquark states

We calculate the mass spectrum of the S-wave fully heavy tetraquark systems $ QQ\bar Q\bar Q~(Q=c,b) $ with both normal $ (J^{PC}=0^{++},1^{+-},2^{++}) $ and exotic $ (J^{PC}=0^{+-},1^{++},2^{+-}) $ C-parities using three different quark potential models (AL1, AP1, BGS). The exotic C-parity systems refer to the ones that cannot be composed of two S-wave ground heavy quarkonia. We incorporate the molecular dimeson and compact diquark-antidiquark spatial correlations simultaneously, thereby discerning the actual configurations of the states. We employ the Gaussian expansion method to solve the four-body Schr\"odinger equation, and the complex scaling method to identify the resonant states. The mass spectra in three different models qualitatively agree with each other. We obtain several resonant states with $ J^{PC} = 0^{++}, 1^{+-}, 2^{++}, 1^{++} $ in the mass region $(6.92,7.30)\, \mathrm{GeV}$, some of which are good candidates of the experimentally observed $X(6900)$ and $X(7200)$. We also obtain several exotic C-parity zero-width states with $ J^{PC}=0^{+-} $ and $ 2^{+-} $. These zero-width states have no corresponding S-wave diquarkonium threshold and can only decay strongly to final states with P-wave quarkonia. With the notation $T_{4Q,J(C)}(M)$, we deduce from the root mean square radii that the $ X(7200) $ candidates $ T_{4c,0(+)}(7173), T_{4c,2(+)}(7214) $ and the state $ T_{4c,1(-)}(7191) $ look like molecular states although most of the resonant and zero-width states are compact states.

The nature of fully charmed tetraquark states remains controversial.Currently, few theoretical works consider compact diquark-antidiquark and molecular diquarkonium spatial configurations simultaneously and perform comprehensive fourbody dynamical calculations.In our previous studies [80][81][82], we incorporated both dimeson and diquark-antidiquark spatial configurations, employing various quark models and few-body methods to conduct benchmark calculations for tetraquark bound states.We have illustrated that the Gaussian expansion method (GEM) is highly efficient in exploring tetraquark states.Our results indicate that there are no bound states in the fully heavy tetraquark systems.In Refs.[70,76], the authors utilized the GEM to conduct four-body dynamic calculations for the fully charmed tetraquark systems.Furthermore, they employed the complex scaling method (CSM) [83][84][85] to distinguish resonant states from dicharmonium scattering states, yielding convincing and intriguing results.However, it is noteworthy that the discussion in Refs.[70,76] is limited to fully charmed tetraquark states, and does not include an exploration of fully bottomed tetraquark states.Furthermore, the discussion excludes the systems with exotic Cparity (J P C = 0 +− , 1 ++ , 2 +− ), which refer to the ones that cannot be composed of two S-wave ground heavy quarkonia.
In this study, we aim to investigate the S-wave fully heavy tetraquark QQ Q Q (Q = b, c) resonances with all possible quantum numbers, employing a framework that has been used to investigate Qsq q states efficiently [86].We employ the GEM to solve the four-body Schrödinger equation, considering both compact diquark-antidiquark and molecular diquarkonium spatial configurations.We utilize the CSM to distinguish resonant states from diquarkonium scattering states.Regarding the (anti)quark-(anti)quark interactions, we compare three quark potential models with well-determined parameters and do not introduce any new free parameters.Additionally, we apply the approach proposed in our previous In the center-of-mass frame, the nonrelativistic Hamiltonian for a tetraquark system reads where m i and p i are the mass and momentum of the (anti)quark i, respectively.V ij represents the two-body interaction between the (anti)quark pair (ij).In this study, we adopt three different quark potential models, namely the AL1 and AP1 potentials proposed in Refs.[87,88] and the potential used to study charmonia in Ref. [89] (denoted as the BGS potential hereafter).These three potentials contain the one-gluon-exchange interaction and quark confinement interaction, and can be written as where r ij is the distance between (anti)quark i and j, λ i and S i are the SU(3) color Gell-Mann matrix and the spin operator acting on (anti)quark i.The first term in the potential is referred to as the color electric term, and the last term is known as the color magnetic term.In the conventional meson and baryon systems, the color factor λ i • λ j is always negative and induces a confining interaction, thus all eigenstates of the Hamiltonian must be bound states.However, scattering states of two color-singlet clusters and possible resonant states are allowed in the tetraquark systems, since they have richer inner color structures than the conventional hadrons, and the color factor might take a zero or positive value.The parameters of the models are taken from Refs.[88,89] and listed in Table I.The parameters for the AL1 and AP1 potential were determined by fitting the meson spectra across all flavor sectors, while those for the BGS potential were determined by the charmonium spectra.The theoretical masses of the charmonia and bottomonia as well as their root-mean-square (rms) radii calculated from the AP1 potential are listed in Table II.It can be seen that all three potential models can give a satisfactory description of the meson spectra.We use the complex scaling method (CSM) to obtain possible bound and resonant states simultaneously, and we apply the Gaussian expansion method (GEM) to solve the complexscaled four-body Schrödinger equation.We also introduce a method to determine the C-parity of the neutral tetraquark states by decomposing the Hilbert space.
In the CSM [83][84][85], the coordinate r and its conjugate momentum p are transformed as Under such a transformation, the complex-scaled Hamiltonian which is no longer hermitian and has complex eigenvalues E(θ).According to the ABC theorem [83,84], the energies of the bound states, resonant states and scattering states can be obtained as the eigenvalues of H(θ) simultaneously.The bound states (zero-width states) lie on the negative real axis in the energy plane and are not changed by the complex scaling.The resonant states with mass M R and width Γ R can be The scattering states line up along beams starting from threshold energies and rotated clockwise by 2θ from the positive real axis.
To solve the complex-scaled Schrödinger equation, we apply the GEM [91] and expand the wave functions of the Swave fully heavy tetraquark states with total angular momentum J and C-parity C as where A is the antisymmetric operator of identical particles.χ J α is the color-spin wave function, given by for all possible combinations of s 1 , s 2 , J. It should be emphasized that one has the flexibility to use either the diquarkantidiquark type color-spin basis in Eq. ( 6) or the dimeson type color-spin basis, which can be written as These two sets of bases are both complete and the transformation between them are shown in Appendix A. It has been demonstrated that the use of different discrete basis functions yields negligible differences once they are complete [82].The S-wave spatial wave function Φ (jac) n1,n2,n3 is written as where (jac) = (a), (b), (c) denotes three sets of spatial configurations (dimeson and diquark-antidiquark) considered in our calculations, and r jac , λ jac , ρ jac are three independent Jacobian coordinates in configuration (jac), as shown in Fig. 1. ϕ ni (r) takes the Gaussian form, where N ni is the normalization factor.Finally, the expansion coefficients C (jac) α,ni (θ) are determined by solving the energy eigenvalue equation, For the fully heavy neutral tetraquark system (QQ Q Q), the eigenstates of the Hamiltonian have definite C-parity C. To determine the C-parity of the obtained states, we decompose the Hilbert space H into H + and H − , where H ± represents the subspace with C = ±1.Under charge conjugation, the basis functions with (jac) = (a) are transformed as In the third line, we exchange xc in the color-spin wave function and the factor (−1) s1+s2−J arises from the Clebsch-Gordan coefficient.In the last line the particle indices 1 ↔ 3, 2 ↔ 4 are exchanged to rewrite the basis function to its original form.The transformation behavior of the basis functions with (jac) = (b), (c) can be obtained similarly, For these two sets of spatial configurations, the indices of Gaussian basis (n FIG. 1.The Jacobian coordinates for two types of spatial configurations: (a), (b) for the dimeson configurations, and (c) for the diquarkantidiquark configuration.

C. Discern between molecular and compact tetraquark states
In contrast to other frameworks, such as those discussed in Refs.[92,93], the quark model does not require a prior assumptions regarding the structures of a multiquark state.Basically, the molecular or compact tetraquark states can be discerned through the analysis of their wave functions in the quark model.The proportions of color components and root-mean-square (rms) radii are commonly used criteria [27,86,94], which reflect the color structure and the spatial structure, respectively.
However, identifying molecular states based on the proportions of color components may be ambiguous and misleading in the systems with identical particles.An ideal loose molecular state consists of two colorless subclusters, which are widely separated.Due to the antisymmetrization of identical particles, its wave function takes the form The two subclusters are widely separated, therefore we have ⟨ψ 1 |ψ 2 ⟩ ≈ 0. We decompose the wave function in Eq. ( 14) in the orthogonal color basis as with For an ideal loose molecular state, we obtain P 1⊗1 : P 8⊗8 = 5 : 4 and P3 ⊗3 : P 6⊗ 6 = 1 : 2. One can see that even for an ideal molecular state, the proportions of χ 8c⊗8c and χ 1c⊗1c are comparable.Therefore, discerning between molecular and compact states solely via the proportions of color components could be misleading in the systems with identical particles.Meanwhile, the naive definition of the rms radius could be misleading when the antisymmetric wave function is required for the identical quarks.For instance, when the mesons (cq) and (sq) form a molecular state, the wave function satisfying the Pauli principle is |ψ A ⟩ = |(cq 1 )(sq 2 )⟩−|(cq 2 )(sq 1 )⟩.One can see that each antiquark q belongs to both mesons simultaneously.Therefore, neither ⟨ψ A |r 2 cq |ψ A ⟩ nor ⟨ψ A |r 2 sq |ψ A ⟩ can reflect the size of the constituent mesons.
In Ref. [86], we proposed a new approach to calculate the rms radii in the Qsq q system, which eliminates the ambiguity arising from the antisymmetrization of identical particles q q.
In the fully heavy tetraquark system, there exist two pairs of identical particles.Here we further extend the definitions of rms radii for the QQ Q Q system.We uniquely decompose the antisymmetric wave function as where s 1 , s 2 sum over spin configurations with total angular momentum J.We denote the non-antisymmetric component of the wave function as (21) where (Q 1 Q3 ) and (Q 2 Q4 ) form color singlets.Instead of using the complete wave function Ψ J (θ), we use |Ψ J nA (θ)⟩ to define the rms radius: This definition discards the contribution of the exchange terms resulting from the antisymmetrization, which may play an important role in the numerical results of the rms radii, especially for compact tetraquark states.However, our primary interest lies in the general clustering behavior of the tetraquark states rather than in specific numerical results.The current definition of the rms radius is useful for investigating the internal spatial structures of the tetraquark states.For example, if the resulting state is a scattering state or a hadronic molecule of η c J/ψ, r c1 c3 and r c2 c4 are respectively expected to be the sizes of J/ψ and η c , and much smaller than the other rms radii.On the other hand, if the resulting state is a compact tetraquark state, all rms radii in the four-body system should be of the same order.One can find detailed examples for comparing two definitions of the rms radius in Appendix C.However, it should be noted that for a hadronic molecule composed of two mesons with the same quantum numbers but different radial excitation, for example ψ(3S)J/ψ, the current definition cannot eliminate the ambiguity arising from the antisymmetrization.As a result, neither r c1 c3 nor r c2 c4 reflects the size of ψ(3S) or J/ψ; instead, they represent the average of the sizes of the two mesons.
It should be emphasized that the inner products in the CSM are defined using the c-product [95], where the square of the wave function rather than the square of its magnitude is used.The rms radius calculated by the cproduct is generally not real, but its real part can still reflect the internal quark clustering behavior if the resonant state is not too broad, as discussed in Ref. [96].

III. RESULTS AND DISCUSSIONS
We investigate the S-wave fully charmed cccc and fully bottomed bb bb tetraquark systems with all possible quantum numbers, including J P C = 0 ++ , 1 +− , 2 ++ , 0 +− , 1 ++ , 2 +− .It should be stressed that the S-wave ground state diquarkonium thresholds exist only in the 0 ++ , 1 +− , 2 ++ systems, namely η c η c and η b η b in the 0 ++ systems, η c J/ψ and η b Υ in the 1 +− systems, J/ψJ/ψ and ΥΥ in the 2 ++ systems.In the following discussions, the 0 ++ , 1 +− , 2 ++ systems are referred to as normal C-parity systems, whereas the 0 +− , 1 ++ , 2 +− systems are referred to as exotic C-parity systems.For convenience, we label the QQ Q Q tetraquark states obtained in the calculations as T 4Q,J(C) (M ), where M is the mass of the state.
A. Fully charmed tetraquark

States with normal C-parity
With the CSM, the complex eigenenergies of the 0 ++ , 1 +− , 2 ++ cccc systems obtained from three different quark potential models are shown in Fig. 2. We choose varying complex scaling angles θ to distinguish resonant states from scattering states.All of the states are above the lowest diquarkonium threshold, so no bound state is obtained.The diquarkonium scattering states rotate along the continuum lines starting from the threshold energies.Moreover, we obtain a series of resonant states whose complex energies are summarized in Table III.For comparison, we also list the results in Ref. [70], where the authors used the BGS potential for the cccc systems.
Qualitatively, the cccc resonances obtained from three different quark potential models are in accordance with each other.Most of the resonances exist in all three models.For a specific resonant state, its width remains consistent across different models, while the mass in the BGS potential is approximately 50-100 MeV larger than those in the AL1 and AP1 potentials.These differences are expected considering that the discrepancies of the predictions of the heavy quarkonium mass spectra from various potentials are up to tens of MeV.
The tetraquark resonant states with different quantum numbers J P C exhibit a similar pattern.A lower resonant state with mass M ≈ 7000 MeV and width Γ ≈ 75 MeV, and a higher resonant state with mass M ≈ 7200 MeV and width Γ ≈ 50 MeV are obtained in the 0 ++ , 1 +− , 2 ++ systems.The lower 0 ++ state can decay into the η c η c , J/ψJ/ψ, η c (2S)η c and ψ(2S)J/ψ channels.Additionally, the higher 0 ++ state can decay into the η c (3S)η c and ψ(3S)J/ψ channels.The lower 2 ++ state can decay into the J/ψJ/ψ and ψ(2S)J/ψ channels.Additionally, the higher 2 ++ state can decay into the ψ(3S)J/ψ channel.Considering that the quark potential models have errors up to tens of MeV and that we have neglected the widths of the quarkonia in our calculations, the lower 0 ++ and 2 ++ states may serve as the candidates for the experimentally observed X(6900) state, while the higher 0 ++ and 2 ++ states may serve as the candidates for the experimentally observed X(7200) state.On the other hand, the lower 1 +− state can decay into the η c J/ψ, η c (2S)J/ψ and η c ψ(2S) channels.Additionally, the higher 1 +− state can decay into the η c ψ(3S) and η c (3S)J/ψ channels.The 1 +− states are not the candidates for X(6900) or X(7200) because they cannot decay into either the J/ψJ/ψ or ψ(2S)J/ψ channels.These states can be searched for in future experiments.
Moreover, we observe several narrow resonant states with different quantum numbers.These states are found in the mass region (6.92, 7.30) GeV.These narrow resonances can be searched for by experiments in the corresponding diquarkonium decay channels.However, we do not observe any signal for resonance in the mass region (6.2, 6.6) GeV, namely no candidate for X(6400) or X(6600) is found.
Comparing the results of resonant states obtained from the BGS potential with those of Ref. [70], our calculations can reproduce the previous results well within numerical uncertainty.In addition, we obtain three extra narrow resonant states in the BGS potential.The reason for these discrepancies might be that a set of complete color-spin basis is used in our calculations while some basis functions are neglected in Ref. [70].For example, in the 0 ++ system, the diquarkonium color configuration ) 8c ] 1c is not included in the previous calculations.These missing basis functions turn out to be crucial for the existence of the extra resonant states.
As mentioned above, the fully charmed tetraquark resonant states obtained from different models qualitatively agree with each other.Therefore, we only choose the results in the AP1 potential to analyze their inner structures.The proportions of the color configurations χ3 c ⊗3c and χ 6c⊗ 6c in the wave functions as well as the rms radii of the fully charmed resonant states are listed in Table IV.For most resonant states, their rms radii are approximately of the same size and less than 1 fm, indicating that they are compact tetraquark states.However, for the X(7200) candidate states T 4c,0(+) (7173), T 4c,2(+) (7214) and the state T 4c,1(−) (7191) , the r rms c1 c3 and r rms c2 c4 are much smaller than the other radii.Compared with the rms radii of the quarkonia listed in Table II, we observe that both r rms c1 c3 and r rms c2 c4 are larger than the rms radii of η c and J/ψ, and smaller than those of η c (3S) and ψ(3S), while r rms c1 c4 , r rms c2 c3 , r rms c1c2 and r rms c3 c4 are much larger than the rms radii of all mesons.These results suggest that these three states might have a molecular configuration.It should be noted that these resonant states have relatively large widths and are located close to the continuum line of the scattering states, therefore the results of their rms radii are less numerically accurate in the CSM and should be considered as qualitative estimates.

States with exotic C-parity
The complex eigenenergies of the fully charmed tetraquark systems with exotic C-parity (J P C = 0 +− , 1 ++ , 2 +− ) ob- tained from three different quark potential models are shown in Fig. 3.We obtain a series of resonant and zero-width states, whose energies are summarized in Table V.Similar to the systems with normal C-parity, the states with exotic C-parity obtained from different models qualitatively agree with each other.The masses of a specific state in the AL1 and AP1 potentials are nearly the same, while the mass in the BGS potential is around 50-100 MeV larger than the former ones.In the following we solely focus on the results in the AP1 potential.The proportions of the color configurations and the rms radii of the resonant and zero-width states are listed in Table VI.The different rms radii of these states are approximately the same and less than 1 fm, indicating that all of these states have compact tetraquark configuration.
For the 0 +− cccc system, there do no exist any S-wave diquarkonium thresholds.Therefore no meson-meson scattering state is observed.All of the states in the system TABLE III.The complex energies E = M − iΓ/2 (in MeV) of the cccc resonant states with normal C-parity from various potential models.The last column lists the results in Ref. [70].The "?" indicates the potential existence of resonant states, which blend into the continuum lines of scattering states and cannot be obtained accurately in the present calculations.The "-" suggests that no corresponding resonance is obtained.
are identified as zero-width states, which lie on the real axis and are not changed by the complex scaling.According to the proportions of color configurations listed in Table VI, the four 0 +− zero-width states can be clearly arranged into two doublets, T 4c,0(−) (6890), T 4c,0(−) (6945) and T 4c,0(−) (7188), T 4c,0(−) (7242) .The color configurations of the two states inside a doublet are orthogonal to each other.The lower state is dominated by the χ 6c⊗ 6c component, while the higher state is dominated by the χ3 c ⊗3c component.In the fully heavy tetraquark system, it is known that the color magnetic term is suppressed by the heavy quark mass, and the dominant color electric interactions between two (anti)quarks are attractive in 6c and [ Q3 Q4 ]6 c configurations.Besides, the color electric term also provides an attractive interaction between the 6 c diquark and 6c antidiquark, which is much stronger than the one between 3c diquark and 3 c antidiquark due to the color SU(3) algebra [18,27,29].The fact that the lower state has the dominant χ 6c⊗ 6c configuration suggests that the strong attraction between two color sextet clusters prevails over the repulsion within the (anti)diquark and contributes to the formation of a deeper state than the χ3 c ⊗3c dominant one, which is consistent with the conclusion in Refs.[18,27].As a result of the interaction mechanism, the rms radii r rms c1c2 and r rms c3 c4 , which characterize the sizes of the diquark and antidiquark, take larger values in the χ 6c⊗ 6c dominant state than in the χ3 c ⊗3c dominant state.On the other hand, the rms radii r rms c1 c3 , r rms c2 c4 , r rms c1 c4 and r rms c2 c3 , which characterize the distance between diquark and antidiquark, take smaller values in the χ 6c⊗ 6c dominant state than in the χ3 c ⊗3c dominant state.It is also worth mentioning that the four rms radii r rms ci cj equal to each other in the 0 +− and 2 +− systems.The reason is that these states consist only of S-wave diquarkantidiquark configuration, as illustrated in Appendix.B, and their spatial wave function is symmetric under the particle exchange c 1 ↔ c 2 or c3 ↔ c4 .The higher doublet states T 4c,0(−) (7188), T 4c,0(−) (7242) have larger rms radii than the lower doublet states T 4c,0(−) (6890), T 4c,0(−) (6945) and can be viewed as the radial excitation of the latter.
Similar to the 0 +− system, S-wave diquarkonium threshold does not exist in the 2 +− cccc system, and all of the states in the system are identified as zero-width states.Due to the restriction of antisymmetrization of identical particles, the only allowed color configuration for these states is χ3 c ⊗3c .The state T 4c,2(−) (7252) is the radial excitation of the ground state T 4c,2(−) (6950).
It should be noted that although S-wave diquarkonium threshold does not exist in the 0 +− and 2 +− cccc systems, diquarkonium thresholds with higher orbital angular momentum do exist in these systems.For example, the P-wave J/ψχ c1 state and the S-wave h c χ c1 state can form the 0 +− or 2 +− cccc system, while the S-wave η c ψ 2 state and the Pwave J/ψχ c2 can form the 2 +− cccc system.These scattering states have the same quantum numbers as the zero-width states obtained in the present calculations.The coupling between them may alter the positions of the zero-width states.Considering the effect of this coupling is beyond the scope of this work.However, if one assumes the coupling effect is small and treats it as perturbation, the positions of the states should not change by much.The zero-width states may obtain nonzero widths and transform into resonant states, which can decay into the diquarkonium channels with lower energies.All of these states lie above the ground state diquarkonium threshold J/ψχ c1 , whose theoretical energy in the AP1 potential is 6593 MeV.Therefore, they may be searched for in the P-wave J/ψχ c1 decay channel in the experiment.

B. Fully bottomed tetraquark
In the cccc system, we observe that the results in different quark potential models qualitatively agree with each other.Therefore, we only choose the AP1 potential to investigate the bb bb system.The complex eigenenergies of the bb bb states with normal and exotic C-parities are shown in Fig. 4.
Similar to the cccc system, we obtain a series of resonant states with normal C-parity in the bb bb system.The complex energies, the proportions of different color configurations and the rms radii of these resonant states are listed in Table VII.The different rms radii of these states are around 0.5 fm, indicating that all of these states have compact tetraquark configuration.We obtain a lower resonant state with mass M ≈ 19780 MeV and width Γ ≈ 60 MeV, and a higher resonant state with mass M ≈ 19950 MeV and width Γ ≈ 40 MeV in the 0 ++ , 1 +− and 2 ++ systems.We also obtain two narrow resonant states T 4b,0(+) (19898) and T 4b,1(−) (19762).These states can decay strongly and be searched for in the corresponding diquarkonium decay channels in the experiment.
For the bb bb system with exotic C-parity, we obtain a series of resonant and zero-width states, whose complex energies, proportions of different color configurations and rms radii are listed in Table VIII.For the 0 +− system, we obtain two zero-width states T 4b,0(−) (19722) and T 4b,0(−) (19754).The lower one is dominated by the χ3 c ⊗3c color configura- tion while the higher one is dominated by the χ 6c⊗ 6c color configuration.The mass hierarchy of the χ3 c ⊗3c dominated state and the χ 6c⊗ 6c dominated state is reversed compared to the cccc system.This suggests that in the bb bb system, the interaction within the (anti)diquark plays a more important role than the interaction between diquark and antidiquark.For the 1 ++ system, we obtain two extremely narrow resonant states T 4b,1(+) (19731) and T 4b,1(+) (19748), whose two-body decay widths are less than 0.1 MeV.For the 2 ++ system, a zerowidth state T 4b,2(−) (19741) is obtained.These states may couple with the P-wave diquarkonium thresholds.Due to the coupling effect, the 0 +− and 2 +− zero-width states may transform into resonant states, which can decay into the P-wave Υχ b1 channel.

C. Resonances with small widths
In our calculations, we obtain several extremely narrow resonant states with quantum numbers J P = 1 + in both the cccc and bb bb systems.Intuitively, one would assume that these states can decay into the diquarkonium channels with the same quantum numbers and lower energies.However, the narrow widths of these states indicate that the two-body decay process is suppressed.The reasons for the suppression are given as follows.
The two-body decay width of the tetraquark state is proportional to the modulus square of the T -matrix element [97], where F denotes the tetraquark state and A, B denote the final mesons.The potential V ij is given in Eq. ( 2), comprising the color magnetic term and spin-independent terms.In the fully heavy system, the color magnetic transition is suppressed by the heavy quark mass.The spin-independent terms contain the color factor λ i • λ j , whose matrix elements are listed in Appendix A. It can be seen that 4 i<j=1 λ i • λ j is proportional to the identity operator and cannot induce transition between different color configurations.Color mixing via spin-independent terms can occur only when the coefficients of λ i • λ j are different.Roughly speaking, the magnitudes of color mixing matrix elements depend on the differences between various coefficients.From Tables IV,VI,VII,VIII, we can see that the different rms radii of the narrow resonances with J P = 1 + are of the same order and the differences between them are rather small.Therefore, the color mixing matrix elements are suppressed in the system.

IV. SUMMARY
In summary, we calculate the mass spectrum of the S-wave fully heavy tetraquark systems with both normal (J P C = 0 ++ , 1 +− , 2 ++ ) and exotic (J P C = 0 +− , 1 ++ , 2 +− ) Cparities using three different quark potential models (AL1, AP1, BGS).The exotic C-parity systems refer to the ones that have no corresponding S-wave ground heavy quarkonia thresholds.We employ the Gaussian expansion method to solve the four-body Schrödinger equation, and the complex scaling method to distinguish resonant states from scattering states.
Our calculations show that the mass spectra in different quark models are in qualitative agreement.We obtain a series of resonant states with J P C = 0 ++ , 1 +− , 2 ++ and 1 ++ .Moreover, we obtain several zero-width states in the 0 +− and 2 +− systems, where S-wave diquarkonium threshold does not exist.For the fully charmed system, we compare the theoretical results in the AP1 and BGS potentials with the experimental results in Fig. 5.We do not display the results in the AL1 potential since they are nearly the same as those in the AP1 potential.We find good candidates for the experimentally observed X(6900) and X(7200) in both the 0 ++ and 2 ++ systems.However, signals for the X(6400) and X(6600) are not seen in our calculations.Several resonant and zero-width fully charmed tetraquark states await to be found in the experiment.For the fully bottomed system, we summarize the theoretical results in the AP1 potential in Fig. 6.Resonant and zero-width fully bottomed tetraquark states are predicted   in the mass region (19.7,20.0) GeV.By investigating the root mean square radii of the states, we find that most of the resonant and zero-width states have compact tetraquark configuration, except that the X(7200) candidates T 4c,0(+) (7173), T 4c,2(+) (7214) and the state T 4c,1(−) (7191) may have molecular configurations.We also study the decay modes of the resonant and zero-width states.The resonant states can decay strongly to S-wave diquarkonium thresholds while the zero-width states can only decay to P-wave quarkonia.Further study that considers the P-wave tetraquark systems is needed to better establish the properties of the 0 +− and 2 +− states. • •   n1,n2,n3 is the S-wave Gaussian spatial wave function.The detail expression for the wave functions can be found in Eqs. ( 6) and (8).
We can see that the S-wave fully heavy tetraquark states with J P C = 0 +− and 2 +− contain only diquark-antidiquark configuration Φ n1,n2,n3 , which arises from the fact that S-wave dimeson systems cannot have such quantum numbers.It should be noted that the minus sign in Eqs.(B7), (B8) and (B20) demands that the radial excitation of the diquark and the antidiquark must be different, therefore the 0 +− or 2 +− tetraquark state is not a particle-antiparticle pair and one cannot simply calculate the C-parity as C = (−1) L+S = +1.

Appendix C: TWO DEFINITIONS OF ROOT MEAN SQUARE RADIUS
In our calculations, we only use the non-antisymmetric component of the wave function to calculate the rms radii.It seems more reasonable to calculate the rms radii of the compact tetraquark states using the complete wave function.However, our primary interest lies in the general clustering behavior of the tetraquark states rather than in specific numerical results of the rms radii, since the latter are not experimentally observables at present.The rms radii calculated using the non-antisymmetric component of the wave function are already capable of distinguishing between compact tetraquark states and loose molecular states.To illustrate this, we compare the results of rms radii calculated using the complete wave function |Ψ J (θ)⟩ and the non-antisymmetric term |Ψ J nA (θ)⟩ in Table X.One compact and one molecular cccc resonant states with J P C = 0 ++ are chosen as examples.For the state T 4c,0(+) (6978), both the results from |Ψ J (θ)⟩ and |Ψ J nA (θ)⟩ indicate that all rms radii are of the same order and support the compact tetraquark configuration.For the state T 4c,0(+) (7173), the results from |Ψ J nA (θ)⟩ can clearly demon-strate the clustering behaviour of a molecular state: (c 1 c3 ) and (c 2 c4 ) form two subclusters, which are widely separated.On the other hand, the results from |Ψ J (θ)⟩ are much less clear due to the antisymmetrization.It should also be noted that all r rms ci cj are the same in the results from |Ψ J (θ)⟩ due to the antisymmetrization of the wave function.
In conclusion, calculating the rms radii using only the nonantisymmetric term can reflect the internal spatial structure of tetraquark states more transparently.TABLE X.The rms radii (in fm) of the cccc resonant states calculated using the complete wave functions |Ψ J (θ)⟩ in Eq. ( 20) and the nonantisymmetric components of the wave functions |Ψ J nA (θ)⟩ in Eq. ( 21), respectively.The last column shows the spatial configurations of the states, where C. and M. represent the compact tetraquark and molecular configurations, respectively.

FIG. 2 .
FIG. 2. The complex energy eigenvalues of the cccc states with normal C-parity in the (a) AL1, (b) AP1 and (c) BGS potential with varying θ in the CSM.The solid lines represent the continuum lines rotating along Arg(E) = −2θ.The resonances do not shift as θ changes and are highlighted by the orange circles.

FIG. 3 .χ3 c ⊗3c χ 6c⊗ 6c r rms b 1 b3 r rms b 2 b4 r rms b 1 b4 = r rms b 2 b3 r rms b 1 b 2
FIG. 3. The complex energy eigenvalues of the cccc with exotic C-parity in the (a) AL1, (b) AP1 and (c) BGS potential with varying θ in the CSM.The solid lines represent the continuum lines rotating along Arg(E) = −2θ.The zero-width states and resonances do not shift as θ changes and are highlighted by the orange circles.

FIG. 4 .
FIG. 4. The energy eigenvalues of the bb bb states with (a) normal and (b) exotic C-parities in the AP1 potential with varying θ in the CSM.The solid lines represent the continuum lines rotating along Arg(E) = −2θ.The zero-width states and resonances do not shift as θ changes and are highlighted by the orange circles.

FIG. 5 .
FIG.5.The mass spectrum of the S-wave cccc states in the (a) AP1, (b) BGS potential.The experimental results reported by LHCb (model I)[13], CMS (noninterference model)[14] and ATLAS (model A and α)[15] are shown in Fig. (c).In Fig.(a) and (b), the green and blue lines represent the theoretical cccc states with widths larger and smaller than 1 MeV, respectively.In Fig. (c), the red lines and the dashed areas represent the central masses and the uncertainties in the experiments, respectively.The dotted lines represent dicharmonium thresholds, whose experimental energies are taken from Ref. [90].

FIG. 6 .
FIG.6.The mass spectrum of the S-wave bb bb states in the AP1 potential.The green and blue lines represent the states with widths larger and smaller than 1 MeV, respectively.The dotted lines represent dicharmonium thresholds, whose experimental energies are taken from Ref.[90].

TABLE I .
The parameters in the AL1, AP1 and BGS quark potential models.

TABLE IV .
The proportions of different color configurations and the rms radii (in fm) of the cccc resonant states with normal C-parity in the AP1 potential.The last column shows the spatial configurations of the states, where C. and M. represent the compact tetraquark and molecular configurations, respectively.

TABLE V .
The complex energies E = M − iΓ/2 (in MeV) of the cccc resonant and zero-width states with exotic C-parity from various potential models.

TABLE VI .
The proportions of different color configurations and the rms radii (in fm) of the cccc resonant and zero-width states with exotic C-parity in the AP1 potential.The last column shows the spatial configurations of the states, where C. and M. represent the compact tetraquark and molecular configurations, respectively.

TABLE VII .
The proportions of different color configurations and the rms radii (in fm) of the bb bb resonant states with normal C-parity in the AP1 potential.The last column shows the spatial configurations of the states, where C. and M. represent the compact tetraquark and molecular configurations, respectively.