Scattering phase shifts and mixing angles for an arbitrary number of coupled channels on the lattice

TL;DR

This paper introduces a lattice method to determine scattering phase shifts and mixing angles for multiple coupled channels, extending previous techniques limited to two channels, applicable to particles with any spin.

Contribution

The paper presents a novel lattice approach for extracting scattering parameters in systems with arbitrary coupled channels, surpassing previous two-channel limitations.

Findings

Method accurately reproduces continuum phase shifts and mixing angles.

Applicable to particles with any spin, demonstrated with spin-1 bosons.

Shows agreement with Schrödinger equation solutions for test potentials.

Abstract

We present a lattice method for determining scattering phase shifts and mixing angles for the case of an arbitrary number of coupled channels. Previous nuclear lattice effective field theory simulations were restricted to mixing of up to two partial waves for scattering of two spin-$1/2$ particles, which is insufficient for analyzing nucleon-nucleus or nucleus-nucleus scattering processes. In the proposed method, the phase shifts and mixing angles are extracted from the radial wave functions obtained by projecting the three-dimensional lattice Hamiltonian onto the partial wave basis. We use a spherical wall potential as a boundary condition along with a channel-mixing auxiliary potential to construct the full-rank $S$ matrix. Our method can be applied to particles with any spin, but we focus here on scattering of two spin-$1$ bosons involving up to four coupled channels. For a considered test potential, the phase shifts and mixing angles extracted on the lattice are shown to agree with the ones calculated by solving the Schr\"odinger equation in the continuum.

Figures

Tables

Table 1. Table 1: Parameters for the lattice calculation depending on the number of coupled scattering channels ( n ch subscript 𝑛 ch n_{\mathrm{ch}} ): lattice length L 𝐿 L , interval [ R I , R O ] subscript 𝑅 𝐼 subscript 𝑅 𝑂 [R_{I},R_{O}] for fitting wave functions, spherical wall radius R W subscript 𝑅 𝑊 R_{W} , coefficient V 0 subscript 𝑉 0 V_{0} of Gaussian auxiliary potential, coefficient U 0 subscript 𝑈 0 U_{0} of mixing potential and number of computed eigenvectors n eig subscript 𝑛 eig n_{\mathrm{eig}} . The lattice spacing has been chosen as a = 1.9733 ​ fm 𝑎 1.9733 fm a=1.9733\;\mathrm{fm} .

$n_{\mathrm{ch}}$	$L$ (units of $a$)	$R_I$ (units of $a$)	$R_O$ (units of $a$)	$R_W$ (units of $a$)	$V_0$ (MeV)	$U_0$ (MeV)	$n_{\mathrm{eig}}$
1	35	9.02	12.02	15.02	0	$-$	10
	41	9.02	12.02	18.02	0	$-$	10
	47	9.02	12.02	21.02	0	$-$	10
2	35	9.02	12.02	15.02	0	20	15
	41	9.02	12.02	18.02	0	20	15
3	35	9.02	12.02	15.02	0	10	70
4	35	9.02	12.02	15.02	0	5	110

Equations

$$r_1,r_2,r_3=0,\dots ,L-1,$$

$$|\mathbf{r}⟩=|\mathbf{r}+L{\hat{\mathbf{e}}}_1⟩=|\mathbf{r}+L{\hat{\mathbf{e}}}_2⟩=|\mathbf{r}+L{\hat{\mathbf{e}}}_3⟩$$

$$H_0|\mathbf{r}⟩$$

$$+\frac{3}{40\mu }\sum _{i=1}^3\left(|\mathbf{r}+2{\hat{\mathbf{e}}}_i⟩+|\mathbf{r}-2{\hat{\mathbf{e}}}_i⟩\right)$$

$$-\frac{1}{180\mu }\sum _{i=1}^3\left(|\mathbf{r}+3{\hat{\mathbf{e}}}_i⟩+|\mathbf{r}-3{\hat{\mathbf{e}}}_i⟩\right)$$

$$V_{\mathrm{wall}}(\mathbf{r})=\Lambda \theta (r-R_W),$$

$$V_{\mathrm{aux}}(\mathbf{r})=V_0\exp [-(r-R_W{)}^2]\theta (R_W-r)$$

$$H=H_0+V+V_{\mathrm{wall}}+V_{\mathrm{aux}}.$$

$${|R⟩}_{s,l,j}$$

$$\times Y_{l,l_z}(\hat{\mathbf{r}}){\delta }_{r,R}|\mathbf{r}⟩\otimes |s_{1,z},s_{2,z}⟩,$$

$${|R⟩}_{\alpha }:={|R⟩}_{s_{\alpha },l_{\alpha },j_{\alpha }}\text{ for }\alpha =1,\dots ,n,$$

$$[H_R(R_1,R_2){]}_{\alpha \beta }$$

$$\times \text{tensor}{[}_{{\alpha }^{\prime }}]⟨R_1|H|R_2⟩{}_{{\beta }^{\prime }},$$

$$[N(R){]}_{\alpha {\alpha }^{\prime }}=\text{tensor}{[}_{\alpha }]\text{innerproduct}RR{}_{{\alpha }^{\prime }}.$$

$$\psi (r)=[{\psi }_1(r),\dots ,{\psi }_n(r){]}^T.$$

$$\psi (r)=A{h}_l^{-}(pr)+B{h}_l^{+}(pr),$$

$$S=B/A=e^{2i{\delta }_l}.$$

$$E(\mathbf{p})$$

$$-\frac{1}{90\mu }\sum _{i=1}^3\cos (3p_i).$$

$$\mathbf{p}=(p\sin \theta \cos \varphi ,p\sin \theta \sin \varphi ,p\cos \theta ).$$

$$E_{s,l,j}(p)=\int d{\Omega }_p\sum _{l_z,s_z}\sum _{l_z^{\prime },s_z^{\prime }}\sum _{s_{1,z}}\sum _{s_{2,z}}{C}_{0,l_z,s_z}^{j,l,s}{C}_{s_z,s_{1,z},s_{2,z}}^{s,s_1,s_2}$$

$$\times C_{0,l_z^{\prime },s_z^{\prime }}^{j,l,s}{C}_{s_z^{\prime },s_{1,z},s_{2,z}}^{s,s_1,s_2}{Y}_{l,l_z}^{\ast }(\hat{\mathbf{p}})Y_{l,l_z^{\prime }}(\hat{\mathbf{p}})E(\mathbf{p}).$$

$$\displaystyle S=\left(\begin{array}[]{ll}v_{1}^{+}&v_{2}^{+}\end{array}\right)\left(\begin{array}[]{ll}v_{1}^{-}&v_{2}^{-}\end{array}\right)^{-1},$$

$$\displaystyle\psi(r)=\left(\begin{array}[]{l}A_{1}h_{l_{1}}^{-}(pr)+B_{1}h_{l_{1}}^{+}(pr)\\ A_{2}h_{l_{1}}^{-}(pr)+B_{2}h_{l_{2}}^{+}(pr)\end{array}\right)$$

$$\displaystyle\psi^{*}(r)=\left(\begin{array}[]{l}A_{1}^{*}h_{l_{1}}^{+}(pr)+B_{1}^{*}h_{l_{1}}^{-}(pr)\\ A_{2}^{*}h_{l_{1}}^{+}(pr)+B_{2}^{*}h_{l_{2}}^{-}(pr)\end{array}\right)$$

$$v_1^{-}=(A_1,A_2{)}^T,$$

$$v_1^{+}=(B_1,B_2{)}^T,$$

$${\psi }^{\ast }=T\psi =\psi \Rightarrow v_1^{\pm }=v_2^{\pm },$$

$$(H_R+U)\psi (r)=E\psi (r)$$

$$\displaystyle U(r)=U_{0}\delta_{r,R_{M}}\left(\begin{array}[]{rr}0&i\\ -i&0\end{array}\right)$$

$$S$$

$$\displaystyle\qquad\times\left(\begin{array}[]{rr}\cos\epsilon&\sin\epsilon\\ -\sin\epsilon&\cos\epsilon\end{array}\right)$$

$$\psi (r)={\left[{\psi }_1(r),{\psi }_2(r)\right]}^T$$

$${\psi }^{\prime }(r)={\left[\mathfrak{R}{\psi }_1(r),\mathfrak{I}{\psi }_1(r),\mathfrak{R}{\psi }_2(r),\mathfrak{I}{\psi }_2(r)\right]}^T.$$

Full text

I Introduction

Lattice simulations provide a powerful computational approach to systems of strongly interacting particles, which is widely used in condensed matter, nuclear, and particle physics. In particular, lattice gauge theory is the only known numerical method that allows one to directly solve QCD in the nonperturbative domain. Here, remarkable progress has been achieved in the recent decades due to the rapid increase of computational power and algorithmic efficiency. In particular, high-precision lattice QCD calculations of hadronic observables, such as the masses and decay constants, are already available for physical values of the quark masses Tanabashi:2018oca . While hadronic reactions and resonance properties can also be addressed in lattice QCD, such calculations appear to be much more challenging and require developing reliable methods for relating the scattering amplitude to discrete finite-volume spectra accessible in lattice simulations, see Refs. Lage:2009zv ; Bernard:2010fp ; Hansen:2012tf ; Briceno:2012rv ; Hammer:2017uqm ; Hammer:2017kms ; Briceno:2017tce for recent work along this line and Briceno:2017max for a review article.

Lattice methods have also proven to be very efficient in describing low-energy nuclear systems in the framework of chiral effective field theory (EFT). Recently, the chiral expansion of the nucleon-nucleon (NN) potential has been pushed to fifth order (N4LO) Epelbaum:2014efa ; Epelbaum:2014sza ; Entem:2014msa ; Entem:2017gor within the continuum formulation. The NN potentials derived in chiral EFT in Ref. Reinert:2017usi allow, for the first time, for a nearly perfect description of the neutron-proton and proton-proton scattering data below the pion production threshold, which is comparable to or even better than that based on the available phenomenological potentials. Three- and four-nucleon forces have been worked out completely up to fourth order (N3LO) of the chiral expansion Ishikawa:2007zz ; Bernard:2007sp ; Bernard:2011zr ; Epelbaum:2007us , see also Refs. Krebs:2012yv ; Krebs:2013kha ; Epelbaum:2014sea ; Girlanda:2011fh for the derivation of selected contributions at N4LO and Refs. Epelbaum:2008ga ; Epelbaum:2012vx ; Machleidt:2011zz for review articles. To apply the interactions derived in chiral EFT to few- and many-nucleon systems, it is necessary to solve the quantum mechanical $A$-body problem, which can be achieved using continuum ab initio methods including Faddeev-Yakubovsky equations Gloeckle:1995jg , the no-core configuration interaction approach Barrett:2013nh , coupled-cluster expansions Hagen:2012fb , the in-medium similarity renormalization group approach Hergert:2012nb , self-consistent Green’s function Soma:2012zd , or quantum Monte Carlo methods Lovato:2013cua . Alternatively, a discretized version of chiral EFT Lee:2008fa ; Lee:2016fhn ; UGMlectures has been successfully applied to a broad range of nuclear systems. This approach has an appealing feature of being well suited for dealing with strongly clustered systems such as the famous Hoyle state in 12C Epelbaum:2011md ; Epelbaum:2012qn ; Epelbaum:2013paa and some of the low-lying states of 16O, which often represent a challenge for continuum methods. See Ref. Freer:2017gip for a recent review on clustering in light nuclei. So far, nuclear lattice simulations have been carried out for light- and medium-mass nuclei and neutron matter up to third order in the chiral expansion Borasoy:2007vk ; Epelbaum:2009zsa ; Lahde:2013uqa . For a recent lattice EFT study of NN scattering at N3LO see Ref. Li:2018ymw . This method was also employed to study the dependence of the triple-$\alpha$ process on the fundamental constants of nature Epelbaum:2012iu ; Epelbaum:2013wla , see Ref. Meissner:2014pma for a related discussion, to investigate the isotopic dependence of nuclear clustering Elhatisari:2017eno and to determine the features of the nuclear force essential for nuclear binding Lu:2018bat . It is important to emphasize that the development of chiral EFT interactions is more difficult on the lattice than in the continuum as it requires establishing efficient techniques for extracting the scattering amplitude from the finite-volume discrete spectra and for dealing with the breaking of rotational Alarcon:2017zcv ; Klein:2018iqa and Galilean invariance Li:2019ldq due to nonzero lattice spacing. Lüscher’s finite-volume method is one possible approach to compute scattering phase shifts on the lattice. For lattice QCD applications, there have been many recent advances on coupled-channel calculations and partial-wave mixing using Lüscher’s formalism Briceno:2014oea ; Moir:2016srx ; Briceno:2017qmb ; Woss:2018irj ; Woss:2019hse . However, for lattice EFT calculations of heavier nuclear systems, the large nuclear binding energies and very small finite-volume scattering energies make it difficult to implement Lüscher’s method with accuracy. For this reason, a more robust approach based on the spherical wall boundary conditions Carlson:1984zz was used in Refs. Borasoy:2007vy ; Lu:2015riz ; Alarcon:2017zcv ; Li:2018ymw . This technique is not only applicable to calculations of NN phase shifts on the lattice, but can also be combined with the adiabatic projection method Rokash:2015hra ; Elhatisari:2016hby , which allows one to access nuclear reactions via lattice simulations, see Ref. Elhatisari:2015iga for the first ab initio study of $\alpha$-$\alpha$ scattering. However, the spherical wall method has so far only been applied to uncoupled partial waves and the cases of two coupled channels, which is insufficient for studying nuclear reactions. The purpose of this paper is to generalize this technique to an arbitrary number of coupled channels.

Our paper is organized as follows. In Sec. II, we introduce the lattice notation, review the method of Ref. Lu:2015riz to compute the scattering parameters with up to two coupled channels, and extend this approach to scattering of particles or nuclear clusters of an arbitrary high spin. As an application, we consider in Sec. III the scattering problem of two spin-$1$ bosons using a toy-model potential, which is similar to the one from Refs. Borasoy:2007vy ; Lu:2015riz . The main results of our study are summarized in Sec. IV.

II Scattering of two particles with arbitrary spin on the lattice

II.1 Calculational setup

We employ a periodic cubic lattice with the length $L$ and spacing $a$, and define orthonormal lattice states $\ket{\mathbf{r}}$ with

[TABLE]

[TABLE]

due to the periodic boundary condition. All quantities in this section are given in dimensionless lattice units, i.e., they must be multiplied by an appropriate power of the lattice spacing $a$ to obtain their physical values. In the following, we briefly review the method to compute the scattering phase shifts and mixing angles for up to two coupled channels introduced in Ref. Lu:2015riz , which will then be generalized to the case of three or more coupled channels.

We consider the general scattering problem of two particles with spins $s_{1},s_{2}$ and masses $m_{1},m_{2}$ interacting with the potential $V(\mathbf{r})$. The free Hamiltonian in the center-of-mass (c.m.) system is discretized as Borasoy:2007vy

[TABLE]

using the reduced mass $\mu=m_{1}m_{2}/(m_{1}+m_{2})$. The above expression corresponds to the $O(a^{4})$-improved free lattice Hamiltonian. To avoid artifacts induced by the periodic boundary conditions, it is convenient to use a spherical wall boundary condition by adding the potential

[TABLE]

where $\theta$ is the Heaviside function, $R_{W}$ is the wall radius, and $\Lambda$ is a large positive constant111Following Ref. Lu:2015riz , we use the value $\Lambda=10^{6}$ (given in dimensionless lattice units) in the numerical calculations. Borasoy:2007vy . Calculating the scattering parameters at low momenta usually requires large lattices, which makes the analysis computationally expensive. It is more convenient to introduce an auxiliary potential outside of the range of $V$, which can be chosen, e.g., of a Gaussian type

[TABLE]

with $V_{0}\leq 0$ in order to control the eigenenergies of the Hamiltonian Lu:2015riz . The complete Hamiltonian including all contributions is then given by

[TABLE]

II.2 Projection onto partial waves

The three-dimensional problem can be reduced to the one-dimensional one by defining radial states for a partial wave ${}^{2s+1}l_{j}$,

[TABLE]

i.e., the lattice sites $\mathbf{r}$ with the same radial distance $R$ are grouped together according to the irreducible representations of the rotational group. Here, $C^{j,l,s}_{j_{z},l_{z},s_{z}}$ and $C^{s,s_{1},s_{2}}_{s_{z},s_{1,z},s_{2,z}}$ are the Clebsch-Gordan coefficients for the spin-orbit and spin-spin couplings, respectively. The spherical harmonics $Y_{l,l_{z}}$ behave like $Y_{l,l_{z}}(\mathbf{0})=\delta_{l,0}/\sqrt{4\pi}$ at the origin. Since the results obtained here do not depend on $j_{z}$ in the continuum limit, we can choose $j_{z}=0$.

Note that the radial states have to be normalized by dividing them by the square root of their norm. States that are not linearly independent or have vanishing norm must be omitted to make the norm matrix invertible. Afterwards, the Hamiltonian can be projected onto the normalized radial states. For $n$ coupled channels with

[TABLE]

one has

[TABLE]

where $N^{-1/2}$ is the inverse square root of the norm matrix

[TABLE]

Multiplying the eigenvectors of the projected Hamiltonian $H_{R}$ by $N^{-1/2}$ from the left yields the radial wave functions

[TABLE]

II.3 Single-channel case

Outside of the range of the potential, the radial wave functions are linear combinations of the spherical Hankel functions $h_{l}^{\pm}(pr)$, where $p$ is the momentum in the center-of-mass system. For a single scattering channel, one obtains a wave function of the form

[TABLE]

which allows one to extract the phase shift $\delta_{l}$ from the $S$ matrix via

[TABLE]

The coefficients $A$, $B$ are computed by fitting the spherical Hankel functions to the wave function in an interval $[R_{I},R_{O}]$ outside of the range of the potential. The momentum $p$ is determined from the eigenenergy of the Hamiltonian using the lattice dispersion relation

[TABLE]

This equation can be expressed in spherical coordinates with

[TABLE]

In order to remove the angular dependence, the dispersion relation must be projected onto partial waves as well:

[TABLE]

The angular integration can be facilitated by Taylor-expanding $E(\mathbf{p})$ up to order $O(p^{\kappa})$, which should be sufficiently high to yield accurate results up to the cutoff momentum $\pi/a$. 222In numerical calculations, we use $\kappa=30$. Afterward, the c.m. system momentum can be computed by solving $E_{s,l,j}(p)$ for $p$.

II.4 Scattering with two coupled channels

For two coupled channels, the $S$ matrix must be constructed as

[TABLE]

where $v_{1}^{\pm},v_{2}^{\pm}$ are linearly independent two-component vectors containing the coefficients in front of the spherical Hankel functions $h_{l_{i}}^{\pm}$. A simple way to obtain these coefficients would be to extract them from a complex wave function

[TABLE]

and its complex conjugate

[TABLE]

using

[TABLE]

However, the Hamiltonian $H_{R}$ commutes with the time-reversal operator $T$ so that it holds

[TABLE]

i.e., the vectors $v_{1}^{\pm},v_{2}^{\pm}$ are linearly dependent. Thus, one runs into the problem of having only one independent solution per lattice energy. In order to circumvent this problem, an auxiliary imaginary potential term that breaks the time-reversal symmetry can be added to the Hamiltonian:

[TABLE]

with

[TABLE]

and $U_{0}\in\mathbb{R}$. The radius $R_{M}$ should lie outside the range of the test potential and can be chosen close to the spherical wall radius, $R_{M}\lesssim R_{W}$. Because the matrix in Eq. (30) mixes the two channels, $U(r)$ will be referred to as the mixing potential. Finally, the $S$ matrix is decomposed according to the Blatt-Biedenharn parametrization Blatt:1952zza

[TABLE]

with the phase shifts $\delta_{1},\delta_{2}$ and the mixing angle $\epsilon$. Since the lattice dispersion relation can yield slightly different momenta $p_{1}$, $p_{2}$ for the two channels at the same energy, we assume that the phase shift $\delta_{\alpha}$ is measured at momentum $p_{\alpha}$ and that the mixing angle is measured at the average momentum $(p_{1}+p_{2})/2$.

II.5 Scattering with an arbitrary number of coupled channels

If $n>2$ coupled channels must be considered, one needs $n$ linearly independent wave functions in each channel. The complex conjugation is not sufficient for this purpose because it can only generate two independent solutions $\psi$ and $\psi^{*}$. In order to find an alternative approach, we first consider the case of two coupled partial waves again. The two-channel wave function

[TABLE]

can be rewritten as

[TABLE]

To reproduce Eqs. (27) and (30), the radial Hamiltonian and the mixing potential must be modified accordingly:

[TABLE]

On the other hand, instead of using Eq. (40), we can regard the wave function vector as having four independent complex components

[TABLE]

A natural extension of Eqs. (50), (51) to three coupled scattering channels is given by introducing

[TABLE]

[TABLE]

A generalization to $n$ channels is straightforward:

[TABLE]

for $\alpha,\alpha^{\prime},\beta,\beta^{\prime}=1,\dots,n$. The wave function vector has the form

[TABLE]

where $\psi_{\beta+(\alpha-1)n}^{\prime}$ denotes the $\beta$th wave function for the $\alpha$th scattering channel with $\alpha,\beta=1,\dots,n$.

More generally, any Hermitian matrix that produces $n$ linearly independent solutions in every channel can be used to define the mixing potential, i.e., it must hold ${U^{\prime}}^{\dagger}=U^{\prime}$ and the matrix $M$ with

[TABLE]

must have rank $n$. The particular choice for the mixing potential in Eq. (71) is consistent with the one employed for two channels in Ref. Lu:2015riz .

Each component of the wave function vector has the form

[TABLE]

with $\alpha,\beta=1,\dots,n$. Since it holds

[TABLE]

for $\beta=1,\dots,n$, one can construct the $S$ matrix as

[TABLE]

The Blatt-Biedenharn parametrization also has to be extended to $n>2$ coupled channels Blatt:1952zza :

[TABLE]

where $O$ is a real orthogonal matrix. (This decomposition is equivalent to computing the eigenvalues and eigenvectors of $S$.) Again, the phase shift $\delta_{\alpha}$ is assigned to the momentum $p_{\alpha}$ in scattering channel $\alpha$. For simplicity, we define the mixing angles as

[TABLE]

for $\alpha,\beta=1,\dots,n$ and $\beta>\alpha$ because a real orthogonal $n\times n$ matrix can be given by $n(n-1)/2$ real parameters.333Note that Eq. (90) allows one to extract the mixing angles from the matrix $O$, which is sufficient for the purpose of this paper. The equation can, however, not be used to reconstruct the $S$ matrix from the phase shifts and mixing angles in a unique way. This problem could be avoided by parametrizing the orthogonal matrix $O$ in terms of generators of the rotation group.

III Test case: scattering of two spin-$1$ particles

The method described in Sec. II allows one to determine scattering phase shifts and mixing angles on the lattice for an arbitrary number of coupled channels and for any type of particles. As a concrete example, we consider the scattering problem of two spin-$1$ bosons having nearly the same mass as the deuteron, $m_{1,2}=2m_{N}=2\times 938.92\;\mathrm{MeV}$. As a test potential, we employ the corresponding generalization of the toy-model potential used for two spin-1/2 fermions in Refs. Borasoy:2007vy ; Lu:2015riz :

[TABLE]

where the spin-dependent part is given by

[TABLE]

Here, $\mathbf{s}_{1}$ and $\mathbf{s}_{2}$ denote the spin matrices for the considered particles. The constants are set to $C=-2\;\mathrm{MeV}$ and $r_{0}=0.02\;\mathrm{MeV}^{-1}$, and the lattice spacing is chosen to be $a=(100\;\mathrm{MeV})^{-1}=1.9733\;\mathrm{fm}$. Notice that by projecting the test potential onto partial waves one obtains up to four coupled scattering channels. We calculate the phase shifts and mixing angles for the following cases:

• (i)

Uncoupled channels: ${}^{3}P_{0}$, ${}^{3}P_{1}$, ${}^{3}D_{2}$, ${}^{3}F_{3}$, ${}^{3}G_{4}$, ${}^{3}H_{5}$, ${}^{5}D_{1}$;

• (ii)

Two coupled channels: ${}^{3}SD_{1}$, ${}^{3}PF_{2}$, ${}^{3}DG_{3}$, ${}^{3}FH_{4}$, ${}^{1}S_{0}/{}^{5}D_{0}$, ${}^{5}PF_{2}$, ${}^{5}DG_{3}$, ${}^{5}FH_{4}$;

• (iii)

Three coupled channels: ${}^{1}P_{1}/{}^{5}PF_{1}$;

• (iv)

Four coupled channels: ${}^{1}D_{2}/{}^{5}SDG_{2}$.

We further emphasize that the considered potential possesses four bound states in the ${}^{3}SD_{1}$, ${}^{1}S_{0}/{}^{5}D_{0}$, ${}^{5}PF_{2}$ and ${}^{1}D_{2}/{}^{5}SDG_{2}$ channels with the binding energies of $-0.258\;\mathrm{MeV}$, $-0.204\;\mathrm{MeV}$, $-0.198\;\mathrm{MeV}$, and $-0.583\;\mathrm{MeV}$, respectively, as visualized in Fig. 1(a). The binding energies have been obtained by computing the eigenvalues of the momentum-space Hamiltonian in the infinite-volume continuum.

On the lattice, the eigenvectors of the Hamiltonian corresponding to the lowest positive eigenenergies have been used as radial wave functions. The parameters used in our calculations, which have partly been adopted from Ref. Lu:2015riz , can be found in Table 1 while Fig. 1(b) visualizes the different contributions to the potential on the lattice. Notice that the Gaussian auxiliary potential can distort the wave function and may even generate additional bound states, see Fig. 2. This may result in the appearance of the outlying points in the calculated phase shifts or mixing angles as will be discussed below. For this reason, instead of imposing the auxiliary potential, we have actually varied the lattice size to generate data at low momenta.

The choice of the mixing potential requires some care, too. In particular, as shown in Table 1, the numerical value of the strength $U_{0}$ has been decreased for three and four channels in order to avoid the appearance of additional bound states as visualized in Fig. 3(a). On the other hand, choosing too small in magnitude values of $U_{0}$ leads to numerical instabilities causing a violation of unitarity in the calculated $S$ matrix, see Fig. 3(b). For large or small values of the coefficient $U_{0}$, outlying data points may appear in the plots of the phase shifts and mixing angles, as illustrated in Fig. 4. However, the outlying points for small coefficients are only caused by numerical round-off errors. We have verified that they can thus be removed using numbers with higher precision. Therefore, any value of $U_{0}$ should be suitable as long as no additional bound states are produced by the mixing potential.

In Figs. 5$-$9, we show the phase shifts and mixing angles for all considered scattering channels, which have been calculated on the lattice using the method presented above. To benchmark our calculations, we have also computed the scattering parameters in the continuum. This has been achieved by solving the radial Schrödinger equation for a fixed energy $E=p^{2}/(2\mu)$. The boundary conditions at $r\simeq 0$ have to be chosen in such a way that one obtains a sufficient number of linearly independent solutions. Finally, the $S$ matrix can again be extracted from the wave function by fitting spherical Hankel functions. As shown in Figs. 5$-$9, the results of the lattice calculations are in essentially a perfect agreement with the ones calculated using the continuum approach in the considered range of c.m. system momenta up to $p=120\;\mathrm{MeV}$. The cutoff momentum associated with the lattice spacing of $a=(100\;\mathrm{MeV})^{-1}=1.9733\;\mathrm{fm}$ employed in our analysis is $\Lambda_{\mathrm{latt}}\sim\pi/a\simeq 314$ MeV. Thus, the lattice and continuum results are expected to agree for momenta well below $\Lambda_{\mathrm{latt}}$. Indeed, for c.m. system momenta higher than $p=120\;\mathrm{MeV}$, the deviations between the continuum and lattice results start to become visible. This observation is in line with the findings of Refs. Borasoy:2007vy ; Lu:2015riz .

We also notice a subtlety in the extraction of phase shifts in multichannel cases due to the fact that the eigenvalues $e^{2i\delta_{1}},\dots,e^{2i\delta_{n}}$ of the $S$ matrix have no predefined ordering Blatt:1952zza . Therefore, the multichannel continuum calculation has been repeated without the off-diagonal elements in the potential matrix, see the gray dashed lines in Figs. 6$-$9. Then, the phase shifts in the coupled channels have been ordered such that they are roughly consistent with the phase shifts obtained without the coupling. The large differences between the solid and dashed lines demonstrate the very important role of channel mixing in the considered toy model. For many-body systems where no continuum calculation is possible, the comparison can be performed with lattice data instead. If the results at very low momenta are available, one may also possibly identify the partial waves from the threshold behavior of the eigenphases. Last but not least, we emphasize that the behavior the phase shifts in the ${}^{3}SD_{1}$, ${}^{1}S_{0}/{}^{5}D_{0}$, ${}^{5}PF_{2}$, and ${}^{1}D_{2}/{}^{5}SDG_{2}$ channels with $\delta(p=0)=\pi$ is consistent with the appearance of a single bound state in each of these channels, see Fig. 1(a), in agreement with Levinson’s theorem.

IV Summary and outlook

In this paper, we considered two-particle scattering by solving the Schrödinger equation on the lattice. A generalization of the method used in Ref. Lu:2015riz for spin-$1/2$ fermions to scattering of particles with any spin and an arbitrary number of coupled scattering channels has been proposed. For the case of two spin-$1$ bosons, the proposed method was benchmarked against the continuum approach and demonstrated to yield accurate and reliable results for phase shifts and mixing angles for momenta well below the lattice cutoff. Our study opens the way to perform ab initio chiral EFT calculations in the four-nucleon continuum and to access nuclear reactions on the lattice using the adiabatic projection method. Work along these lines is in progress.

Acknowledgements.

We are grateful to Ning Li and Bing-Nan Lu for helpful discussions and to Ulf-G. Meißner and Xiu-Lei Ren for useful comments on the manuscript. We also thank Ning Li for sharing his nucleon-nucleon scattering code. This work was supported by DFG (SFB/TR 110, “Symmetries and the Emergence of Structure in QCD”), the BMBF (Grant No. 05P2015) and the U.S. Department of Energy (Grants No. DE-SC0018638 and No. DE-AC52-06NA25396).