Most Strange Dibaryon from Lattice QCD

The $\Omega\Omega$ system in the $^1S_0$ channel (the most strange dibaryon) is studied on the basis of the (2+1)-flavor lattice QCD simulations with a large volume (8.1 fm)$^3$ and nearly physical pion mass $m_{\pi}\simeq 146$ MeV at a lattice spacing $a\simeq 0.0846$ fm. We show that lattice QCD data analysis by the HAL QCD method leads to the scattering length $a_0 = 4.6 (6)(^{+1.2}_{-0.5}) {\rm fm}$, the effective range $r_{\rm eff} = 1.27 (3)(^{+0.06}_{-0.03}) {\rm fm}$ and the binding energy $B_{\Omega \Omega} = 1.6 (6) (^{+0.7}_{-0.6}) {\rm MeV}$. These results indicate that the $\Omega\Omega$ system has an overall attraction and is located near the unitary regime. Such a system can be best searched experimentally by the pair-momentum correlation in relativistic heavy-ion collisions.

Introduction: Dibaryon is defined as a baryon number B = 2 system (equivalently a 6-quark system) in quantum chromodynamics [1][2][3]. So far, only one stable dibaryon, the deuteron, has been observed: It is a loosely bound system of the proton and the neutron in spin-triplet and isospin-singlet channel. In recent years, there are renewed experimental interests in the dibaryons due to exclusive measurements in hadron reactions [4] as well as the direct measurement in relativistic heavyion collisions [5]. Also from the theoretical side, (2+1)flavor lattice QCD simulations of the 6-quark system in a large box with nearly physical quark masses became possible recently [6]. The main aim of this Letter is to report the first result and physical implication of ΩΩ, the strangeness S = −6 dibaryon ("most strange dibaryon"), in full QCD simulations with the lattice volume (8.1 fm) 3 and the pion mass m π 146 MeV at a lattice spacing a 0.0846 fm [7].
Before entering the detailed discussions of our study, we first introduce the reason why such an exotic channel (S = −6) is of interest in QCD. Let us consider octet 8 and decuplet 10 baryons in the flavor SU(3) classification. All the members of 8 are stable under strong decay. This is why the forces between octet baryons in 8 ⊗ 8 are most relevant in the physics of hypernuclei and of neutron stars. Also, the elusive H-dibaryon (a combination of ΛΛ, N Ξ and ΣΣ) is in this representation [8][9][10] and does not suffer from the Pauli exclusion principle in the flavor-SU(3) limit.
On the other hand, only Ω in 10 is stable under strong decay. Therefore, in the 8 ⊗ 10 representation, the most promising candidate of stable dibaryon is N Ω [11]. The Pauli exclusion principle does not work in this case too, so that there is a possibility to have a bound state in the S-wave and total-spin 2 channel [12]. Such a system is indeed studied by the two-particle momentum correlation in high-energy heavy-ion collisions both theoretically and experimentally [13].
In the decuplet-decuplet channnel, we have where "sym." and "anti-sym." stand for the flavor symmetry under the exchange of two baryons. Only possible stable state under strong decay is the ΩΩ system in the symmetric 28 representation. Again, the quark Pauli principle does not operate in this channel [14]. Note that the celebrated ABC resonance (∆∆ in the spin-3 and isospin-0 channel) [4,15] belongs to the anti-symmetric 10 * representation, while ∆∆ in the spin-0 and isospin-3 channel is in the same multiplet with ΩΩ.
The ΩΩ interaction at low energies has been investigated so far by using phenomenological quark models or by using lattice QCD simulations with heavy quark masses. Very recently, the chiral effective field theory has also been applied to the scattering of the Ω baryons [16]. In some quark models, strong attraction is reported [17,18], while other models show weak repulsion [19,20]. A (2+1)-flavor lattice QCD study with m π 390 MeV by using the finite volume method shows weak repulsion [21], while a study with m π 700 MeV by using the HAL QCD method shows moderate attraction [22]. Under such a controversial situation, it is most important to carry out first-principles lattice QCD simulations in a large volume with the pion mass close to the physical point. HAL QCD method: In the HAL QCD method [23][24][25][26], the observables such as the binding energy and phase shifts are obtained from the equal-time Nambu-Bethe-Salpeter (NBS) wave function ψ(r) and associated twobaryon irreducible kernel U (r, r ). The traditional finite volume method with the plateau fitting [27] is challenging in practice for B = 2 systems in large volumes because of the difficulty in differentiating each scattering state [28][29][30][31]. On the other hand, the time-dependent HAL QCD method [25] can avoid such a problem since all the elastic scattering states are dictated by the same kernel U (r, r ) arXiv:1709.00654v2 [hep-lat] 24 May 2018 and there is no need to identify each scattering state in a finite box.
The equal-time NBS wave function ψ(r) has the property that its asymptotic behavior at large distance reproduces the phase shifts, which can be proven from the unitarity of S-matrix in quantum field theories [23]. Moreover, it is related to the following reduced four-point (4-pt) function, R(r, t > 0) = 0|Ω(r, t)Ω(0, t)J (0)|0 /e −2m Ω t = n a n ψ n (r)e −(δWn)t + O(e −(∆E * )t ). (1) Here J (0) is a source operator creating the (B, S) = (2, −6) system at Euclidean time 0, and a n is the matrix element defined by n|J (0)|0 with |n representing the elastic states in a finite volume. The energy is represented as δW n = 2 m 2 Ω + k 2 n − 2m Ω with the Ω baryon mass m Ω and the relative momentum k n . Typical excitation energy of a single Ω-baryon is denoted by ∆E * , so that the last term in Eq.(1) is exponentially suppressed as long as t (∆E * ) −1 ∼ Λ −1 QCD [25] with Λ QCD ∼ 200 − 300 MeV being the QCD scale parameter. A local interpolating operator for the Ω baryon has a general form with a, b and c being color indices, γ k being the Dirac matrix, α being the spinor index, and C ≡ γ 4 γ 2 being charge conjugation. An appropriate spin projection is necessary from this operator to single out a particular spin state as mentioned later. The reduced 4-pt function R has been shown to satisfy the following master equation [25], which is valid as long as t (∆E * ) −1 . We emphasize that we do not need to isolate each scattering state with the energy δW n , so that only the moderate values of t up to 1.5-2 fm are sufficient for a reliable extraction of the kernel U . This is a crucial difference from the finite volume method which requires t (δW 1 ) −1 ∼ 10 fm (for the lattice volume as large as 8 fm) to identify each δW n . (See a recent summary [31] and references therein.) At low energies, one can make the derivative expansion with respect to the non-locality of the kernel [24,32], Then, we introduce an "effective" leading-order potential V (r), which provides a good leading-order approximation of U (r, r ) to obtain physical observables at low energies, as long as its t-dependence is sufficiently small.
Interpolating operator: The interpolating operator Ω sz (x) for the Ω baryon with spin- 3 2 and the z-component s z = ± 3 2 , ± 1 2 can be readily constructed by the appropriate spin projection of the upper two components of Eq.(2) as shown in [22]. The asymptotic ΩΩ system can now be characterized by 2s+1 L J with the total spin (s), the orbital angular momentum (L) and the total angular momentum (J). We consider L = 0 where maximum attraction is expected at low energies. Then, the Fermi statistics leads s to be even (either s=0 or 2). Here we consider an s = 0 system with the interpolating operator For J (0), we use the wall-type quark source with the s = 0 projection given above. With this source, the total momentum of the system is automatically zero. Also, it has good overlap with the scattering states where |r| in Eq. (3) is larger than the typical baryon size. To extract the L = 0 and s = 0 state at t on the lattice, we employ Eq.(5) for the sink operator together with the projection to the A 1 representation of the cubic-group, where R i is an element of the cubic group acting on the relative distance r. Note here that R(r, t) and U (r, r ) depend on the choice of interpolating operators, while observables calculated from these quantities are independent of the choice thanks to the Nishijima-Zimmermann-Haag theorem [24]. Lattice Setup: By using the 11PFlops supercomputer K at RIKEN Center for Computational Science, (2 + 1)flavor gauge configurations on the 96 4 lattice are generated with the Iwasaki gauge action at β = 1.82 and nonperturbatively O(a)-improved Wilson quark action with stout smearing. The lattice spacing is a 0.0846 fm (a −1 2.333 GeV) [7] and the pion mass, the kaon mass and the nucleon masses are m π 146 MeV, m K 525 MeV and m N 964 MeV, respectively. (These masses are higher than the physical values by about 8%, 6 % and 3 %, respectively, due to slightly larger quark masses at the simulation point.) The lattice size, La 8.1 fm, is sufficiently large to accommodate two baryons in a box.
We employ the wall quark source with the Coulomb gauge fixing, and the periodic (Dirichlet) boundary condition is used for spatial (temporal) directions. Forward and backward propagations are averaged to reduce the statistical fluctuations. We pick 1 configuration per each 5 trajectories, and make use of the rotation symmetry and the translational invariance for the source position to increase the statistics. The total statistics in this Letter amounts to 400 configurations × 2 (forward/backward) [ [34][35][36][37] and the correlation functions are calculated using the unified contraction algorithm [38].
The Ω-baryon mass extracted from the effective mass m eff (t) ≡ ln G(t)/G(t + a) with G(t) being the baryonic two-point function is m Ω = 1712 ± 1 MeV (from the plateau in t/a = 17 − 22) and m Ω = 1713 ± 1 MeV (from t/a = 18 − 25) with the statistical errors. These numbers are about 2% higher than the physical value, 1672 MeV. We take the former number in the following analysis. Numerical results: The 1 S 0 potential V (r) obtained from Eq.(4) with the lattice measurement of R(r, t) is shown in Fig. 1 for t/a = 16, 17, and 18. Here Laplacian and the time-derivative in Eq.(4) are approximated by the central (symmetric) difference. The statistical errors for V (r) at each r are estimated by the jackknife method with the bin size of 40 configurations. A comparison with the bin size of 20 configurations shows that the bin size dependence is small. The particular region t/a = 17 ± 1 in Fig. 1 is chosen to suppress contamination from excited states in the single Ω propagator at smaller t and simultaneously to avoid large statistical errors at larger t. We observe that the potentials at t/a = 16, 17, and 18 are nearly identical within statistical errors as expected from the time-dependent HAL QCD method [25].
The ΩΩ potential V (r) has qualitative features similar to the central potential of the nucleon-nucleon (N N ) interaction, i.e., the short range repulsion and the intermediate range attraction [6]. There are, however, two quantitative differences: (i) the short range repulsion is much weaker in the ΩΩ case possibly due to the absence of quark Pauli exclusion effect, and (ii) the attractive part is much short-ranged due to the absence of pion exchanges.
For the purpose of converting the potential to the physical observables such as the scattering phase shifts and the binding energy, we fit V (r) in Fig. 1 in the range r = 0 − 6 fm by three Gaussians, V fit (r) = The finite volume effect on the potential is expected to be small due to the large lattice size. The naive order estimate of the finite a effect for the physical observables is also small ((Λa) 2 ≤ 1 %) thanks to the non-perturbative O(a) improvement, but an explicit confirmation would be desirable in the future. The ΩΩ scattering phase shifts δ(k) in the 1 S 0 channel obtained from V fit (r) are shown in Fig.2 for t/a = 16, 17, and 18 as a function of the kinetic energy in the center of mass frame, E CM = 2 k 2 + m 2 Ω − 2m Ω . The error bands reflect the statistical uncertainty of the potential in Fig.  1. All three cases show that δ(0) starts from 180 • , which indicates the existence of a bound ΩΩ system.
The scattering length a 0 and the effective range r eff in the 1 S 0 channel is extracted from δ(k) through the effective range expansion (ERE), k cot δ(k) = − 1 a0 + 1 2 r eff k 2 + · · · , with the sign convention of nuclear and atomic physics: The central values and the statistical errors in the first parentheses are extracted from δ(k) at t/a = 17, and the systematic errors in the second parentheses are estimated from the results at t/a = 16 and 18. The origin of this systematic error is the truncation of the higher-derivatives of the non-local potential as well as the contaminations from inelastic states. To get a feel for the magnitude of these values, we recapitulate here the experimental values of a 0 and r eff in the N N systems; (a 0 , r eff ) spin-triplet = (5.4112 (15) FIG. 3. The dimensionless ratio of the effective range r eff and the scattering length a0 as a function of r eff for the ΩΩ system in the 1 S0 channel as well as for the spin-triplet N N system (the deuteron channel) and for the spin-singlet N N system (the neutron-neutron channel). The error bar for ΩΩ is the quadrature of the statistical and systematic errors in Eqs. (7) and (8). ing parameters in the N N systems and those in the ΩΩ system should be similar. Nevertheless, it is remarkable that they are all close to the unitary region where r eff /a 0 is substantially smaller than 1 as shown in Fig.3. Shown in Fig.4 are the bound state energy given by the opposite sign of the binding energy, −B ΩΩ , and the rootmean-square distance ( r 2 ) of the ΩΩ bound state obtained from the potential. The blue diamond is taken from the data at t/a = 17 without the Coulomb repulsion. The blue solid and dashed lines are the statistical error at t/a = 17 and the systematic error estimated from the data at t/a = 17 ± 1, respectively: B (QCD) ΩΩ = 1.6(6)( +0.7 −0.6 ) MeV. (9) As an alternative estimate, the truncation error of the derivative expansion on the binding energy is evaluated perturbatively, and is found to be less than 20% even if the magnitude of dimensionless next-to-leading order potential is the same order as that of the effective leadingorder potential. It is an important future subject to determine higher-order potentials explicitly by using the method of multiple quark sources [33]. The binding energy is consistent with the value obtained from the general formula for loosely bound states [40] with (7) and (8); Associated with this small binding energy, r 2 is as large as 3-4 fm which is consistent with the expectation, r 2 ∼ a 0 , for loosely bound states. The Coulomb repulsion can be evaluated by adding α/r with α = e 2 /4π to the potential obtained from lattice QCD, i.e. V (QCD+Coulomb) ≡ V (QCD) + α/r. This reduces the above binding energy by a factor of two, B It is in order here to remark that there are three energy scales in the present problem: 2m Ω 3400 MeV |V (r 0.5fm)| ∼ 50 MeV B ΩΩ ∼ 1 MeV. Since only the relative difference between the interacting and noninteracting two-Ω systems matters, the absolute magnitude of the uncertainty of 2m Ω is not reflected directly in V (r). This is why we could extract V (r) rather accurately as shown in Fig.1 despite of the large scale difference between 2m Ω and V (r). Then the small binding energy B as well as the large scattering length a 0 are the natural consequence of the large cancellation between the long-range attraction and the short-range repulsion of V (r), a situation common in nuclear and atomic physics. Although V (r) is not a direct observable, it provides an important intermediate step to link the QCD scale (GeV) to nuclear physics scale (MeV), since it is difficult to measure the binding energy directly from lattice QCD using the finite volume method for large lattice volumes and physical quark masses (see the critical review [31]).
Finally let us estimate the effect of slightly heavy quark masses in our simulation. First of all, the attractive part of the ΩΩ potential would be slightly longranged at the physical point due to the kaon mass, m Summary and Discussions: In this Letter, we presented a first realistic calculation on the most strange dibaryon, ΩΩ, in the 1 S 0 channel on the basis of the (2+1)-flavor lattice QCD simulations with a large volume and nearly physical quark masses. The scattering length, effective range and the binding energy obtained by the HAL QCD method strongly indicate that the ΩΩ system in the 1 S 0 channel has an overall attraction and is located in the vicinity of the unitary regime. From the phenomenological point of view, such a system can be best searched by the measurement of pair-momentum correlation C(Q) with Q being the relative momentum between two baryons produced in relativistic heavy-ion collisions [5]. Experimentally, each Ω can be identified through a successive weak decay, Ω − → Λ + K − → p + π − + K − . Note that a large scattering length (not the existence of a bound state) is the important element for C(Q) to have characteristic enhancement at small relative momentum Q. Moreover, the effect of the Coulomb interaction can be effectively eliminated by taking a ratio of C(Q) between small and large collision systems [13].
We are currently underway to increase the statistics of the lattice data with the same lattice setup. These results together with the detailed examination of the spectrum analysis in a finite lattice volume and the effective range expansion will be reported.