Partial wave decomposition on the lattice and its applications to the HAL QCD method

TL;DR

This paper applies Misner's partial wave decomposition to lattice QCD data within the HAL QCD method, successfully extracting S-wave components and reducing artifacts, thereby improving the analysis of hadron-hadron interactions.

Contribution

It demonstrates the effectiveness of Misner's partial wave decomposition in lattice QCD calculations, enabling cleaner extraction of partial wave components in the HAL QCD method.

Findings

Misner's method successfully extracts the l=0 component from lattice data.

The method reduces unphysical structures in the potential caused by higher partial waves.

Results from Misner's method agree with conventional approaches on scattering phase shifts.

Abstract

The approximated partial wave decomposition method to the discrete data on a cubic lattice, developed by C. W. Misner, is applied to the calculation of $S$-wave hadron-hadron scatterings by the HAL QCD method in lattice QCD. We consider the Nambu-Bethe-Salpeter (NBS) wave function for the spin-singlet $\Lambda_c N$ system calculated in the $(2+1)$-flavor QCD on a $(32a~\mathrm{fm})^3$ lattice at the lattice spacing $a\simeq0.0907$ fm and $m_\pi \simeq 700$ MeV. We find that the $l=0$ component can be successfully extracted by Misner's method from the NBS wave function projected to $A_1^+$ representation of the cubic group, which contains small $l\ge 4$ components. Furthermore, while the higher partial wave components are enhanced so as to produce significant comb-like structures in the conventional HAL QCD potential if the Laplacian approximated by the usual second order difference is applied to the NBS wave function, such structures are found to be absent in the potential extracted by Misner's method, where the Laplacian can be evaluated analytically for each partial wave component. Despite the difference in the potentials, two methods give almost identical results on the central values and on the magnitude of statistical errors for the fits of the potentials, and consequently on the scattering phase shifts. This indicates not only that Misner's method works well in lattice QCD with the HAL QCD method but also that the contaminations from higher partial waves in the study of $S$-wave scatterings are well under control even in the conventional HAL QCD method. It will be of interest to study interactions in higher partial wave channels in the HAL QCD method with Misner's decomposition, where the utility of this new technique may become clearer.

Figures

Equations

$${\psi }^{A_1^{+}}(\vec{x})\equiv P^{A_1^{+}}\psi (\vec{x})=\frac{1}{48}\sum _{g\in O_h}\psi (g^{-1}\vec{x}),$$

$$\psi (\vec{x})=\sum _{l=0}^{\infty }\sum _{m=-l}^l{g}_{lm}(r)Y_{lm}(\theta ,\varphi ),$$

$$\left(E_k-H_0\right){\psi }^{(W_k)}(\vec{x})=\int d^3{x}^{\prime }U(\vec{x},\overrightarrow{x^{\prime }}){\psi }^{(W_k)}(\overrightarrow{x^{\prime }}),\left(E_k=\frac{|\vec{k}{|}^2}{2\mu },H_0=-\frac{{\vec{\nabla }}^2}{2\mu }\right),$$

$${\psi }^{(W_k)}(\vec{x})e^{-W_{k}t}=\frac{1}{\sqrt{Z_1}\sqrt{Z_2}}\sum _{\vec{y}}⟨0|B_1(\vec{x}+\vec{y},t)B_2(\vec{y},t)|2B;W_k⟩.$$

$$G(\vec{x},t-t_0)$$

$$G(\vec{x},t-t_0)\to \sqrt{Z_1}\sqrt{Z_2}{\psi }^{(W_0)}(\vec{x})e^{-W_0(t-t_0)}{A}_0+\mathcal{O}(e^{-W_1(t-t_0)}).$$

$$R(\vec{x},t-t_0)$$

$$\left[\left(\frac{1+3{\delta }^2}{8\mu }\right)\frac{{\partial }^2}{\partial t^2}-\frac{\partial }{\partial t}-H_0+\mathcal{O}(\frac{{\delta }^2}{2m_{B_1}{m}_{B_2}}\frac{{\partial }^3}{\partial t^3})\right]R(\vec{x},t-t_0)=\int d^3{x}^{\prime }U(\vec{x},\overrightarrow{x^{\prime }})R(\overrightarrow{x^{\prime }},t-t_0),$$

$$U(\vec{x},\overrightarrow{x^{\prime }})=V(\vec{x},\vec{\nabla }){\delta }^{(3)}(\vec{x}-\overrightarrow{x^{\prime }}),$$

$$V_{\mathrm{LO}}^{{}^1{S}_0}(r)=\frac{\left[\left(\frac{1+3{\delta }^2}{8\mu }\right)\frac{{\partial }^2}{\partial t^2}-\frac{\partial }{\partial t}-H_0\right]R^{{}^1{S}_0}(\vec{x},t-t_0)}{R^{{}^1{S}_0}(\vec{x},t-t_0)},$$

$$R^{{}^1{S}_0}(\vec{x},t-t_0)\equiv P^{A_1^{+}}{P}^{S=0}R(\vec{x},t-t_0).$$

$$g_{lm}(R)={\int }_{S}d\Omega \overline{Y_{lm}(\theta ,\varphi )}\psi (\vec{x};r=R),$$

$${\int }_{S}d\Omega \overline{Y_{lm}(\theta ,\varphi )}{Y}_{l^{\prime }{m}^{\prime }}(\theta ,\varphi )={\delta }_{l{l}^{\prime }}{\delta }_{m{m}^{\prime }},$$

$$\mathcal{G}\equiv \left\{{\mathcal{G}}_{lm,l^{\prime }{m}^{\prime }}\equiv \frac{1}{N_R}\sum _{\vec{x}\in \{\vec{x}|r=R\}}\overline{Y_{lm}(\theta ,\varphi )}{Y}_{l^{\prime }{m}^{\prime }}(\theta ,\varphi )||m|\le l,|m^{\prime }|\le l^{\prime }\right\},$$

$${\int }_{R-\Delta }^{R+\Delta }dr{r}^2\overline{G_n^{R,\Delta }(r)}{G}_m^{R,\Delta }(r)={\delta }_{nm}.$$

$$G_n^{R,\Delta }(r)=\frac{1}{r}\sqrt{\frac{2n+1}{2\Delta }}{P}_n\left(\frac{r-R}{\Delta }\right),$$

$$S_{R,\Delta }\equiv \left\{\vec{x}|R-\Delta \le r\le R+\Delta \right\},$$

$${\int }_{S_{R,\Delta }}{d}^{3}x\overline{{\mathcal{Y}}_{nlm}^{R,\Delta }(r,\theta ,\varphi )}{\mathcal{Y}}_{n^{\prime }{l}^{\prime }{m}^{\prime }}^{R,\Delta }(r,\theta ,\varphi )={\delta }_{n{n}^{\prime }}{\delta }_{l{l}^{\prime }}{\delta }_{m{m}^{\prime }},$$

$${\int }_{S_{R,\Delta }}{d}^{3}x\equiv {\int }_{R-\Delta }^{R+\Delta }{r}^{2}dr{\int }_{S}d\Omega .$$

$$\psi (\vec{x})=\sum _{n=0}^{\infty }\sum _{l=0}^{\infty }\sum _{m=-l}^l{c}_{nlm}^{R,\Delta }{\mathcal{Y}}_{nlm}^{R,\Delta }(r,\theta ,\varphi )$$

$$c_{nlm}^{R,\Delta }={\int }_{S_{R,\Delta }}{d}^{3}x\overline{{\mathcal{Y}}_{nlm}^{R,\Delta }(r,\theta ,\varphi )}\psi (\vec{x}).$$

$$g_{lm}(r)=\sum _{n=0}^{\infty }{c}_{nlm}^{R,\Delta }{G}_n^{R,\Delta }(r).$$

$${\int }_{S_{R,\Delta }}{d}^{3}x⟹\sum _{\vec{x}}{\omega }^{R,\Delta }(\vec{x}),$$

$${\omega }^{R,\Delta }(\vec{x})=\{\begin{array}{ll}{a}^3& \text{for ∣r−R∣<Δ−21a},\\ 0& \text{for ∣r−R∣>Δ+21a},\\ a^2\left(\Delta +\frac{1}{2}a-|R-r|\right)& \text{otherwise,}\end{array}$$

$$⟨f|g{⟩}_{S_{R,\Delta }}$$

$${\tilde{\mathcal{Y}}}_A^{R,\Delta }(\vec{x})\equiv {\sum _B}^{\prime }{\mathcal{Y}}_B^{R,\Delta }(\vec{x}){\mathcal{G}}_{BA}^{-1},$$

$$⟨{\tilde{\mathcal{Y}}}_A^{R,\Delta }|{\mathcal{Y}}_B^{R,\Delta }{⟩}_{S_{R,\Delta }}$$

$$\psi (\vec{x})\simeq \sum _{n=0}^{n_{\max }}\sum _{l=0}^{l_{\max }}\sum _{m=-l}^l{c}_{nlm}^{R,\Delta }{\mathcal{Y}}_{nlm}^{R,\Delta }(r,\theta ,\varphi )$$

$$g_{lm}(R)$$

$$\mathbf{\Psi }$$

Full text

1 Introduction

A determination of hadron-hadron interactions from the first-principle is one of the ultimate goals in both particle and nuclear physics. In a lattice QCD calculation of a two-hadron system, the quantum numbers of the system are specified by the source and/or sink operators in the corresponding correlation function. Among quantum numbers, the partial wave can be specified by the angular dependence in terms of the relative coordinate of two hadron operator, in principle. In lattice QCD, however, the extraction of the designated partial wave becomes non-trivial, because the rotational symmetry is broken to the cubic symmetry due to the finite volume (IR-effect) as well as the finite lattice spacing (UV-effect), which introduces the mixing between different partial waves [1, 2]. In addition, full account of the angular dependence cannot be obtained since the data are available only on discretized spacial coordinates.

In Ref. [3], a general method to (approximately) obtain a radial function in the particular partial wave on the cubic lattice was proposed by C. W. Misner, and it has been applied to the analyses of gravitational waves simulated on the grid points [4]. Applying this method to lattice QCD enables us to verify how much the partial wave mixing is induced by the breaking of the rotational symmetry and evaluate the corresponding systematics in the calculations of hadron-hadron interactions. In addition, this method could open a new possibility to extract interactions in higher partial waves which are otherwise difficult to be studied in lattice QCD calculations.

The HAL QCD method is a promising method to calculate the hadron-hadron interactions in lattice QCD, in which we construct the hadron-hadron “potential” from the Nambu-Bethe-Salpeter (NBS) wave function and the physical observables such as scattering phase shifts are calculated by solving the Schrödinger equation with the potential in the infinite volume [5, 6, 7, 8]. One of the unique features in the procedure of the HAL QCD method is that the spatial correlation of the NBS wave function is calculated at all discrete points on a cubic lattice, and thus the hadron-hadron interactions from the HAL QCD method are expected to be a good application of Misner’s method (or Misner’s decomposition/extraction). In fact, noticing that not all discrete points at a given radial coordinate $r$ are necessarily transformed to each other by the cubic rotation, it is realized that there is a room to develop a better methodology for the partial wave decomposition than the standard projection method based on the irreducible representation of the cubic group.

In this paper, we apply Misner’s method for the first time to the lattice QCD study for hadron-hadron scatterings in the framework of the HAL QCD method. We extract the potential from the $l=0$ ($S$-wave) component of the NBS wave function by Misner’s method, which is then compared with the conventional HAL QCD potential.

In the conventional method, the $S$-wave projection of the NBS wave function is approximated on the lattice by the $A_{1}^{+}$ projection as

[TABLE]

where the cubic group $O_{h}$ consists of cubic rotations and the parity. The $A_{1}^{+}$ representation contains not only $l=0$ component but also higher partial waves with $l\geq 4$. If $l\geq 4$ components of the NBS wave function were absent, the NBS wave function (and also the potential) would be isotropic and thus depend only on the radial coordinate $r=|\vec{x}|$. However, we often observe comb-like structures in the potential in terms of $r$ (for example, see Fig. 2 in Ref. [9]), which represent the anisotropy of the potential. This observation indicates the existence of non-negligible $l\geq 4$ components. We often observe that these comb-like structures lead to superficial fluctuations, whose magnitudes are larger than those of the genuine statistical fluctuations.

Throughout this paper, we write the NBS wave function as $\psi(\vec{x})$ with $\vec{x}=(x,y,z)$ or $(r,\theta,\phi)$, which is expanded in term of the spherical harmonics $Y_{lm}(\theta,\phi)$ as

[TABLE]

where we call $g_{lm}(r)$ “spherical harmonics amplitude” for a $(l,m)$ component.

This paper is organized as follows. After briefly reviewing the HAL QCD method in Sec. 2, we explain Misner’s method in detail in Sec. 3, together with some remarks on its application to the HAL QCD method. Our main results are given in Sec. 4, where we consider the spin-singlet $\Lambda_{c}N$ system as a representative example. In Sec. 4.1, we extract the $l=0$ component of the NBS wave function by Misner’s method, where the comb-like structures indeed disappear. We however found that contaminations from $l\geq 4$ partial waves are small in the $A_{1}^{+}$ projected NBS wave function. In Sec. 4.2, we analyze the Laplacian of the NBS wave function by Misner’s method. We found that $l\geq 4$ components are enhanced by applying Laplacian to the NBS wave function in the conventional HAL QCD method, while such a problem is absent in Misner’s method, where the Laplacian is calculated analytically for each partial wave component. In Sec. 4.3, we investigate parameter dependencies of the potential in Misner’s method. In Sec. 5, we calculate the scattering phases shifts from the HAL QCD potentials with the conventional $A_{1}^{+}$ projection and with Misner’s $S$-wave extraction. We found that not only the central values but also statistical errors agree in both cases. We briefly discuss a reason for this agreement. Summary and conclusion are presented in Sec. 6. In appendix A, a simpler but less general method is considered to extract the $l=0$ component from the $A_{1}^{+}$ projected NBS wave function.

2 HAL QCD method

In the HAL QCD method [5, 6, 7, 8], a non-local but energy-independent potential is defined through the Schrödinger equation as

[TABLE]

where the NBS wave function in the center-of-mass frame is given by

[TABLE]

Here we consider a two-baryon system as a representative case where $B_{i}(\vec{x},t)~{}(i=1,2)$ is the local interpolating operator for a baryon $B_{i}$ with its renormalization factor $\sqrt{Z_{i}}$. The state $|2B;W_{k}\rangle$ stands for a QCD eigenstate for the two-baryon system at the total energy $W_{k}=\sqrt{|\vec{k}|^{2}+m_{B_{1}}^{2}}+\sqrt{|\vec{k}|^{2}+m_{B_{2}}^{2}}$ with baryon masses $m_{B_{1,2}}$ and a relative momentum $\vec{k}$, and $\mu$ denotes a reduced mass. Since the asymptotic behavior of the NBS wave function at large $r=|\vec{x}|$ is identical to that of the scattering wave in quantum mechanics, whose phase shift is the phase of QCD $S$-matrix [10, 11], the potential defined from the NBS wave functions reproduces the scattering phase shifts in QCD. Note that the non-local potential is constructed so as to be energy-independent below the inelastic threshold [6, 8].

In terms of the NBS wave functions, the two-baryon four-point correlation function is expressed as

[TABLE]

where $\mathcal{J}^{(J^{P})}(t_{0})$ stands for a source operator defined so as to create two-baryon states at $t=t_{0}$ with the total angular momentum $J$ and the parity $P$, and the ellipses represent contributions from inelastic states. For a sufficiently large $t-t_{0}$, the NBS wave function for the ground state is extracted from the four-point function as

[TABLE]

In practice, however, this extraction of the NBS wave function from the ground state saturation is very difficult due to increasing statistical noises at large $t-t_{0}$ [12, 13, 14]. Therefore, in Ref. [7], an improved method was proposed to extract the potential without the requirement of the ground state saturation. We consider the normalized baryon four-point correlation function ($R$-correlator) defined by

[TABLE]

where $\Delta W_{n}=W_{n}-(m_{B_{1}}+m_{B_{2}})$. If contributions from inelastic states are negligible (“elastic state saturation”), this $R$-correlator satisfies

[TABLE]

where $\delta=(m_{B_{1}}-m_{B_{2}})/(m_{B_{1}}+m_{B_{2}})$. Since the elastic state saturation can be generally achieved at much smaller $t-t_{0}$ than the case of the ground state saturation, we can obtain reliable results from this “time-dependent HAL QCD method” [12, 13, 14, 15].

In order to handle the non-locality of the potential, we introduce the derivative expansion as

[TABLE]

where $V(\vec{x},\vec{\nabla})$ is expanded in terms of $\vec{\nabla}$ [16]. At low energies, since the leading order (LO) potential of the derivative expansion dominates [15], the interaction for the $S$-wave spin singlet system is well approximated by the LO central potential given as

[TABLE]

where the $R$-correlator for the ${}^{1}S_{0}$ state is defined as

[TABLE]

Here $P^{S=0}$ represents the projection to the state with the total spin $S=0$, while $P^{A_{1}^{+}}$ stands for the projection to the $A_{1}^{+}$ representation in Eq. (1). As we have explained before, the $A_{1}^{+}$ projection contains not only $l=0$ component but also $l\geq 4$ components, which produces angular dependencies of the $R$-correlator as well as the potential.

3 Approximated partial wave decomposition

3.1 Misner’s method

Let us consider the extraction $g_{lm}(r)$ for a given $r\equiv|\vec{x}|=R$ from the NBS wave function. In the continuum space, we can obtain $g_{lm}(R)$ by taking the spherical surface integral at $r=R$ as

[TABLE]

because of the orthogonality of the spherical harmonics,

[TABLE]

where the overline represents its complex conjugation. This is even true on a finite $L_{s}^{3}$ box in the continuum space for $0<R\leq L_{s}/2$.

In the discrete space such as the cubic lattice, however, we can not obtain $g_{lm}(R)$ for any $R$, since $N_{R}$, the number of points which satisfy $r=R$, is finite, so that the infinite dimensional matrix,

[TABLE]

is non-invertible due to its zero eigenvalues. If one knows that contributions from higher partial waves are negligible, one can introduce some approximation in order to obtain $g_{lm}(R)$ for small $l$ and some restricted $R$. For example, we may impose the condition that $l,l^{\prime}\leq l_{\rm max}$, where $l_{\rm max}$ is chosen so that the finite dimensional matrix ${\cal G}^{l_{\rm max}}\equiv\{{\cal G}_{lm,l^{\prime}m^{\prime}}|l,l^{\prime}\leq l_{\rm max}\}$ becomes invertible. In appendix A, we consider this approximation in detail.

In Ref. [3], a more general and sophisticated approximation was proposed. To explain the method, let us start from the continuum case. We first introduce the basis functions in the radial coordinate $G^{R,\Delta}_{n}(r)$ ($n=0,\cdots,\infty$) which are orthonormal in the radial interval $[R-\Delta,R+\Delta]$ as

[TABLE]

One of the candidates for $G^{R,\Delta}_{n}(r)$ is given by

[TABLE]

where $P_{n}(r)$ is the Legendre polynomial, which obviously satisfies Eq. (15). When we consider a spherical shell $S_{R,\Delta}$ with thickness $2\Delta$ surrounding the sphere surface $r=R$ defined by

[TABLE]

an orthonormal basis function $\mathcal{Y}^{R,\Delta}_{nlm}(r,\theta,\phi)\equiv G^{R,\Delta}_{n}(r)Y_{lm}(\theta,\phi)$ obeys

[TABLE]

where the integral over $S_{R,\Delta}$ is defined as

[TABLE]

The NBS wave function in a spherical shell $S_{R,\Delta}$ is expanded in terms of the orthonormal basis functions as

[TABLE]

with coefficients $c^{R,\Delta}_{nlm}$, which can be determined by

[TABLE]

We finally obtain the spherical harmonics amplitude $g_{lm}(r)$ for $R-\Delta\leq r\leq R+\Delta$ as

[TABLE]

Ref. [3] employed ${\cal Y}^{R,\Delta}_{n,l,m}(r,\theta,\phi)$ as the basis function for the approximation in the case of the discrete space on a cubic lattice. In this case, the volume integral is replaced by the discrete sum as

[TABLE]

where $\omega^{R,\Delta}({\vec{x}})$ is a weight factor, which corresponds to a volume of the overlapped region between the shell $S_{R,\Delta}$ and a unit cube around the point $\vec{x}$. For example, if the unit cube lies entirely inside the shell $S_{R,\Delta}$, $\omega^{R,\Delta}({\vec{x}})=a^{3}$ with a lattice spacing $a$, while $\omega^{R,\Delta}({\vec{x}})=0$ when the unit cube lies entirely outside the shell. Since $\omega^{R,\Delta}({\vec{x}})$ in the general cases is rather complicated, it is approximated in Ref. [3] as

[TABLE]

which corresponds to the overlapped volume of a unit cube parallel to the radial direction. Using this, we define an inner product of functions $f(\vec{x})$ and $g(\vec{x})$ in the shell $S_{R,\Delta}$ as

[TABLE]

Let us consider ${\cal G}_{AA^{\prime}}\equiv\langle{\cal Y}_{A}^{R,\Delta}|{\cal Y}_{A^{\prime}}^{R,\Delta}\rangle_{S_{R,\Delta}}$ with a shorthand notation $A=(n,l,m)$. The finite dimensional Hermitian matrix ${\cal G}$ constructed from ${\cal G}_{AA^{\prime}}$ with a restriction that $n,n^{\prime}\leq n_{\rm max}$ and $l,l^{\prime}\leq l_{\rm max}$ becomes invertible if one properly chooses $n_{\rm max}$ and $l_{\rm max}$. Using ${\cal G}$ (whose dependencies on $n_{\rm max}$ and $l_{\rm max}$ are implicit here), one can define the dual basis functions $\tilde{\mathcal{Y}}^{R,\Delta}_{A}$ as

[TABLE]

which satisfies

[TABLE]

where the prime in the summation indicates the upper bounds $n_{\rm max}$ and $l_{\rm max}$.

Assuming that $c^{R,\Delta}_{nlm}$ is negligibly small for $l>l_{\rm max}$ or $n>n_{\rm max}$, Eq. (20) is approximately written as

[TABLE]

with the coefficient $c^{R,\Delta}_{nlm}=\langle\tilde{\mathcal{Y}}^{R,\Delta}_{nlm}|\psi\rangle_{S_{R,\Delta}}$. Finally the spherical harmonics amplitude $g_{lm}(R)$ can be approximated as

[TABLE]

3.2 Misner’s method as a minimization

Misner’s method can be also understood from the viewpoint of the least square minimization and we here give the explicit correspondence following Ref. [17].

Let us denote a $N_{R,\Delta}$ component vector $\bm{\Psi}$ of the NBS wave function as

[TABLE]

where $N_{R,\Delta}$ is the number of points in the shell $S_{R,\Delta}$, equivalently, the number of data with non-zero $\omega^{R,\Delta}(\vec{x})$. Similarly, we define a $N_{R,\Delta}\times M$ rectangular matrix $Y$ of the basis functions, whose components are defined by

[TABLE]

where $M$ is the number of the basis functions and is given by $M=(n_{\rm max}+1)(l_{\rm max}+1)^{2}$. For $n_{\mathrm{max}}=2$ and $l_{\mathrm{max}}=2$, for example, $Y$ becomes

[TABLE]

where the number of columns is $M=27$. Defining a $N_{R,\Delta}\times N_{R,\Delta}$ diagonal matrix $W$ for the non-zero weight $\omega^{R,\Delta}(\vec{x})$, the $M\times M$ matrix $\bm{\mathcal{G}}\equiv\{\mathcal{G}_{AB}\}$ is simply written as $\bm{\mathcal{G}}=$ $Y$†$WY$.

Using these notations, we introduce a trial $N_{R,\Delta}$-component vector function \bm{\tilde{\Psi}}\equiv$$Y\tilde{C} as

[TABLE]

where a $M$-component vector $\tilde{C}$ corresponds to fit parameters that should minimize

[TABLE]

Since dF(\mbox{\boldmath\tilde{C}})/d\mbox{\boldmath\tilde{C}}=0 at the minimum \mbox{\boldmath\tilde{C}}=\mbox{\boldmath\tilde{C}}_{\rm min}, we obtain

[TABLE]

so that

[TABLE]

which agrees with Misner’s method, Eq. (28). Therefore, Misner’s method is equivalent to finding a solution of $\tilde{\psi}(\vec{x}_{i})$ which minimizes the difference between the data $\psi(\vec{x}_{i})$ and the fit function $\tilde{\psi}(\vec{x}_{i})$ defined by the norm $\langle\tilde{\psi}-\psi|\tilde{\psi}-\psi\rangle_{S_{R,\Delta}}$.

3.3 Remarks

The calculation of the potential in the HAL QCD method requires the Laplacian applied to the NBS wave function, which is conventionally approximated by a finite difference, and thus contains discretization errors. In the application of Misner’s decomposition to the HAL QCD method, we can instead employ an analytical expression for the Laplacian, which operates on the (approximately obtained) partial wave $g_{lm}(r)Y_{lm}(\theta,\phi)$ as

[TABLE]

where

[TABLE]

and $P^{\prime\prime}_{n}$ is the second-order derivative for the Legendre polynomial. Unlike the conventional HAL QCD method in which the difference operator for the Laplacian is applied to (all partial wave components of) the NBS wave function, it is clear that the analytic derivative in Misner’s method does not induce contributions from other partial waves than the targeted one. A comparison between two implementations for the Laplacian operator will be given in Sec. 4.2.

In Misner’s method, it is practically important to choose $n_{\rm max}$, $l_{\rm max}$ and $\Delta$ appropriately. While larger $n_{\rm max},l_{\rm max}$ and smaller $\Delta$ gives a better approximation of the NBS wave function, it leads to a small $N_{R,\Delta}-M$ that may cause some numerical instability due to small eigenvalues of $\mathcal{G}_{AB}$111If $N_{R,\Delta}-M<0$ (and there is no symmetry), $\mathcal{G}_{AB}$ has zero eigenvalues. or may give an over-fitting, where $N_{R,\Delta}$ and $M$ correspond to the numbers of data and fit parameters, respectively. For example, if the variation of the NBS wave function in the radial coordinate is large, one should increase $n_{\mathrm{max}}$ to approximate the spherical harmonics amplitude $g_{lm}(r)$ better, but not too much so as to avoid the instability or the over-fitting. In Ref. [18], the scaling of the discretization error in Misner’s method is discussed. It is found that the error depends on $\Delta$ as ${\cal O}(\Delta^{n_{\mathrm{max}}+2})$, which also indicates that the choice of $\Delta={\cal O}(a)$ is preferable. In addition, the volume integral in the shell with the approximated weight (Eq. (24)) gives ${\cal O}(a^{2})$ error with the choice of $\Delta={\cal O}(a)$. Therefore, the choice of the parameters $n_{\mathrm{max}}=2$ with $\Delta={\cal O}(a)$ is found to be good for the second order accuracy.

The author [18] also suggests $\Delta=3a/4$ as a rule of thumb by numerical investigations, but we have to examine whether the results are stable against the change of parameters $\Delta,n_{\mathrm{max}},l_{\mathrm{max}}$ case by case, as will be presented in the next section.

In practice, the most costly calculation in Misner’s method is the construction of the matrix $\mathcal{G}_{AB}$. Once we calculate the matrix, however, we can use it for different lattice data (e.g. NBS wave functions calculated on different gauge samples). Therefore, it is better to calculate the dual basis functions (Eq. (26)) once before the calculation of the spherical harmonics amplitude from NBS wave functions and use them for the later analyses. One possible obstacle in this procedure is that the dual basis functions $\tilde{\mathcal{Y}}^{R,\Delta}_{nlm}(\vec{x})$ consume large amount of memory to store, $L^{3}(n_{\rm max}+1)(l_{\rm max}+1)^{2}\times 16$ Bytes for each given value of $R$. For example, the required memory for $\tilde{\mathcal{Y}}^{R,\Delta}_{nlm}(\vec{x})$ with $L=32$, $n_{\mathrm{max}}=4$ and $l_{\mathrm{max}}=6$ becomes $32^{3}\times 5\times 7^{2}\times 16$ (Bytes) $=122.5$ MB. In order to reduce the memory usage by a factor of $(n_{\rm max}+1)$, we instead store

[TABLE]

which needs only $24.5$ MB. Using this function, the spherical harmonics amplitude $g_{lm}(R)$ can be calculated directly as

[TABLE]

Furthermore, it is sufficient to store $F^{R,\Delta}_{lm}(\vec{x})$ only at point $\vec{x}$ included in the shell $S_{R,\Delta}$ where the weight function $\omega^{R,\Delta}({\vec{x}})$ is non-zero, which leads to extra large reduction for the memory usage.

4 HAL QCD potentials with Misner’s method

For the numerical calculation, we consider the spin-singlet $\Lambda_{c}N$ system in the $(2+1)$-flavor full lattice QCD with the renormalization-group improved Iwasaki gluon action and a nonperturbatively ${\cal O}(a)$ improved Wilson-clover quark action on a $(32a~{}\mathrm{fm})^{3}\times(64a~{}\mathrm{fm})$ volume with the lattice spacing $a\simeq 0.0907$ fm at $m_{\pi}\simeq 700$ MeV. We apply Misner’s method to the same data of NBS wave function calculated in Ref. [19], where the results in the conventional HAL QCD method are given. For the source operator, we employ the wall-type source operator and thus the $\Lambda_{c}N$ system belongs to the $A_{1}^{+}$ representation of the cubic group. In order to reduce the statistical fluctuations, we also impose the $A_{1}^{+}$ projection on the sink operator as given in Eq. (1). Recall that the $A_{1}^{+}$ representation contains the partial waves $l=0,4,6,\cdots$. In this study, we consider the $A_{1}^{+}$ projected $R$-correlator taken at $t-t_{0}=13a$. The total number of configuration is 399, and the statistical errors are estimated by the jackknife method with a bin size of 57 configurations (the total number of bins is 7). For more details on the simulation setup, see Ref. [19].

4.1 NBS wave function for the spin-singlet $\Lambda_{c}N$ system

Fig. 1 shows the results of the $R$-correlator defined in Eq. (7) for spin-singlet $\Lambda_{c}N$ system. The gray points represent the $A_{1}^{+}$ projected $R$-correlator divided by $Y_{00}$, which is actually used in Ref. [19] to construct the $\Lambda_{c}N$ potential. The red points correspond to the spherical harmonics amplitude $g_{lm}(r)$ for $l=m=0$ component calculated by Misner’s method for the radial coordinate $r$ from $2a$ to $16a$ with the interval $\Delta r=0.2a$, so that some data are used several times. We do not perform Misner’s method for $r<2a$ and $r>L_{s}/2=16a$: In the former case, the number of data points in the spherical shell $S_{r,\Delta}$ is too small, whereas the spherical shell is not contained in the $L_{s}^{3}$ cubic lattice for the latter. We here employ $\Delta=a$, $n_{\mathrm{max}}=2$ and $l_{\mathrm{max}}=4$ as the parameters in Misner’s method, and we found that $g_{00}(r)$ has a weak parameter dependence.

Fig. 1 shows small comb-like structures in the $A_{1}^{+}$ projected $R$-correlator, which however do not appear in the $l=0$ component extracted by Misner’s method. This observation indicates that $l\geq 4$ components exist in the $A_{1}^{+}$ projected $R$-correlator and their angular dependencies become manifest as the comb-like structures in the radial-coordinate. Such higher partial wave components can be explicitly extracted by Misner’s method as shown in Fig. 2 for the $l=4$ component ($g_{40}(r)=g_{4\>\pm 4}(r)$). Note that $g_{4m}(r)=0$ for $m\neq 0,\pm 4$ for the $A_{1}^{+}$ representation. We find that $l=4$ component indeed exists while its magnitude is small (by a factor of ${\cal O}(10^{-3})$ compared to that of the $l=0$ component). On the other hand, the absence of comb-like structures in $l=0,4$ components obtained by Misner’s method with $l_{\mathrm{max}}=4$ indicates that $l\leq 4$ components are sufficient to explain the $A_{1}^{+}$ projected $R$-correlator, which is explicitly confirmed by observing that $l=6$ component extracted by Misner’s method with $l_{\mathrm{max}}=6$ is actually negligible.

The mixing of $l=4$ component is most likely induced due to the rotational symmetry breaking by the finite volume (IR-effect), except for $r\lesssim a$ where there could also exist the effect by the finite lattice spacing (UV-effect).

4.2 Laplacian and HAL QCD potential

We then study the effect of Laplacian applied to the NBS wave function. Note that the term containing the Laplacian, i.e., the third term in the rhs of Eq. (10), is known to give the dominant contribution for the potential. In the conventional HAL QCD method, the Laplacian approximated by a finite second-order difference is applied to (all partial wave components of) the NBS wave function, while it is analytically calculated as Eq. (42) in Misner’s method for the designated partial wave component. In Fig. 3, we compare the Laplacian applied to the NBS wave functions between Misner’s method (red points) and the conventional HAL QCD method (gray points). In the case of Misner’s method, the Laplacian applied to $g_{00}(r)$ analytically does not exhibit the comb-like structure. In the case of the conventional HAL QCD method, on the other hand, we find that the comb-like structures are much larger than those of the $R$-correlator itself, indicating that the $l\geq 4$ components are larger for the conventional Laplacian. The partial wave decomposition of the conventional Laplacian term reveals that the $l=4$ component is indeed larger than the case of the $R$-correlator as shown in Fig. 4. The origin of these $l\geq 4$ components in the conventional Laplacian is most likely the $l\geq 4$ components in the $R$-correlator enhanced by the Laplacian, rather than the discretization error in the conventional Laplacian operator itself. In fact, in the case of Misner’s method, the difference between the results with the conventional Laplacian operator applied to $g_{00}(r)$ and those with the analytic Laplacian operator is found to be marginal.

Shown in Fig. 5 are the HAL QCD potentials for the spin-singlet $\Lambda_{c}N$ system from Misner’s $S$-wave extraction (red points) and the conventional $A_{1}^{+}$ projection (gray points). In the case of Misner’s extraction, the potential is found to be free from comb-like structures. In the case of the conventional projection, however, the potential has large comb-like structures. The main origin is attributed to $l\geq 4$ components in the Laplacian term for the potential, which are enhanced by applying the Laplacian to the $R$-correlator, even though $l\geq 4$ components are small in the $R$-correlator itself.

4.3 Parameter dependencies for potentials in Misner’s method

We here discuss dependencies of potentials on various parameters in Misner’s method, $\Delta,n_{\mathrm{max}}$ and $l_{\mathrm{max}}$, which correspond to the thickness of the spherical shell, the maximum number of bases for the radial function and that for the spherical harmonics, respectively. Throughout this section, the potential is constructed from the $l=0$ component of the NBS wave function.

We first show the $n_{\mathrm{max}}$ dependence with other parameters fixed to $\Delta=a$ and $l_{\mathrm{max}}=4$, in upper two figures of Fig. 6, where we vary the value of $n_{\mathrm{max}}$ from $1$ to $5$. In Fig. 6, we plot only the central values of the potentials without statistical errors in order to make it easier to see the dependence on parameters in Misner’s method. Note that the magnitude of statistical errors are found to be stable against changing parameters. For $n_{\mathrm{max}}=1$, the potential does not reproduce the correct behavior. The reason is that the second derivative of the Legendre polynomial $P_{n}^{\prime\prime}(x)$ is zero for $n=0$ and $1$ so that the Laplacian term in the potential from the spherical harmonics amplitude vanishes. The small contribution to the potential at $n_{\mathrm{max}}=1$ comes from the time-derivative terms (the first and second terms in Eq. (10)). Thus it is necessary to take $n_{\mathrm{max}}\geq 2$, for which we find that the potentials are almost stable against the change of $n_{\mathrm{max}}$. While we observe small oscillations of the potential for $n_{\mathrm{max}}\geq 4$, they are probably due to the numerical instabilities associated with small eigenvalues of $\mathcal{G}_{AB}$ caused by a large number of $n_{\mathrm{max}}$, as discussed in Sec. 3.3. Recall also that the discretization errors of the radial orthonormal basis function $G^{R,\Delta}_{n}(r)$ is known to be $\mathcal{O}(\Delta^{n_{\mathrm{max}}+2})$ [18]. Since the discretization errors in our lattice QCD action [19] is $\mathcal{O}(a^{2})$, a choice of the parameters $n_{\mathrm{max}}=2$ with $\Delta=\mathcal{O}(a)$ is reasonable.

We next present the $\Delta$ dependence of the potential in middle two figures of Fig. 6, where we take $\Delta=0.2a,0.5a,a,1.5a$ and $2a$ with $n_{\mathrm{max}}=2$ and $l_{\mathrm{max}}=4$ fixed. For $\Delta=0.2a$ and $0.5a$, we find the comb-like structures even in the potential constructed from the $l=0$ component, while such structures are absent for the potential with $\Delta\geq a$. This can be understood from the fact that the number of the data points in the spherical shell becomes too small for small $\Delta$ to reproduce the spherical harmonics amplitude accurately. For larger $\Delta$, on the other hands, there appear small deviations from the potential with $\Delta=a$ at short distances, while the potentials at long distances are stable against the change of $\Delta$. This may be explained by the fact that variations of the $R$-correlator in the spherical shell become sizable for larger $\Delta$, so that we need to enlarge $n_{\mathrm{max}}$ accordingly to approximate the spherical harmonics amplitudes precisely. In fact, by taking larger value of $n_{\mathrm{max}}$ in the case of $\Delta=1.5a,2a$, we find that the results tend to converge to that with $\Delta=a,n_{\mathrm{max}}=2$. From these observations, we take $\Delta=a$ and $n_{\rm max}=2$ in this paper.

Finally, we investigate the $l_{\mathrm{max}}$ dependence with $\Delta$ and $n_{\mathrm{max}}$ fixed, as shown in lower two figures of Fig. 6. We take $l_{\mathrm{max}}=0,4,6$ and $8$, while keeping $\Delta=a$ and $n_{\mathrm{max}}=2$. The potentials are stable against the change of $l_{\mathrm{max}}$. While the potential is rather reasonable even in the case of $l_{\mathrm{max}}=0$, small oscillations are observed in the potential for this case. Such oscillations are absent for $l_{\mathrm{max}}\geq 4$, indicating that $l_{\mathrm{max}}=4$ component has small but non-negligible contributions in the $R$-correlator, while $l_{\mathrm{max}}>4$ components are sufficiently small. The $l=6$ component in the $R$-correlator obtained with $l_{\mathrm{max}}\geq 6$ is actually found to be negligible. Therefore we take $l_{\mathrm{max}}=4$ in this paper, as a conservative choice to avoid numerical instabilities for larger $l_{\mathrm{max}}$.

5 Phase shifts for the spin-singlet $\Lambda_{c}N$ system

We here compare the scattering phase shifts calculated from the HAL QCD potential obtained by Misner’s $S$-wave extraction with those by the conventional $A_{1}^{+}$ projection, in order to estimate effects from $l\geq 4$ partial waves to physical observables. For this purpose, we fit both potentials using the fit function

[TABLE]

with the fitting range $r\in[2a,16a]$, where both fit and original data are shown in Fig. 7. While the conventional HAL QCD potential has large comb-like structures, the fit parameters $a_{i}$ are almost identical for both potentials, and, more surprisingly, the magnitude of statistical errors are also found to be similar. This observation indicates that $l\geq 4$ contributions in the conventional HAL QCD potential hardly affect the fit of the potential. The agreement for the fit parameters between two cases is most likely attributed to that the fit function for the potential is taken to be isotropy (See Eq. (46)). In the fit, we employ the uncorrelated fit and a more systematic study with the correlated fit is left for future studies.

By solving the Schrödinger equation numerically with the fitted potentials, we extract the scattering phase shifts, which are shown in Fig. 8. As expected from the fit results of the potentials, not only the central values of the scattering phase shifts but also their statistical errors are almost identical between two methods.

In the analysis with the conventional HAL QCD method, while the potential is affected by contaminations from higher partial waves with $l\geq 4$, it is confirmed that the results of the scattering phase shifts are not affected by such systematics for both central values and the magnitude of errors. The conventional HAL QCD potential sometimes shows large fluctuations, in particular for $NN$ channels at lighter pion masses. The results in this paper tell us that these fluctuations are mainly systematic ones due to contaminations of higher partial waves enhanced by the second difference approximation of the Laplacian. By removing higher partial wave components from the $R$-correlator and calculating the Laplacian analytically, the analysis with Misner’s method reveals genuine statistical errors of the potential. An agreement in errors of the scattering phase shifts (or equivalently the fitted potential) between two analyses with the conventional method and Misner’s method provides a valuable check of the results.

6 Summary and conclusion

In this paper, we have performed the approximated partial wave decomposition by Misner’s method to lattice QCD data in order to extract the $l=0$ component of the $A_{1}^{+}$ projected NBS wave functions and its Laplacian for the $\Lambda_{c}N$ system in the spin-singlet channel, calculated in the $(2+1)$-flavor QCD on $(32a~{}\mathrm{fm})^{3}\times(64a~{}\mathrm{fm})$ at $m_{\pi}\simeq 700$ MeV [19].

We obtain the following results. While the $A_{1}^{+}$ projected NBS wave functions contain small contaminations from $l\geq 4$ partial waves, such contaminations are enhanced if the Laplacian approximated by the second order difference is applied to the NBS wave function, which cause large fluctuations in the conventional HAL QCD potential. With the use of Misner’s method, since the Laplacian can be calculated analytically for the designated ($l=0$) partial wave component in the NBS wave function, the potential is free from contaminations from higher partial waves, and thus its fluctuations become much smaller. Therefore, Misner’s method is very useful to reduce superficial fluctuations of the potential. If we fit the potentials as a function of $r$, not only the central values but also statistical errors of the fit parameters are almost independent of whether the conventional HAL QCD potential or the potential extracted with Misner’s method are used as input. Consequently, the scattering phase shifts agree between two methods. This agreement demonstrates not only that Misner’s method works well in the HAL QCD method but also the contaminations from higher partial waves in the study of $S$-wave scatterings are well under control even in the conventional HAL QCD method.

Since one can approximately obtain the spherical harmonics amplitude for an arbitrary $l$ component by Misner’s method, it is interesting to apply the method to extract the potentials from higher partial wave channels. For example, in order to extract the tensor potential, one needs to obtain $l=0$ and $l=2$ components of the NBS wave function separately. In the conventional HAL QCD method, one extracts the $l=2$ component from the NBS wave function by the projection, $(1-P^{A_{1}^{+}})$, which however also contains $l\geq 4$ components. By employing Misner’s method, on the other hand, one can extract $l=0$ and $l=2$ components separately without contaminations from higher partial waves, so as to obtain the tensor potential as well as the central potential without comb-like structures.

——————————————–

Acknowledgments

Numerical data used in this study were obtained by the KEK supercomputer system (BG/Q) [Project number: 14/15-21, 15/16-12]. This work is supported in part by JSPS Grant-in-Aid for Scientific Research, No. JP19K03879, JP18H05236, JP18H05407, JP16H03978, JP15K17667, by a priority issue (Elucidation of the fundamental laws and evolution of the universe) to be tackled by using Post “K” Computer, and by Joint Institute for Computational Fundamental Science (JICFuS). The authors thank all the members of the HAL QCD Collaboration for discussion.

Appendix

Appendix A Approximated partial wave decomposition at fixed $r$

In this appendix, we propose a simpler method to extract $g_{00}(r)$ from the discrete data, which is compared with Misner’s method. We here consider the $A_{1}^{+}$ projected NBS wave function defined by

[TABLE]

where $P^{A_{1}^{+}}$ is the projection operator to the $A_{1}^{+}$ representation, $Y_{00}^{A_{1}^{+}}(\theta,\phi)$ and $Y_{4m}^{A_{1}^{+}}(\theta,\phi)$ stand for the $A_{1}^{+}$ projected spherical harmonics, given by

[TABLE]

and the ellipsis denotes higher angular momentum components such as $l=6,8,\cdots$. In the $A_{1}^{+}$ representation, the $l=4$ components become non-zero only for $m=0,\pm 4$. Since $Y_{40}^{A_{1}^{+}}$ (also $Y_{4,\pm 4}^{A_{1}^{+}}$) has an angular dependence, the NBS wave function at given $r$ is multi-valued, which make comb-like structures in the potential on the radial coordinate.

Let us assume that there are $N$ points $\vec{x}_{i}$ $(i=1,\cdots,N)$ which satisfy $|\vec{x_{i}}|=r$ but can not be transformed each other by the cubic rotation. Neglecting components with $l\geq 6$, we have

[TABLE]

where we omit the arguments for the constant $Y_{00}^{A_{1}^{+}}$, and $g_{4}(r)\equiv g_{40}(r)+\sqrt{\frac{5}{14}}(g_{44}(r)+g_{4-4}(r))$. Eq. (51) can be compactly written as

[TABLE]

where the matrix in the right-hands side is an $N\times 2$ rectangular matrix in general. If $N>2$, we solve Eq. (60) by using the Singular Value Decomposition (SVD) for the rectangular matrix. The SVD for a $M\times N~{}(M>N)$ rectangular matrix $A$ is denoted as $A=U\Sigma V^{\dagger}$, where $U,V$ are unitary matrices and $\Sigma$ is a diagonal matrix for the singular values. Then the generalized inverse matrix is defined as $A^{-1}=V\Sigma^{-1}U^{\dagger}$, where $\Sigma$ and $\Sigma^{-1}$ are given by

[TABLE]

where \mbox{\boldmath0}_{M\times N} represents the $M\times N$ zero matrix. We here assume that all singular values are non-zero, otherwise the generalized inverse matrix cannot be defined.

An extension to higher angular momentum components than $l=4$ is straightforward. Including the $l=6$ component, for instance, we can extract the radial functions for $l=0,4,$ and $6$ by solving the equation

[TABLE]

Since the matrix in the right-hands side is a $N\times 3$ rectangular matrix, we need at least 3 points which satisfy $|\vec{x}_{i}|=r$ but are not transformed by the cubic rotation.

We next consider the extraction of the Laplacian for the radial function such as $g_{00}(r)$ in order to construct the potentials. Neglecting components with $l\geq 6$ again, the Laplacian of the NBS wave function in Eq. (51) becomes

[TABLE]

where the second term is evaluated in the continuum relation as

[TABLE]

Combining Eqs. (51) with (76), we have

[TABLE]

which can be solved by SVD. Note that the Laplacian in the left-hand-side of Eq. (76) is approximated by the second order difference, which has $\mathcal{O}(a^{2})$ discretized errors.

An advantage of this method over Misner’s method is that we need no parameter to extract the $l=0$ component of the NBS wave function. Applying this method, we extract the $l=0$ component of the $A_{1}^{+}$ projected $R$-correlator for the spin-singlet $\Lambda_{c}N$ system. In Fig. 9, $g_{00}(r)$ and its Laplacian are compared with those extracted by Misner’s method. The spherical harmonic amplitude $g_{00}(r)$ extracted by this method does not show comb-like structures and is consistent with $g_{00}(r)$ in Misner’s method, while its Laplacian has comb-like structures, which probably originate from discretized errors due to the Laplacian in the left-hand-side of Eq. (76). Consequently, the potential constructed from the $l=0$ component in this method, shown in Fig. 10, also has comb-like structures. Therefore Misner’s method works better for the HAL QCD potential than the method in this section, which however may be used to extract $g_{00}(r)$ only.222On this point, the method in this section essentially corresponds to Misner’s method with the thickness of the spherical shell $\Delta\rightarrow 0$. With the generalization of $\Delta\neq 0$, Misner’s method can utilize more data points and thus obtain $g_{00}(r)$ at more points of $r$ (red) than the method in this section (blue) at short distances.