Hyperon-nucleon force from lattice QCD

We calculate potentials between a proton and a $\Xi^0$ (hyperon with strangeness -2) through the equal-time Bethe-Salpeter wave function, employing quenched lattice QCD simulations with the plaquette gauge action and the Wilson quark action on (4.5 fm)^4 lattice at the lattice spacing $a \simeq 0.14$ fm. The ud quark mass in our study corresponds to $m_{\pi}\simeq 0.37$ and 0.51 GeV, while the s quark mass corresponds to the physical value of $m_K$. The central $p \Xi^0$ potential has a strong (weak) repulsive core in the $^1S_0$ ($^3S_1$) channel for $r \lsim 0.6$ fm, while the potential has attractive well at the medium and long distances (0.6 fm $\lsim r \lsim 1.2$ fm) in both channels. The sign of the $p \Xi^0$ scattering length and its quark mass dependence indicate a net attraction in both channels at low energies.


Introduction
Modern nucleon-nucleon (N N ) potentials [1] give successful description of the N N scattering data at low energies and have been used for precision calculations in light nuclei [2]. Similarly, hyperon-nucleon (Y N ) and hyperon-hyperon (Y Y ) potentials in the coordinate space are known to be quite useful and important for studying the properties of hypernuclei [3] as well as the equation of state for hyperonic matter in neutron stars [4]. However, there are still large experimental uncertainties in Y N and Y Y interactions at present, because the short life-time of the hyperons makes the scattering experiments difficult. Accordingly, phenomenological Y N and Y Y potentials are not well constrained from data even under some theoretical guides [5,6,7,8,9,10].
In such a situation, it would be desirable to analyse hyperon interactions on the basis of the first principle lattice QCD simulations. Studies along this line for nucleon interactions was initiated in Ref. [11], where the N N scattering lengths were extracted from the quenched simulations by using the Lüscher's finite volume method [12]. The same method was applied later to the (2+1)-flavor simulations with the mixed action [13]. The scattering length of ΛN system was first examined in [14] in quenched simulations on a small lattice box. Subsequently, ΛN and ΣN scattering lengths were studied in (2+1)-flavor simulations with the mixed action [15]. Such low-energy scattering parameters calculated in lattice QCD would be valuable inputs to construct phenomenological Y N and Y Y potentials to be used for studying hypernuclei and hyperonic matter.
Recently, an alternative but closely related approach to [12] has been proposed to define the N N potential from lattice QCD [16,17,18]. Since the potential is not a direct physical observable, one cannot make quantitative comparison of this lattice N N potential with the phenomenological N N potentials. However, they are both designed to reproduce the correct scattering phase shifts, and their spatial structures are found to have common features [16]; the attraction at long and intermediate distances and the strong repulsion at short distance in the S-wave channel [19].
The purpose of the present paper is to report our first attempt to apply the above approach to the hyperonnucleon systems. For the N N potential, only two channels (isovector and isoscalar) exist in the flavor SU(2) space, i.e., 2 ⊗ 2 = 3 ⊕ 1. Including the strange quark extends this to 8 ⊗ 8 = 27 ⊕ 10 * ⊕ 1 ⊕ 8 ⊕ 10 ⊕ 8 in the flavor SU(3) space. The isovector (isoscalar) channel in the N N sector is assigned to be a subset in the 27-plet (10 * -plet) representation. The potentials in newly arising channels are hardly determined from experiments so far.
In this paper, we focus on the N Ξ potential in the isovector (I = 1) channel as a first step. (A preliminary account of this system has been reported in [20].) There are two main reasons for picking up this channel: (i) Theoretically, it is the simplest generalization of the N N system; pΞ 0 is obtained from pn by replacing the d-quarks in the neutron by the s-quarks. Also it is a channel which do not have strong decay into other Y N systems. 2 (ii) Experimentally, not much information has been available on the N Ξ interaction except for a few studies; a recent report gives the upper limit of elastic and inelastic cross sections [21], and earlier publications suggest weakly attractive Ξ-nucleus potential [22]. The Ξ-nucleus interaction will be soon studied as one of the day-one experiments at J-PARC [23] via (K − , K + ) reaction with nuclear target. This paper is organized as follows. In Section 2, we describe the basic formulation to derive the pΞ 0 potential through the Bethe-Salpeter amplitude measured in the lattice QCD simulations. Our lattice setup is explained in Section 3. In Section 4, we show numerical results of the potentials in the spin-singlet and spin-triplet channels and their quark mass dependence. We also show the estimate of the scattering lengths in these channels. Section 5 is devoted to summary and concluding remarks.

Basic formulation
Our methodology to obtain the baryon potentials is along the lines of [16,17,18]. (See also [12,24] for the seminal attempts to introduce similar notion of the potential.) We consider the low-energy N -Ξ scattering and start with an effective Schrödinger equation for the equal-time Bethe-Salpeter (BS) wave function φ( r) obtained from the Lippmann-Schwinger equation [18,25]; Here µ = m N m Ξ /(m N + m Ξ ) and E ≡ k 2 /(2µ) are the reduced mass of the N Ξ system and the nonrelativistic energy in the center-of-mass frame, respectively. Eq. (1) is derived from QCD by adopting local three-quark operator as the nucleon interpolating operator. It can be shown that the non-local potential U ( r, r ′ ) is energy independent. Also, it is unique given the interpolating operator as long as the BS-amplitude for all energies is known [18]. If we take other nucleon interpolating operator, both the BS wave function and the non-local potential change their spatial structure without affecting the observables such as the binding energy and the phase shift. At low energies, the nonlocality of the potential can be expanded as [26] U ( r, where the dimensionless expansion parameter is ∇/M with M being the typical scale of the strong interaction such as the excitation energy of a single baryon. The general expression of the potential V N Ξ is known to be [27] Here S 12 = 3( σ 1 · n)( σ 2 · n) − σ 1 · σ 2 is the tensor operator with n = r/| r|, S ± = ( σ 1 ± σ 2 )/2 are symmetric (+) and antisymmetric (−) spin operators, L = −i r × ∇ is the orbital angular momentum operator, and τ N ( τ Ξ ) is isospin operator for N = (p, n) T (Ξ = (Ξ 0 , Ξ − ) T ). We note that the antisymmetric spin-orbit forces (V ALS and V ALSτ ) do not arise in the N N case because of the identical nature of the nucleon within the isospin symmetry. According to the above expansion, the wave function should be classified by the total isospin I, the total angular momentum J = L + S + and the parity. A particular spin (isospin) projection is made in terms of σ N · σ Ξ ( τ N · τ Ξ ), e.g., for the isospin projection we have P (I=0) = (1 − τ N · τ Ξ )/4 and P (I=1) = (3 + τ N · τ Ξ )/4. The equal-time BS wave function for (I, I z ) = (1, 1) and L = 0 (S-wave) on the lattice is obtained by where α and β denote the Dirac indices, a, b and c the color indices, and C = γ 4 γ 2 the charge conjugation matrix. The summation over R ∈ O is taken for cubic transformation group to project out the S-wave state. 3 The summation over x is to select the state with zero total momentum. Here we take local field operator p α (x) and Ξ 0 β (y) for the proton and Ξ 0 . The wave function and the potential (or equivalently the off-shell behavior of the scattering amplitude) depend on the choice of interpolating operators. This is the situation common to any field theories. In this paper, we focus exclusively on the local operators as introduced above and leave further discussions on the operator dependence to [25]. We take the upper components of the Dirac indices to construct the spin singlet (triplet) channel by P The BS wave function φ( r) is understood as a probability amplitude to find "nucleon-like" three-quarks located at point x + r and "Ξ-like" three-quarks located at point x. Our BS wave function has information not only of the elastic amplitude N Ξ → N Ξ but also of the inelastic amplitudes: The inelastic effects are localized in coordinate space at low energies and do not affect the asymptotic behavior of φ( r) and hence the phase shift.
In the actual simulations, the BS wave function is obtained from the four-point correlator, Here . The eigen-energy and the eigen-state of the six quark system are denoted by E n and |E n , respectively, with the matrix element A n (t 0 ) = E n |J pΞ 0 (t 0 )|0 . For t − t 0 ≫ 1, the F pΞ 0 and hence the wave function φ are dominated by the lowest energy state.
The lowest energy state |E 0 created by the wall source J pΞ 0 (t 0 ) contains not only the S-wave N Ξ component but also the components which can mix by QCD dynamics, such as the D-wave due to the tensor force and ΛΣ state due to the rearrangement of the quarks inside baryons. In principle, these components can be disentangled by preparing appropriate operator sets for the sink. Study along this line to extract the mixing between the S-wave and the D-wave in low energy N N interaction was put forward recently in [28]. In the present paper, instead of making such decomposition, we define an effective central potential V C (r) for the S-wave component according to our previous works [16,17,18]: Such a potential in the I = 1 sector has only the spin dependence as V where effects of the tensor force and coupled-channel effects are implicitly taken into account inṼ 0,σ .
It is in order here to make some remarks on the potential we have defined.
(i) The description by the energy-independent non-local potential U ( r, r ′ ) is equivalent to that by the momentum-dependent local potential, V ( r, ∇). The potential V C (r) in Eq. (9) is the leading order term of the derivative expansion of V ( r, ∇) and can be determined by the BS wave function at E ≃ 0. Higher order terms can be extracted successively by the BS wave functions with different E. It was recently reported that momentum-dependence of V ( r, ∇) is small at least up to the center of mass momentum p = 250 MeV for the N N system [29]. (ii) The non-local potential U ( r, r ′ ) has one-to-one correspondence to the baryon interpolating operator adopted in defining the Bethe-Salpeter amplitude. Different choices of the interpolating operator would give different BS wave function and the baryon-baryon potential, although the phase shifts and binding energies are unchanged. An advantage of working directly in lattice QCD is that we can unambiguously trace the relation between the baryon interpolating operator and the BS wave function or equivalently the associated potential.

Lattice setup and hadron masses
We use the plaquette gauge action and the Wilson fermion action with the gauge coupling β = 5.7 on the 32 3 ×32 lattice. The heatbath algorithm combined with the overrelaxation is used to generate quenched gauge configurations. After skipping 3000 sweeps for thermalization, measurement is made for every 200 sweep. The Dirichlet (periodic) boundary condition is imposed for quarks in the temporal (spatial) direction. The wall source is placed at t 0 = 5 with the Coulomb gauge fixing. These setup are same as our previous calculation for the N N potential [16].
Masses of pseudo-scalar and vector mesons obtained at hopping parameter (κ 1 , κ 2 ) are fitted as for κ i = 0.1678, 0.1665, 0.1640 (i = 1, 2). The fit gives κ c = 0.16930(1) for the critical hopping parameter. From m π = 135 MeV and m ρ = 770 MeV, we determine the hopping parameter for physical ud quarks and the lattice spacing as κ phys = 0.16910(1) and a = 0.1416(9) fm (1/a = 1.393(9) GeV), respectively, while the hopping parameter for the strange quark mass is given by κ s = 0.16432(6) from m K = 494 MeV. A corresponding lattice volume is (4.5fm) 4 , which is large enough to accommodate two baryons. In order to check the thresholds of two baryon systems with strangeness S = −2 (ΛΛ,N Ξ,ΛΣ and ΣΣ), we calculate the two-point correlator, C(t; t 0 ) = x 0|B α ( x, t)J Bα (t 0 )|0 , for the octet baryons (B = N, Ξ, Λ, Σ), and J Bα (t 0 ) is the wall-source for B. The interpolating fields for Λ and Σ + employed in this work are Σ + β (y) = −ε abc (u a (y)Cγ 5 s b (y)) u cβ (y).   from totally 1300 gauge configurations for the average, while no such exceptional configuration appears at κ ud = 0.1665 (m π ≃ 511 MeV). As seen in Table 1, the present results for the baryon masses are consistent with the experimentally observed ordering of the two-baryon thresholds in the strangeness −2 sector: E th (ΛΛ) < E th (N Ξ) < E th (ΛΣ) < E th (ΣΣ). This warrants that N Ξ in the I = 1 channel treated in this paper is indeed the lowest energy scattering state.

BS wave function
The left panel of Fig. 1 shows the wave functions obtained at time slice t − t 0 = 6 with t 0 = 5. The open circle (triangle) corresponds to the case of the 1 S 0 ( 3 S 1 ) channel. They are normalized to be unity at the spatial boundary r = (32/2, 0, 0). For r < ∼ 0.7 fm, the wave function is calculated for all possible values of r, while, for the outer region, it is calculated only on the x, y, z axis and their nearest neighbors to reduce the numerical cost.
As seen in the Figure, the wave functions are suppressed at short distance and are slightly enhanced at medium distance in both 1 S 0 and 3 S 1 channels. There is a sizable difference of the suppression pattern at short distance between 1 S 0 and 3 S 1 ; the repulsion is stronger in the 1 S 0 channel.

central potential
To derive the N Ξ potential through Eq. (9), we need to know the non-relativistic energy E = k 2 /(2µ). It has been shown in Ref. [24] that k 2 and thus E can be accurately determined by fitting the wave function in the asymptotic region in terms of the lattice Green's function  Table 1 Hadron masses in the unit of MeV calculated for the hopping parameters κ ud and κs, with the number of gauge configurations N conf .  Fig.1 (right). The results indicated by [6 − 10] are obtained by a constant fit to the data for t − t 0 = 6 − 10 with a single-elimination of the jackknife method.
which is the solution of (△ + k 2 )G( r, k 2 ) = −δ L ( r) with δ L ( r) being the periodic delta function [12,24]. Results of the fit in the range (12 ≤ x ≤ 16, 0 ≤ y ≤ 1, z = 0) are shown in the left panel of Fig. 1. The fitting range is determined so that it is outside the range of the interaction (see the right panel of Fig. 1). The non-relativistic energies thus obtained in the 1 S 0 ( 3 S 1 ) channel becomes E = −0.4(2) MeV (E = −0.8(2) MeV). We have checked that the results of the fit in different ranges 11 ≤ x ≤ 16 and 13 ≤ x ≤ 16 introduce systematic errors only less than half of the statistical errors for E. Note that E can be negative for the scattering state in a finite box if the interaction is attractive. The effective central potentials for pΞ 0 system obtained from Eq. (9) with the wave function φ and the energy E at t − t 0 = 6 are shown in the right panel of Fig. 1 for the 1 S 0 and 3 S 1 channels at m π ≃ 368 MeV. In order to check whether the ground state saturation of the pΞ 0 system is achieved in the present results, we plot the t-dependence of the potential at distances r = 0.14, 0.20, 0.71, 1.42 and 2.27 fm in Fig. 2. The fact that the potential is stable against t for t − t 0 ≥ 6 within errors indicates that the ground state saturation is indeed achieved at t − t 0 = 6. This is the reason why we adopted the values at t − t 0 = 6 in the right panel of Fig.1. Note that a constant fit to data for t − t 0 = 6 − 10 with a single elimination of the jackknife method does not introduce appreciable change of the final values and the errors of the potential as shown in Fig.2.
In the right panel of Fig. 1, the potential in the N Ξ system shows a repulsive core at r < ∼ 0.5 fm surrounded by an attractive well, similar to the N N system [16,17]. In contrast to the N N case, however, one finds that the repulsive core of the pΞ 0 potential in the 1 S 0 channel is substantially stronger than that in the 3 S 1 channel. Such a large spin dependence is also suggested by the quark cluster model [6]. The relatively weak attraction in the medium to long distance region that (0.6 fm < ∼ r < ∼ 1.2 fm) is similar in both 1 S 0 and 3 S 1 channels. As the energies determined above and scattering lengths given later indicate, the present potentials are weakly attractive on the whole in both spin channels in spite of the repulsive core at short distance. Also, the attraction 3 S 1 channel seems a little stronger than that in the 1 S 0 channel. the region that (10 ≤ x ≤ 16, 0 ≤ y ≤ 1, z = 0) is adopted to fit the wave function by the Green's function. We have checked that the results of the fit in different ranges 9 ≤ x ≤ 16 and 11 ≤ x ≤ 16 introduce systematic errors only less than half of the statistical errors for E.

quark mass dependence
The height of the repulsive core increases as the ud quark mass decreases, while the significant difference is not seen in the medium to long distances within the error bars.

one-pion exchange
The exchange of a single π 0 would conduct the interaction at long distance in the pΞ 0 system: The pseudo-vector πN N coupling f πN N and the πΞΞ coupling f πΞΞ are related as f πΞΞ = −f πN N (1 − 2α) with the parameter α = F/(F + D) [5]. Also we define g πN N ≡ f πN N mπ 2mN . The solid lines in Fig.3 is the one pion exchange potential (OPEP) obtained from Eq. (14) with m π ≃ 368 MeV, m N ≃ 1167 MeV (corresponding to κ ud = 0.1678) and the empirical values, α ≃ 0.36 [30] and g 2 πN N /(4π) ≃ 14.0 [1]. Unlike the N N potential in the S-wave, the OPEP in the present case has opposite sign between the spin-singlet channel and the spin-triplet channel. Also, the absolute magnitude of OPEP is weak due to the factor 1 − 2α. As is seen from Fig.3, we do not find clear signature of OPEP at long distance (r > 1.2 fm) in our potential within statistical errors. On the other hand, there is a clear departure from OPEP at medium distance (0.6fm < r < 1.2fm) in both 1 S 0 and 3 S 1 channels. These observations may indicate a mechanism of state-independent attraction such as the correlated two pion exchange.
It should be mentioned here that there is in principle a quenched artifact to the baryon-baryon potentials from the flavor singlet hairpin diagram (the ghost exchange) [33]. Its contribution to the central potential has a spin-dependent exponential tail V η C (r) ∝ σ N · σ Ξ exp(−m π r), which dominates over the Yukawa potential at large distances. Its significance can be estimated by comparing the sign and the magnitude of e mπr V C (r) in the spin-singlet and spin-triplet channels. Our present data at lightest quark mass m π = 368 MeV shows no evidence of the ghost contribution at large distances within errors. This may indicate the weak coupling of the η to N and Ξ.

scattering length
Using the energy E = k 2 /(2µ) obtained from the asymptotic behavior of the wave function in Sec.4.2, the N Ξ scattering lengths can be deduced from the Lüscher's formula [12,24], where Z 00 (1; q 2 ) with q = kL 2π is obtained by the analytic continuation of the generalized zeta-function Z 00 (s; q 2 ) = 1 √ 4π n∈Z 3 (n 2 − q 2 ) −s defined for Re s > 3/2. The sign of the S-wave scattering length a 0 is defined to be positive for weak attraction. The results are plotted in Fig.4 for two different values of the pion masses, m π ≃ 368 MeV and m π ≃ 511MeV. Both a 0t (the scattering length in the 3 S 1 channel) and a 0s (the scattering length in the 1 S 0 channel) slightly increase as the quark mass decreases: a 0s = 0.05(1) → 0.10(6) fm and a 0t = 0.09(1) → 0.21(6) fm. The positive sign indicates that the pΞ 0 interaction is attractive on the whole in both 1 S 0 and 3 S 1 channels. Due to the possible non-linear dependence of the scattering length as a function of the quark mass as is well-known in the case of the N N system [31] (see also [13,32] and references therein), we are not able to obtain a reliable estimate of the scattering length of the pΞ 0 system at the physical ud quark mass. In addition, a possible contamination from the ghost exchange may arise at small quark masses in the quenched QCD simulations [33]. Therefore, it is necessary in the future to carry out the full QCD calculation toward lighter quark masses to predict the precise value of a 0 .

Summary
We study the pΞ 0 interaction in the 1 S 0 and 3 S 1 channels through the equal-time Bethe-Salpeter amplitude measured by the quenched lattice QCD simulations on a (4.5 fm) 4 lattice with the quark masses corresponding to m π /m ρ ≃ 0.45 and 0.59. We adopt specific choice of the baryon interpolating operators as given in Eqs. (5) and (6) and extract the potential associated with these operators. The effective central potential deduced from the Bethe-Salpeter wave function has repulsive core at short distance surrounded by attractive well at the medium and long distances. The scattering lengths determined by fitting the asymptotic wave function with the solution of the Helmholtz equation indicates that the pΞ 0 interaction is attractive on the whole in both 1 S 0 and 3 S 1 channels. There is a slight tendency that the 3 S 1 interaction is more attractive than that in the 1 S 0 channel.
To reduce the uncertainties due to the lattice discretization at short distances, we need to carry out the simulations with smaller lattice spacing and/or with the improved lattice action. Systematic studies in various other channels such as ΛN , ΣN and ΛΛ are not only interesting by themselves but also important for studying the structure of hypernuclei and the interior of neutron stars. Also, (2+1)-flavor QCD simulations with realistic quark masses will be the ultimate goal to make qualitative predictions of the N N and Y N interactions: Studies along this line with the PACS-CS gauge configurations [34] are now under way.