Scattering processes and resonances from lattice QCD

TL;DR

This paper reviews how lattice QCD techniques are used to study hadron resonances, providing insights into their properties and the underlying strong interaction dynamics in a non-perturbative regime.

Contribution

It presents the leading lattice QCD approach for extracting scattering amplitudes and resonance properties from finite-volume spectra, highlighting recent progress and challenges.

Findings

Lattice QCD can determine resonance masses and widths.

Finite-volume spectra encode scattering information.

Progress addresses current methodological challenges.

Abstract

The vast majority of hadrons observed in nature are not stable under the strong interaction, rather they are resonances whose existence is deduced from enhancements in the energy dependence of scattering amplitudes. The study of hadron resonances offers a window into the workings of quantum chromodynamics (QCD) in the low-energy non-perturbative region, and in addition, many probes of the limits of the electroweak sector of the Standard Model consider processes which feature hadron resonances. From a theoretical standpoint, this is a challenging field: the same dynamics that binds quarks and gluons into hadron resonances also controls their decay into lighter hadrons, so a complete approach to QCD is required. Presently, lattice QCD is the only available tool that provides the required non-perturbative evaluation of hadron observables. In this article, we review progress in the study of few-hadron reactions in which resonances and bound-states appear using lattice QCD techniques. We describe the leading approach which takes advantage of the periodic finite spatial volume used in lattice QCD calculations to extract scattering amplitudes from the discrete spectrum of QCD eigenstates in a box. We explain how from explicit lattice QCD calculations, one can rigorously garner information about a variety of resonance properties, including their masses, widths, decay couplings, and form factors. The challenges which currently limit the field are discussed along with the steps being taken to resolve them.

Tables

Table 1. Table 1: The lowest partial waves, ℓ ℓ \ell , that are subduced into the Λ Λ \Lambda irrep for the cubic group for bosonic systems – the dimension of each irrep is given in parentheses. An irrep of dimension d 𝑑 d has d 𝑑 d equivalent “rows” – analogous to the 2 ​ ℓ + 1 2 ℓ 1 2\ell+1 m 𝑚 m -values in the rotationally invariant case.

$\Lambda$(dim)	$\ell$
$A_1\left(1\right)$	$0,4,\dots$
$T_1\left(3\right)$	$1,3,4,\dots$
$T_2\left(3\right)$	$2,3,4,\dots$
$E\left(2\right)$	$2,4,\dots$
$A_2\left(1\right)$	$3,\dots$

Equations

$$\sqrt{s}=\sqrt{E^2-{\mathbf{P}}^2}=E^{\star }=\sqrt{m_1^2+{\mathbf{q}}^{\star 2}}+\sqrt{m_2^2+{\mathbf{q}}^{\star 2}},$$

$$q^{\star}=\frac{1}{2}\left[s-2(m_{1}^{2}+m_{2}^{2})+\frac{\big{(}m_{2}^{2}-m_{1}^{2}\big{)}^{2}}{s}\right]^{1/2},$$

$$\mathrm{Im}\frac{1}{{\mathcal{M}}_{\ell }}=-\frac{1}{16\pi }\frac{2q^{\star }}{E^{\star }}\Theta (E^{\star }-E^{thr.}),$$

$$\mathrm{Re}\frac{1}{{\mathcal{M}}_{\ell }}=\frac{1}{16\pi }\frac{2q^{\star }}{E^{\star }}\cot {\delta }_{\ell }(E^{\star }).$$

$$\mathrm{Im}{\mathcal{M}}_{\ell }^{\ast }\propto |{\mathcal{M}}_{\ell }{|}^2,$$

$$\mathrm{Im}{\left({\mathcal{M}}_{\ell }^{-1}(E^{\star })\right)}_{ab}=-{\delta }_{ab}\frac{1}{16\pi }\frac{2q_a^{\star }}{E^{\star }}\Theta (E^{\star }-E_a^{thr.}),$$

$${\left({{\mathcal{M}}_{\ell }}^{-1}\right)}_{ab}={\left({{\mathcal{K}}_{\ell }}^{-1}\right)}_{ab}-i{\delta }_{ab}\frac{1}{16\pi }\frac{2q_a^{\star }}{E^{\star }},$$

$${\mathcal{Z}}_E=\int \mathcal{D}U\mathcal{D}\psi \mathcal{D}\bar{\psi }{e}^{-S_E(\psi ,\bar{\psi },U)},$$

$$\frac{1}{{\mathcal{Z}}_E}\int \mathcal{D}U\mathcal{D}\psi \mathcal{D}\bar{\psi }f(\psi ,\bar{\psi },U)e^{-S_E(\psi ,\bar{\psi },U)},$$

$${\mathcal{Z}}_E=\int \mathcal{D}U\det Q(U)e^{-S_g(U)}=\int \mathcal{D}U{e}^{-{\tilde{S}}_E(U)},$$

$$\det\big{[}F^{-1}(E,\mathbf{P};L)+\mathcal{M}(E)\big{]}=0.$$

$$p=\frac{2\pi }{L}n-\frac{2}{L}\delta (p).$$

$$C_L(x_4-y_4,\mathbf{P})\equiv {\int }_L$$

$$\displaystyle\times\Big{[}\langle 0|T\mathcal{A}(x)\mathcal{B}^{\dagger}(y)|0\rangle\Big{]}_{L}\,,$$

$$L^{6}\sum_{n}e^{-E_{n}(x_{4}-y_{4})}\big{\langle}0\big{|}\mathcal{A}(0)\big{|}E_{n},\mathbf{P};L\big{\rangle}\,\big{\langle}E_{n},\mathbf{P};L\big{|}\mathcal{B}^{\dagger}(0)\big{|}0\big{\rangle},$$

$$\int \frac{d\text{k}}{(2\pi {)}^3}⟶$$

$$=\sum _{\mathbf{k}}\frac{\Delta n_x\Delta n_y\Delta n_z}{L^3}=\frac{1}{L^3}\sum _{\mathbf{k}}.$$

$${\mathcal{F}}_L$$

$$\times \mathcal{L}(P-k,k)\Delta (k)\Delta (P-k){\mathcal{R}}^{†}(P-k,k),$$

$${\mathcal{F}}_L$$

$$\displaystyle\times\mathcal{L}(P-k,k)\frac{1}{E-\omega_{k}-\omega_{Pk}+i\epsilon}\mathcal{R}^{\dagger}(P-k,k)\bigg{|}_{k_{4}=i\omega_{k}}$$

$${\mathcal{L}}_{\mathrm{on}}(P,{\text{k}}^{\star })$$

$${\mathcal{R}}_{\mathrm{on}}^{†}(P,{\text{k}}^{\star })$$

$${\mathcal{F}}_L=-{\mathcal{L}}_{\ell m}(P)F_{\ell m;{\ell }^{\prime }{m}^{\prime }}(P,L){\mathcal{R}}_{{\ell }^{\prime }{m}^{\prime }}^{†}(P),$$

$$\displaystyle F_{\ell m;\ell^{\prime}m^{\prime}}(P,L)\equiv\bigg{[}\frac{1}{L^{3}}\sum_{\mathbf{k}}-\int\!\!\frac{d\textbf{k}}{(2\pi)^{3}}\bigg{]}\quad\frac{4\pi\,Y_{\ell m}(\hat{\textbf{k}}^{\star})\,Y^{*}_{\ell^{\prime}m^{\prime}}(\hat{\textbf{k}}^{\star})}{2\omega_{k}\,2\omega_{Pk}\,(E-\omega_{k}-\omega_{Pk}+i\epsilon)}\left(\frac{k^{\star}}{q^{\star}}\right)^{\ell+\ell^{\prime}}\,.$$

$$C_L(x_4-y_4,\text{P})$$

$$\times L^3{A}^{\star }(E_n,\text{P})\mathcal{R}(E_n,\text{P})B^{\star }(E_n,\text{P})$$

$$\mathcal{R}(E_n,\text{P})$$

$$\displaystyle\lim_{P_{4}\rightarrow iE_{n}}\left(-(iP_{4}+E_{n})\Big{[}F^{-1}(P,L)+\mathcal{M}(P)\Big{]}^{-1}\right)\,.$$

$$\det\big{[}F^{-1}(E,\mathbf{P};L)+\mathcal{M}(E)\big{]}=0,$$

$$\cot {\delta }_{\ell }(E)=-\cot {\varphi }_{\ell }^{\Lambda }(\mathbf{P},L,E)$$

$$\mathcal{M}=\left(\begin{array}[]{cc}\mathcal{M}_{\pi\pi;\pi\pi}&\mathcal{M}_{\pi\pi;K\overline{K}}\\ \mathcal{M}_{\pi\pi;K\overline{K}}&\mathcal{M}_{K\overline{K};K\overline{K}}\end{array}\right),$$

$$\det\left[\left(\begin{array}[]{cc}F^{-1}_{\pi\pi}&0\\ 0&F^{-1}_{K\overline{K}}\end{array}\right)+\left(\begin{array}[]{cc}\mathcal{M}_{\pi\pi;\pi\pi}&\mathcal{M}_{\pi\pi;K\overline{K}}\\ \mathcal{M}_{\pi\pi;K\overline{K}}&\mathcal{M}_{K\overline{K};K\overline{K}}\end{array}\right)\right]=0,$$

Full text

1. I Introduction

2. II Resonances, composite particles, and scattering amplitudes

3. II.1 Pole singularities

4. II.2 Coupled-channel scattering

5. II.3 Diagrammatic representation

6. III Lattice QCD

7. IV Scattering in a finite-volume

8. IV.1 Scattering in non-relativistic quantum mechanics in one space dimension

9. IV.2 Scattering in a periodic cubic volume

10. IV.3 Relating scattering amplitudes to finite-volume spectra

11. IV.3.1 Dominance of the lowest partial-wave

12. IV.3.2 Coupled-channel scattering and parameterization of scattering amplitudes

13. IV.3.3 Examples of finite-volume spectra for simple scattering amplitudes

14. V Determining the finite-volume spectrum

15. V.1 Variational analysis of correlation matrices

16. V.2 Operator construction

17. V.2.1 The importance of “multi-hadron” operators

18. VI Examples of resonance determination

19. VI.1 Elastic resonances in $\pi\pi$ scattering

20. VI.2 Resonances in coupled-channel meson-meson scattering

21. VII Other approaches to resonance determination

22. VII.1 Resonances in the Lüscher formalism in the narrow-width approximation

23. VII.2 Resonances and ‘naive’ level counting

24. VII.3 Finite volume EFT Hamiltonian approach

25. VII.4 Unitarized chiral perturbation theory and chiral extrapolations

26. VII.5 Other approaches

27. VIII Coupling resonances to external currents

28. VIII.1 Determining matrix elements in lattice QCD

29. VIII.2 Lellouch-Lüscher formalism and its generalizations

30. VIII.3 Applications

31. IX Contemporary Extensions

32. IX.1 Particles with nonzero intrinsic spin

33. IX.2 Three-particle systems

34. IX.3 Elastic form factors of resonances

35. X Outlook

36. XI Numerical evaluation of $F$

I Introduction

In the last couple of decades, the field of hadron spectroscopy has undergone a renaissance, spurred in part by the observation in experiments of a large number of states whose presence had not been anticipated. Many of these, known colloquially as the “X,Y,Z” states, have been found in the heavy quark sector (see Liu (2014); Chen et al. (2016) for recent reviews on the topic), and a complete understanding of them remains elusive, with piecemeal proposals ranging from tetraquark constructions, to meson-meson molecules, to quark-gluon hybrid states etc… All these ‘pictures’ strongly depend on the model chosen to locally approximate the low-energy behavior of the fundamental theory of the strong interaction, quantum chromodynamics (QCD). Ultimately it is QCD which builds these states, and it is to QCD that we should turn to understand them.

In practice, one does not have to go high in the experimental hadron spectrum to find states that put our understanding of QCD to a test. An iconic example is the $\sigma/f_{0}(500)$ resonance, which couples to the $\pi\pi$ scattering channel, and is both extremely light and unnaturally broad, having a decay width larger than its mass. Until recently even its existence was not a settled issue Pelaez (2016); Garcia-Martin et al. (2011); Caprini et al. (2006), and at the QCD level, the origin of the $\sigma$ remains a mystery. On the other hand, there is also a spectrum of ‘conventional’, mainly narrow, resonances, the lightest of which is the $\rho$, whose mass pattern is fairly well described by models that consider mesons to be constructed from a quark and an antiquark, having a relatively weak coupling to their decay channels Shepherd et al. (2016).

Motivated in part by the excitement in contemporary experimental exploration of the spectrum of hadrons, there is demand for a parallel theoretical program that can shed light on these states. Models, either at the quark level, or the hadron level, can help us to gain insight into how these states are constructed, but eventually we need to anchor our understanding in first-principles calculation within QCD. Historically, this has been extremely challenging, with one relevant complication being that the vast majority of hadrons observed in nature are unstable under the strong interactions. We observe all but the lightest hadrons through their decay products as resonances, whose lifetimes are of the same order as the time scale of QCD interactions. This forces us to consider theoretical techniques which can consistently encode the physics of both binding and decay. Presently, the only rigorous theoretical tool which has been shown to meet these demands is lattice QCD.

Lattice QCD considers quark and gluon fields on a discrete grid of points of finite size, and by sampling possible configurations of these fields, with a probability dictated by the lagrangian of QCD, hadronic observables can be estimated, along with a measure of statistical uncertainty that can be reduced with increased computer time. Systematic errors arise from the choice of lattice spacing, the size of the box, and in many practical calculations, the values chosen for the quark masses, which for computational cost reasons may not be as low as the experimentally measured quark masses. These approximations are controlled, and the corresponding uncertainties can be systematically decreased, in principle.

The vast majority of lattice QCD calculations to date have focused on the properties of states that are stable within QCD, and these calculations have matured to the level where they can be considered realistic, with up, down, strange and charm quark masses tuned to their correct physical values, and even some including the relatively small effects due to QED Duncan et al. (1996); Blum et al. (2007); Borsanyi et al. (2015); Horsley et al. (2016); Blum et al. (2010); Aoki et al. (2012b); de Divitiis et al. (2013).

The study of hadron resonances using lattice QCD, and consequently the study of the QCD spectrum, is at a much earlier stage of its development, but rapid progress is being made that we will review here. An approach that relates the discrete spectrum of QCD eigenstates in the finite-volume defined by the lattice, to scattering amplitudes which may contain resonances, has been a powerful tool. Generically known as the “Lüscher method”, it allows us access to hadron scattering amplitudes computed from first principles in QCD. By applying techniques similar to those used in the analysis of experimental scattering data, we are able to infer the resonance content of these amplitudes by analytically continuing into the complex energy plane, where resonances appear as pole singularities.

From a practical point of view, a first challenge has been to develop the algorithmic and computational techniques to make possible the determination of a tower of excited states, and recent years have seen tremendous advances in lattice QCD methods, such as those described in Luscher (1986b); Luscher and Wolff (1990); Foley et al. (2005); Peardon et al. (2009); Blossier et al. (2009); Dudek et al. (2009); Morningstar et al. (2011); Dudek et al. (2010, 2013a); Liu et al. (2012); Dudek and Edwards (2012); Dudek (2011); Edwards et al. (2011); Dudek et al. (2011a). It is now not unusual to observe lattice calculations determining as many as two dozen states in a single quantum number channel, with many methods being sufficiently flexible as to place almost no restriction on this number in future calculations.

The detailed discrete spectra of states extracted in lattice calculations have been used to constrain two-body scattering amplitudes in both elastic and coupled-channel cases, using a rigorous formalism that has been made rather general in derivations presented over the past few years. In this review we will consider all aspects of resonance physics within lattice QCD, including the coupling of resonances to external currents, and will outline challenges which still remain and possible avenues to overcome them.

II Resonances, composite particles, and scattering amplitudes

The asymptotic states of QCD are QCD-stable hadrons, built from quarks and gluons, like the pion or the proton. Our interest in this review is in composite particles, which we will consider through their dual identities as dynamical enhancements in hadron scattering amplitudes, and as objects built from quarks and gluons. Two categories of composite particle are resonances, which are unstable, decaying into multiple stable hadrons, and bound-states, which are stable against decay by virtue of being lighter than the relevant decay threshold. Bound-states can therefore be added to the list of possible asymptotic states of QCD.

Familiar bound-states in hadron physics are atomic nuclei, with the simplest being the deuteron, a composite particle having the quantum numbers of a proton and a neutron in a $S$-wave, with small fractional coupling to $D$-wave. A detailed spectrum of hadron resonances has been observed experimentally Patrignani et al. (2016), which includes both baryons and mesons with a wide range of angular momenta. Resonances with relatively long lifetimes, narrow resonances, can often be observed as ‘bump-like’ enhancements in the invariant mass distribution of their decay products.

Bound-states and resonances may be considered more rigorously as being associated with pole singularities of scattering amplitudes. We can illustrate this in the context of elastic scattering of two spinless particles of masses $m_{1}$, $m_{2}$. The total energy and momentum in any frame, $(E,\mathbf{P})$, the center-of-momentum frame energy, $E^{\star}$, and cm-frame momentum, $\mathbf{q}^{\star}$, are related to Mandelstam $s$ by

[TABLE]

such that the magnitude of $\mathbf{q}^{\star}$ is determined,

[TABLE]

but the direction $\hat{\mathbf{q}}^{\star}$ is not, reflecting the possible angular dependence of the scattering amplitude, which can be expressed in terms of the angle $\theta_{\hat{\mathbf{q}}^{\star}}$ in the scattering plane, or via Mandelstam $t$.

The elastic scattering amplitude, $\mathcal{M}(s,t)$, can be decomposed in terms of partial-wave amplitudes, ${\mathcal{M}=\frac{1}{4\pi}\sum_{\ell}P_{\ell}(\cos\theta_{\hat{\mathbf{q}}^{\star}})\,\mathcal{M}_{\ell}(s)}$, and hadron resonances of definite angular momentum are expected to contribute to just one of the infinite set of partial-wave amplitudes. Conservation of probability above the kinematic threshold ($E^{\star}>E^{\mathrm{thr.}}=m_{1}+m_{2}$) is enforced by the elastic unitarity condition,

[TABLE]

and with the imaginary part specified by unitarity, it is common to express the real part of the elastic scattering amplitude as a function of a real variable, the phase-shift, $\delta_{\ell}(E^{\star})$,

[TABLE]

The presence of $q^{\star}$ in the unitarity relation indicates an important property of partial-wave amplitudes when they are considered to be functions of a complex value of Mandelstam $s$. The square-root in Eq.(1) means that $\mathcal{M}_{\ell}(s)$ features a branch cut beginning at threshold {s=\big{(}m_{1}+m_{2}\big{)}^{2}}, and because of this there are two Riemann sheets. The first, or physical sheet, has $\mathrm{Im}(q^{\star})>0$, and is so named because it contains the real energy axis ($s+i\epsilon$ with $\epsilon\to 0^{+}$), where physical scattering occurs. The second, or unphysical sheet, has $\mathrm{Im}(q^{\star})<0$ and can be reached by moving down through the branch cut from the real axis.

II.1 Pole singularities

The existence of a composite particle of angular momentum $\ell$ is indicated by the presence in $\mathcal{M}_{\ell}(s)$ of a pole singularity — in the vicinity of a pole at $s_{0}$, the elastic partial wave amplitude takes the form $\mathcal{M}_{\ell}\sim\frac{g^{2}}{s_{0}-s}$. Causality forbids there to be poles off the real axis on the physical sheet Gribov (2008), but poles on the real axis below threshold on the physical sheet are allowed and are identified with bound states, with the bound-state mass being $\sqrt{s_{0}}$. Pole above threshold on the real axis violate unitarity, which one case see from Eq. 2, which implies that,

[TABLE]

which cannot be satisfied at a real-valued pole.

Poles off the real axis are allowed to appear if they are on the unphysical sheet — they appear as complex conjugate pairs and can be identified with resonances. The real and imaginary parts of their position are often associated with a mass and a width for the resonance, $\sqrt{s_{0}}=m_{R}\pm i\frac{1}{2}\Gamma_{R}$. The proximity of the lower half-plane of the unphysical sheet to the real energy axis means that it is usually the pole at $\sqrt{s_{0}}=m_{R}-i\frac{1}{2}\Gamma_{R}$ which has the dominant effect on the measured scattering amplitude. Such a pole with a small value of $\Gamma_{R}$ lies close to the region of physical scattering and will give rise to a prominent ‘bump’ in the energy region around $E^{\star}\approx m_{R}$. Resonance poles lying further from the real axis, or poles which lie close to the opening of a new threshold, do not necessarily have a simple ‘bump-like’ signature on the real energy axis.

There is another possibility we have not yet considered: that a scattering amplitude has a pole singularity on the real energy axis below threshold, but on the unphysical sheet. This case is known as a virtual bound-state, and while the singularity will produce an enhancement at threshold, there is no asymptotic state possible. Virtual bound states can arise in cases when interactions are attractive, but not attractive enough to form a bound-state, with a famous experimental example being the di-neutron (spin-singlet $NN$ scattering).

We should be careful not to think of bound-states, virtual bound-states and resonances as necessarily having fundamentally different origins. For example, it has been observed in lattice QCD calculations that as the masses of the light $u,d$ quarks are increased from their physical values, the lightest resonance in $\pi\pi$ $P$-wave scattering, the $\rho$, has a width which decreases, until at a certain point the width becomes zero and the $\rho$ becomes a stable bound-state  DeGrand (1991); Bernard et al. (1993); Leinweber and Cohen (1994); Allton et al. (1999); Ali Khan et al. (2002); Bernard et al. (2001); Dudek et al. (2013a); Lin et al. (2009); Feng et al. (2011, 2015); Dudek et al. (2013b); Wilson et al. (2015a) (illustrated in Fig. 1). Another example is the di-neutron which, as presently understood, evolves from a virtual bound-state when the $u,d$ quark masses take their physical value, to a true bound-state when the quarks are somewhat heavier Berkowitz et al. (2017); Beane et al. (2012b, 2013); Yamazaki et al. (2012, 2015).

II.2 Coupled-channel scattering

Generalizing to the case where there are multiple two-body scattering channels kinematically accessible is relatively straightforward if we introduce the scattering matrix $\mathcal{M}$, with matrix elements, $\mathcal{M}_{ab}$, where $a,b$ label channels (e.g. $\pi\pi$, $K\overline{K}$, $\eta\eta$ …). The symmetry of QCD under time-reversal ensures that this matrix is symmetric, and a constraint on the imaginary part of $\mathcal{M}$ is furnished by unitarity. The unitarity condition is most compactly expressed in terms of the matrix inverse of $\mathcal{M}$ in the partial wave basis,

[TABLE]

where we see that there are only imaginary pieces along the diagonal and only above the relevant kinematic threshold. The presence of the channel momentum, $q^{\star}_{a}$, indicates that $\mathcal{M}$ now has a branch cut starting at each kinematic threshold.111This complicates the Riemann sheet structure for complex $s$, leading to $2^{n_{\mathrm{chan}}}$ sheets for $n_{\mathrm{chan}}$ open channels. Typically only a small number of those sheets are close to physical scattering.

A resonance pole will appear at the same position in each element of the matrix, $\left(\mathcal{M}_{\ell}(s)\right)_{ab}\sim\frac{r_{ab}}{s_{0}-s}$. The residue can be factorized Gribov (2008), $r_{ab}=g_{a}\cdot g_{b}$, such that we obtain couplings describing the resonance’s connection to each scattering channel.

Empirically one can measure only quantities which lead to estimates of $\left(\mathcal{M}_{\ell}(s)\right)_{ab}$ for real values of $s=E^{\star 2}$, and the behavior in the complex plane must be obtained by analytically continuing parameterized functions of $s$. Unitarity provides a strong constraint on the possible forms of such parameterizations, and a convenient way to implement this is by utilizing a $K$-matrix, where we write

[TABLE]

and where by choosing $\mathcal{K}$ to be a matrix of real functions of $s$ for real energies above kinematic threshold, we ensure the unitarity condition is satisfied.

Analyticity constraints on scattering amplitudes in principle constrain the set of allowed forms for $\mathcal{K}$, but often these are not manifestly enforced, with forms typically chosen for their simplicity used in fits to experimental data. For example, choosing the $K$-matrix to contain just a simple pole, $\mathcal{K}_{ab}=\frac{\gamma_{a}\gamma_{b}}{m^{2}-s}$, leads to a scattering matrix of the Flatté form Flatte (1976b, a), or in the single-channel case, a Breit-Wigner.

With an explicit choice of parameterization made, its free parameters can be varied to try to best describe experimental data for real energies. The amplitudes may then be analytically continued into the complex plane and examined for pole singularities. If they are found, the pole position can be interpreted in terms of a mass and width, and the pole residue can be factorized into couplings.

II.3 Diagrammatic representation

We will later explore the behavior of scattering amplitudes when our theory is placed in a finite box, and the derivations that we will present are most easily performed using a diagrammatic representation. While QCD is a theory whose fundamental degrees of freedom are quark and gluons, these fields are confined, and the relevant asymptotic states for scattering processes feature only color-singlet hadrons. The Feynman diagrams we will be discussing in the following feature only these asymptotically allowed hadrons.222One could always envision that QCD can be mapped onto an all-encompassing low-energy effective field theory, which dictates all of the dynamics of hadrons.

The diagrammatic representation of the scattering amplitude describing $n$ incoming and $n^{\prime}$ outgoing hadrons, $\mathcal{M}(n\to n^{\prime})$, is the sum over all diagrams with $n$ incoming and $n^{\prime}$ outgoing legs that have been amputated and put on-shell. All intermediate propagators are evaluated using the $i\epsilon$-prescription and all intermediate loop momenta are integrated. We make the fairly standard $i\epsilon$ prescription and momentum integration explicit in this statement because later we will have cause to adjust them.

For the simple case of elastic scattering of spinless hadrons, the sum of diagrams for $2\to 2$ amplitudes can be presented immediately, as in Fig. 2, where the diagrams can be drawn in a remarkably simple form by introducing the Bethe-Salpeter kernel, Fig. 2, and the fully-dressed single-particle propagator, Fig. 2. The use of these two objects, which are themselves infinite sums, ensures that all diagrams, including those involving intermediate multi-particle states, are included in the definition of the scattering amplitude. This seemingly simple representation is in fact exact to all orders in the perturbative expansion. We will return to this diagrammatic representation of the scattering amplitude later when we consider correlation functions evaluated in a finite spatial volume.

III Lattice QCD

Before proceeding to a discussion of the determination of resonance properties, we first present a basic overview of the numerical approach known as lattice QCD (for a more detailed introduction to lattice QCD we point the reader to Gattringer and Lang (2009)).

Lattice QCD is a non-perturbative approach to QCD, where the quark and gluon fields are quantized on a discrete grid of spacetime points of finite size Wilson (1974). By transforming to Euclidean time, the path integral comes to feature a factor which can be treated as a probability according to which gauge-field configurations may be drawn in a Monte-Carlo approach. The gauge fields in lattice QCD are expressed in terms of $SU(3)$ matrices, $U$, one for each link of the lattice, and in terms of these and the quark fields, $\psi,\bar{\psi}$, the Euclidean partition function can be expressed as

[TABLE]

where the integral is to be thought of as being over all possible configurations of the gauge and quark fields, and where the action $S_{E}$ is any suitable discretization of the QCD action. Correlation functions which feature fields at various space-time points can be similarly defined,

[TABLE]

where for a suitable choice of function we may be able to relate the value of the correlation function to a physical observable.

Because the QCD lagrangian is bilinear in $\psi$, $\bar{\psi}$, the fermionic part of the integral can be done exactly leaving

[TABLE]

where $Q(U)$ is the “Dirac matrix” which appears between $\bar{\psi}$ and $\psi$ in the discretization of the QCD lagrangian. In the case of correlation functions, all possible Wick contractions of $\psi$ and $\bar{\psi}$ in the function are replaced with “propagators”, $Q^{-1}(U)$, when the fermion fields are integrated out.

The quantity $e^{-\widetilde{S}_{E}}$ is positive definite, and can be treated as a probabilistic measure of the importance of a given gauge field configuration, to be used in a Monte-Carlo generation of a finite number of possible gauge field configurations. With these configurations in hand, the value of a correlation function can be computed on each configuration, and by averaging over the ensemble of configurations, a statistical estimate obtained.

In this approach one discretizes and truncates spacetime, which introduces respectively an ultraviolet cutoff (via the lattice spacing, $a$) and an infrared cutoff (via the lattice volume, $L$), and QCD in its usual form can in principle be recovered as $L\to\infty$ and $a\to 0$. In this document we will not be overly concerned with the behavior of lattice QCD under changes in $a$, but the dependence on $L$ will prove to be precisely the tool needed to investigate scattering of hadrons in lattice QCD.

The introduction of a space-time boundary requires us to specify boundary conditions for the fields, and the most common choice is to have fermion(gauge) fields be anti-periodic(periodic) in the temporal direction, and to have all fields be periodic in spatial directions. Periodic boundary conditions on a cube, $L\times L\times L$, mean that free particles can only have three-momenta ${\mathbf{p}=\frac{2\pi}{L}(n_{x},n_{y},n_{z})}$ for integer $n_{i}$. The finite size in the temporal direction, usually denoted $T$, effectively puts the system at a finite temperature, but provided $m_{\pi}T\gg 1$, where $m_{\pi}$ is the mass of the lightest asymptotic particle in QCD (the pion), this effect is usually negligibly small, being suppressed exponentially.

Truncation to a finite periodic spatial volume has several impacts on the theory, some of which we will take advantage of to determine scattering amplitudes — these will be discussed in the next section. Other impacts are not so useful. An example is the change in the properties of stable hadrons due to them being able to ‘see’ themselves around the periodic volume. These ‘polarization’ effects can be shown to be exponentially suppressed with increasing volume Luscher (1986a), scaling like $e^{-m_{\pi}L}$, and by working in large enough volumes, $m_{\pi}L\gg 1$, they can be reduced to a negligible level.

The parameters which must be specified in order to carry out a computation in lattice QCD include the masses of the quarks. Historically, most lattice QCD calculations have been carried out with mass values for the light up and down quarks somewhat larger than their physical masses, owing to the large computational cost of working with very light fermion fields. An obvious effect of this is that the masses of hadrons computed in these lattice QCD calculations come out larger than they would in the physical version of QCD (or in experiment). In particular it is usual to characterize the quark mass used by quoting the pion mass calculated using that quark mass. In recent years we have seen an increasing number of lattice QCD calculations (of relatively simple quantities) that use a quark mass value that is rather close to the physical value.

IV Scattering in a finite-volume

In a finite volume, strictly speaking, we cannot introduce the asymptotic states we require to define a scattering system – when we separate particles by large distances they begin to feel significant effects from the boundary, such that they are not truly free or asymptotic. Furthermore, the spectrum of eigenstates of QCD in a finite volume is qualitatively different from that in infinite volume. The continuous spectrum of multihadron states we observe in scattering experiments is present because the volume of the system is effectively infinite, allowing a particle to have any continuous value of momentum. In a finite volume, no such continuous distribution can exist, as application of the boundary conditions will quantize momenta, leading to a discrete spectrum of states.

In this chapter we will show that information about scattering amplitudes can be obtained from the discrete spectrum in a finite volume, and in particular from its dependence on the volume. We will illustrate the basic idea using the simplest possible system, non-relativistic quantum mechanics in one space dimension, before moving to the desired case, relevant to lattice QCD calculations, of scattering in quantum field theory in three space dimensions when the boundary is a cube with periodic boundary conditions. In the latter case, while conceptually the physics is the same, complications arise from the mismatch between the partial-wave expansion, which relies upon continuous rotational invariance, and the cubic geometry of the lattice boundary.

The final result, which we will refer to as the Lüscher quantization condition 333The idea of extracting scattering information from finite-volume spectra predates DeWitt (1956); Huang and Yang (1957) Lüscher’s seminal work, but he was the first to find a non-perturbative relation , is of the following form:

[TABLE]

This equation features the determinant of a sum of two complex, energy-dependent matrices, and is quite general, with a version of it being applicable to all possible $2\to 2$ scattering processes, be they elastic, coupled-channel, featuring spinless hadrons or hadrons with spin  Luscher (1986b, 1991); Rummukainen and Gottlieb (1995); Feng et al. (2004); He et al. (2005); Bedaque (2004); Liu et al. (2006); Kim et al. (2005); Christ et al. (2005); Lage et al. (2009); Bernard et al. (2011); Fu (2012); Leskovec and Prelovsek (2012); Briceno and Davoudi (2013a); Hansen and Sharpe (2012); Guo et al. (2013); Li and Liu (2013); Briceno et al. (2014); Briceno (2014). The components are $\mathcal{M}$, which is a matrix, diagonal in total angular momentum, built out of the infinite volume scattering matrices introduced in Section II.2, and $F^{-1}$ which encodes the ‘kinematics’ of the finite-volume. $F^{-1}$ is in general not diagonal in angular momentum, but is diagonal in the space of dynamically coupled channels, and it differs depending on the value of the total momentum of the two-body system, $\mathbf{P}$.

We can interpret Eq.(10) in the following way – for a given total momentum, $\mathbf{P}$, and a set of scattering amplitudes, $\mathcal{M}_{\ell}(E)$ (which might be matrices in the space of kinematically open channels), in an $L\times L\times L$ volume, the discrete spectrum of states having a specified quantum number 444we will discuss later what these ‘quantum numbers’ are with a cubic symmetry, given that the broken rotational symmetry indicates that they will not be values of angular momentum , $E_{n}(\mathbf{P},L)$, is given by all energies for which the determinant evaluates to zero. Techniques to determine the scattering amplitudes from values of $E_{n}(\mathbf{P},L)$ computed using lattice QCD will be discussed later, after we first motivate the finite-volume approach, and then sketch a derivation of Eq.( 10).

IV.1 Scattering in non-relativistic quantum mechanics in one space dimension

We can illustrate the essential relationship between the discrete spectrum of states in a finite periodic volume and the infinite volume scattering amplitudes using the simple case of elastic scattering of two identical spinless bosons in one spatial dimension in non-relativistic quantum mechanics DeGrand and DeTar (2006).

We begin by considering the infinite-volume system where we suppose that the bosons, separated by a distance $|x|$, interact through a finite-range potential, $V(|x|)$, with $V(|x|>R)=0$. Outside the potential, the wavefunction of the two-boson system will be of the form \psi_{p}(|x|)\sim\cos\big{(}p|x|+\delta(p)\big{)}, where all values of the momentum, $p\geq 0$, are allowed, and where the elastic scattering phase-shift, $\delta(p)$, describes the scattering amplitude. In principle, if the potential is specified, one can solve for the positive energy eigenfunctions inside the well (thus accounting for the dynamics of the system), and match the wavefunctions at $|x|=R$ to determine the phase-shift.

Now consider putting the system in a periodic box (a circle of circumference $L>R$). A single boson in this box would have a simple momentum spectrum, $p_{n}=\frac{2\pi}{L}n$ for integer $n$, which follows from applying the periodic boundary conditions to free-particle wavefunctions, $e^{ipx}$. For the interacting boson pair, we apply the periodic boundary conditions (at $x=\pm L/2$) to the wavefunction, \cos\big{(}p|x|+\delta(p)\big{)}, and its derivative, which leads to the following condition on the momentum:

[TABLE]

This simple result illustrates most of the important features of the spectrum in a finite-volume: there will be only discrete values of $p$ (and hence $E$) which solve this equation; if there is no scattering ($\delta(p)=0$), we recover the free-particle spectrum; when the particles interact, the discrete spectrum depends on the infinite-volume scattering amplitude (via $\delta(p)$) and the volume of the ‘box’ ($L$).

The result for quantum field theory in a three-dimensional periodic cubic box (a torus) shares all these qualitative features – let us now sketch a derivation.

IV.2 Scattering in a periodic cubic volume

To arrive at Eqn. 10, we equate two different but equivalent representations of the two-point correlation function in a finite volume: the dispersive representation which expresses the correlation function in terms of the discrete spectrum of eigenstates, and an all-orders diagrammatic representation. The derivation sketched here follows the approach first presented in Kim et al. (2005); Hansen and Sharpe (2012), using a notation similar to that in Briceno and Hansen (2016).

A two-point correlation function in a finite Euclidean spacetime of spatial dimension $L\times L\times L$ can be written555on a lattice, the integral $\int_{L}\!d\mathbf{x}$ would be replaced by a finite sum over lattice sites,

[TABLE]

where $B^{\dagger}$ and $\mathcal{A}$ are creation and annihilation operators having the quantum numbers of the hadron-hadron scattering channel we wish to study, and where we have projected into definite three-momentum, {\mathbf{P}=\frac{2\pi}{L}\big{(}n_{x},n_{y},n_{z}\big{)}}. The dispersive representation follows by inserting a complete set of discrete eigenstates of the finite-volume Hamiltonian, $C_{L}(x_{4}-y_{4},\mathbf{P})=$

[TABLE]

where the eigenstates in a finite volume are normalized according to {\big{\langle}E_{n}^{\prime},\mathbf{P}^{\prime};L\big{|}E_{n},\mathbf{P};L\big{\rangle}=\delta_{n,n^{\prime}}\,\delta_{\mathbf{P}^{\prime},\mathbf{P}}}.

The diagrammatic representation of the same finite-volume correlation function, for energies below any three-particle threshold, can be constructed according to the approach laid out in Kim et al. (2005). As illustrated in Fig. 3 those diagrams where an intermediate two-particle state can go on-shell play the dominant role in determining the dependence of the correlation function on the finite-volume. Qualitatively this can be understood by recognizing that on-shell particles can propagate over arbitrary distances and hence sample the boundaries of the volume, and the quantitative manifestation of this will be pole singularities at energies corresponding to allowed free two-particle states. Diagrams in which the intermediate two-particle state cannot go on shell can be shown to contribute at a level which is exponentially suppressed $\sim e^{-m_{\pi}L}$, and these can be neglected for volumes $L\gg m_{\pi}^{-1}$.

The core object in the diagrammatic representation is the difference between the two-particle loop in finite and infinite volume as shown in Fig. 3. We will illustrate this in the simplest case of two identical spinless particles – the extension to cases of unequal masses and particles having non-zero spin can be found in the literature Briceno and Davoudi (2013a); Hansen and Sharpe (2012); Briceno (2014); Briceno and Hansen (2016). Inside the loop, one particle carries momentum k and the other $\textbf{P}-\textbf{k}$, and their on-shell energies are ${\omega_{k}=\sqrt{k^{2}+m^{2}}}$ and ${\omega_{Pk}=\sqrt{(\textbf{P}-\textbf{k})^{2}+m^{2}}}$ respectively. In a finite-volume, the integral over continuous particle momentum is replaced by a sum over the allowed discrete momenta $\textbf{k}=(2\pi/L)\textbf{n}$ with $\textbf{n}\in\mathbb{Z}^{3}$,

[TABLE]

The difference between finite and infinite-volume loop functions, having arbitrary smooth functions $\mathcal{L}(P-k,k)$ and $\mathcal{R}^{\dagger}(P-k,k)$ at the vertices, takes the form 666For identical particles, one should include a symmetry factor of $\frac{1}{2}$. Through this review we ignore such factors, both in finite- and infinite-volume quantities, in such a way that the relationship between these is consistent.

[TABLE]

with $\Delta$ being the fully-dressed Euclidean single-particle propagator, defined with unit residue at the poles. Performing the contour integration over $k_{4}$, dropping terms which do not have singularities for physical $E$ (i.e. terms which vanish exponentially in the volume), and approximating the non-singular part of the function by the value at the pole, we arrive at

[TABLE]

The presence of the pole at $E=\omega_{k}-\omega_{Pk}+i\epsilon$ will ensure the dominance of the on-shell values of $\mathcal{L}$ and $\mathcal{R}^{\dagger}$. These have a familiar decomposition when expressed in the cm-frame. The on-shell condition written in terms of the cm-frame relative momentum, $\mathbf{q}^{\star}$, is $E^{\star}=2\omega_{q}^{\star}$, which determines the magnitude of $\mathbf{q}^{\star}$, but not its direction, which coincides with $\hat{\textbf{k}}$ boosted into the cm-frame. We can decompose $\mathcal{L}$ and $\mathcal{R}^{\dagger}$ in terms of spherical harmonics, as

[TABLE]

and the factors $\mathcal{L}_{\ell m},\mathcal{R}^{\dagger}_{\ell m}$, being independent of the loop momentum, can be taken outside the integral and sum to give

[TABLE]

where

[TABLE]

The factors \bigg{(}\frac{k^{\star}}{q^{\star}}\bigg{)}^{\ell} are introduced to cancel the ambiguity in defining the spherical harmonics at threshold.

In summary, the difference between a finite-volume two-particle loop and its infinite-volume counterpart can be written as a product of matrices in angular momentum space featuring only on-shell quantities. We may apply this result to the diagrams appearing in Fig. 3 — there the functions appearing inside the loops are either Bethe-Salpeter kernels or overlaps with the source and sink operators. For each diagram we replace it with its infinite volume expression plus a finite-volume correction which is functionally identical to Eq.(14) with $\mathcal{L}$ and $\mathcal{R}^{\dagger}$ replaced with the appropriate functions. The sum of all $2\to 2$ on-shell amputated diagrams is replaced with the scattering amplitude 777whose partial-wave expansion we generalize to ${\mathcal{M}(P,\hat{\mathbf{q}}^{\prime\star},\hat{\mathbf{q}}^{\star})=Y^{*}_{\ell m}(\hat{\mathbf{q}}^{\star})\,\mathcal{M}_{\ell m;\ell^{\prime}m^{\prime}}(E^{\star})\,Y^{*}_{\ell^{\prime}m^{\prime}}(\hat{\mathbf{q}}^{\prime\star})}$ for arbitrary directions of initial and final-state particles, realizing that rotational invariance in infinite volume will make the matrix $\mathcal{M}_{\ell m;\ell^{\prime}m^{\prime}}$ diagonal in angular momentum $\mathcal{M}$. After a substantial amount of algebra, one can show that the correlation function can be written in terms of the infinite-volume correlation function, $C_{\infty}$, and a geometric series in $-\mathcal{M}\,F$ which appears between the fully dressed overlaps to the source and sink operators $B^{\star}$ and $A^{\star}$,

[TABLE]

In order to complete the derivation by comparing to the dispersive representation, Eq.(12), we must evaluate the integral over $P_{4}$, which can be done by considering the analytic properties of the integrand. In order for Eq.(12) to result, the integrand must include a series of poles of the form $\frac{1}{P_{4}-iE_{n}(L)}$. Since the position of the poles corresponds to the spectrum in a finite volume, they are surely volume dependent, and as such they must arise from the factor \big{[}F^{-1}(P,L)+\mathcal{M}(P)\big{]}^{-1} in Eq.(16). This factor will diverge appropriately if the matrix \big{[}F^{-1}(P,L)+\mathcal{M}(P)\big{]} is singular, that is if \det\big{[}F^{-1}(P,L)+\mathcal{M}(P)\big{]}=0, which we recognize as Eq.(10).

The relationship between the finite-volume spectrum and the infinite-volume scattering amplitudes, Eq.(10) has thus been derived, but one additional useful result concerning matrix-elements follows if we evaluate the $dP_{4}$ integral, giving,

[TABLE]

where $\mathcal{R}$ is the residue of the finite-volume two-particle propagator,

[TABLE]

We will consider the importance of this result later when we investigate the computation of resonant matrix elements in QCD.

IV.3 Relating scattering amplitudes to finite-volume spectra

Eq.(10), which we repeat here,

[TABLE]

describes the relationship between the finite-volume energy spectrum, $E_{n}(\mathbf{P},L)$, and the infinite-volume scattering amplitudes. As the derivation in the previous section indicates, the matrices $F^{-1}$ and $\mathcal{M}$ are, even in the simplest elastic spinless case, formally of infinite size, featuring all possible integer values of $\ell\geq 0$. How then can we hope to use this equation in practice? To make progress we first note that at low energies, higher partial-wave amplitudes typically reduce rapidly in magnitude as $\ell$ is increased. This is expected on the grounds of angular momentum conservation which leads to a behavior at threshold which much fall at least as fast as $\mathcal{M}_{\ell}\sim(q^{\star})^{\ell}$. It follows that, in practice, the contribution of higher partial waves to Eq.(10) is numerically negligible at low energies, and one is well justified truncating the angular momentum space to the lowest few $\ell$ values.

An additional simplification comes from considering the fact that, although the use of a cubic symmetry for the lattice boundary has broken the continuous rotational symmetry needed for angular momentum to be a good quantum number, there is a smaller residual symmetry group still present. In the case of the rest frame ($\mathbf{P}=\mathbf{0}$), this is the group of rotations which leave a cube invariant, the cubic or octahedral group. This group has a finite number of irreducible representations, or irreps. For example, for bosonic systems the irreps are labelled $A_{1},A_{2},T_{1},T_{2}$ and $E$ (plus an additional parity label), and these can be though of as the ‘quantum numbers’ carried by finite-volume eigenstates. When projected into one of these irreps (subduced is the term-of-art used), the dense matrix $F^{-1}$ (whose rows and columns are labelled by $\ell,m$ values) becomes block diagonal, reflecting the fact that only certain $\ell$ values subduce into each irrep, see Table 1.

For non-zero momenta $\mathbf{P}$, the symmetry is reduced still further, to the appropriate little group, each of which has its own set of irreducible representations. We will not go into detail here (see e.g. Thomas et al. (2012)) except to note that in general, these in-flight irreps are typically more dense in $\ell$ content and usually feature both parities.

Although subduction simplifies the behavior of $F$ in angular momentum space somewhat, the dependence on energy remains complicated, featuring many singularities — some illustration is provided in Appendix XI.

IV.3.1 Dominance of the lowest partial-wave

An extremely important case, the one most commonly considered in the literature to date, has only the lowest partial-wave subduced into an irrep having a non-negligible amplitude, with all higher partial-waves assumed small enough to ignore. In the case of elastic scattering this assumption leads to $F^{-1}$ and $\mathcal{M}$ being $1\times 1$ matrices, and Eq.(10) reducing to a simple equation, $\mathcal{M}=-F^{-1}$, which is often expressed in terms of the elastic scattering phase-shift as

[TABLE]

where the right-hand-side is related to the finite-volume function we have previously discussed, ${\cot\phi^{\Lambda}_{\ell}(\mathbf{P},L,E)=16\pi\,\frac{E^{\star}}{2q^{\star}}\,\mathrm{Re}\,F^{\Lambda}_{\ell}}$, where the energy dependence of $F$ is illustrated in Appendix XI.

Under these conditions, determining the scattering amplitude from values of $E_{n}(\mathbf{P},L)$ computed using lattice QCD is straightforward – one simply inserts the lattice energy value into the function on the right-hand-side of Eq.(19), to yield a value of $\delta_{\ell}$ at that energy. The more energy levels are determined, the more discrete points one will have on the phase-shift curve.

One might wonder how we can know, in any particular calculation, whether we are justified in neglecting higher partial waves. The simplest approach is to directly determine the corresponding higher partial-wave phase-shift using a relevant irrep. For example, if we are considering $\pi\pi$ scattering with $I=1$, we might be concerned about the contribution of the $\ell=3$ partial-wave to the $T_{1}^{-}$ irrep, which should be dominated by $\ell=1$ at low energies. We can use energy levels in the $T_{2}^{-}$ and $A_{2}^{-}$ irreps, which have $\ell=3$ as their lowest partial wave888 $\ell=2$ cannot contribute to $I=1$ in this case due to Bose symmetry. We can also make use of levels in moving-frame irreps which have $\ell=3$ as their leading contribution., assuming that still higher partial-waves are negligible.

Going beyond this simplest case of a single partial wave dominating, Eq.(10), suitably subduced and truncated to include only the lowest few relevant partial-waves, in principle has each finite-volume energy level as a function of the phase-shift for all included partial waves. Of course, level–by–level this is just one equation in multiple unknowns, and it cannot be solved. However, in practice, we will determine many energy levels with each one providing a constraint. If the scattering amplitudes are smooth functions of the energy, they can be suitably parameterized, and we can attempt to describe the spectrum as a whole, by varying a small number of parameters. A similar situation arises in the case of coupled-channel scattering, and indeed is more acute there, as there is usually no sense in which one channel can be considered ‘weak’ relative to a ‘dominant’ channel, so we will illustrate the parameterization of amplitudes approach there.

IV.3.2 Coupled-channel scattering and parameterization of scattering amplitudes

We can illustrate the approach most easily if we initially assume a system dominated by $S$-wave scattering of two coupled channels of spinless particles that we will label $\pi\pi$ and $K\overline{K}$. If $m_{\pi}<m_{K}$ there is a region of elastic $\pi\pi$ scattering before the $K\overline{K}$ threshold opens. $S$-wave scattering is described by the matrix of scattering amplitudes,

[TABLE]

where each entry is a complex number at each value of cm-energy. Multichannel unitarity, Eq.(5), specifies the imaginary parts so that just three real numbers are needed at each value of energy to completely specify the scattering. This matrix appears in a version of Eq.(10),

[TABLE]

where $F^{-1}_{\pi\pi}$ and $F^{-1}_{K\overline{K}}$ are known functions of $E$ and $L$ whose imaginary parts are such that this equation is real for matrices $\mathcal{M}$ satisfying the multichannel unitarity condition (see Appendix XI). It follows that for any particular finite-volume energy level $E_{n}(L)$, we have one real equation and three unknowns.

An approach that has proven successful in a number of explicit calculations (which will be discussed in Section VI) was outlined in Guo et al. (2013). In it the energy dependence of the scattering matrix is parameterized by some explicit analytic form featuring a small number of free parameters, $\{a_{i}\}$. For any particular set of values of these free parameters, the explicit parameterized form for $\mathcal{M}$ can be substituted into Eq.(20) and that equation solved for a discrete spectrum of finite-volume energies in an $L\times L\times L$ box, $E^{\mathrm{par.}}_{n}(L;\{a_{i}\})$. This spectrum can be compared to the spectrum obtained in a lattice calculation, and by varying $\{a_{i}\}$ a $\chi^{2}$ can be minimized to find a best description.

Of course this approach requires us to propose particular parameterization forms for $\mathcal{M}$. Since they must satisfy multichannel unitarity if Eq.(20) is to have solutions, the question of how to parametrize $\mathcal{M}$ is equivalent to to choosing a parametrization for the $K$-matrix introduced in Eq.(6). One might be correctly concerned that by having choosen a particular parametrization for the $K$-matrix one might introduce possible systematic bias. This can and should be tested by varying the choice made in the parametrization. It is empirically observed in explicit calculations presented in Section VI that in cases where a large density of energy levels are determined, strongly constraining the free parameters in the parameterization, any sufficiently flexible form for the $K$-matrix will give compatible amplitudes. Furthermore, if the channel contains a fairly narrow resonance, the analytic continuation of all parametrization forms into the complex plane will show pole singularities at consistent locations.

IV.3.3 Examples of finite-volume spectra for simple scattering amplitudes

In this subsection we will present some illustrative examples of the finite-volume spectra obtained by solving Eq.(10) for some simple scattering amplitudes. The very simplest case is that of no scattering at all, and in this case the finite-volume spectrum will simply be the discrete spectrum of non-interacting hadron pairs. For example for a pair of spinless particles, which we might call $\pi\pi$, in the rest frame, $E(L)=2\sqrt{m_{\pi}^{2}+k^{2}}$ with $\mathbf{k}=\frac{2\pi}{L}(n_{x},n_{y},n_{z})$.

The simplest non-trivial case we will consider is weakly interacting elastic $S$-wave scattering, described by $\tan\delta(q^{\star})=a\,q^{\star}$, with scattering length $a<0$ being repulsive and $a>0$ being attractive. Figure 4 shows the finite volume spectra in the rest frame $\mathbf{P}=[000]$ ($A_{1}^{+}$ irrep) and, for illustration, one possible moving frame, $\mathbf{P}=[110]$ ($A_{1}$ irrep). We observe that in the weak attractive case ($a>0$, green curves), the energies lie close to, but systematically lower than, the non-interacting energy levels (black dashed curves), while in the repulsive case ($a<0$, red curves) they lie systematically higher.

A more interesting case is presented in Figure 5, where the elastic amplitude is a relativistic $S$-wave Breit-Wigner,

[TABLE]

with $m$ chosen to be 1182 MeV and three increasing values of coupling, $g$. For the narrowest resonance ($g=1.0$, upper panel), it is clear that at every value of $L$ there is an ‘extra’ level (beyond those expected in the non-interacting case) in the vicinity of $1182$ MeV, and in those locations where a non-interacting curve (black dashed curves) crosses this energy, there is an ‘avoided level crossing’. Levels at energies far from $1182$ MeV are observed to lie very close to the non-interacting curves, as we might expect given that at those energies, $\delta_{0}\approx 0^{\circ}$ or $180^{\circ}$, either of which correspond to no scattering. As the coupling is increased (lower panels), the avoided level crossings become broader and the effect of the resonance is effectively spread over a larger energy region, where it is no long possible identify any single level as ‘belonging’ to the resonance.

Figure 6 shows an example of an $S$-wave two-coupled-channel process constructed to feature a simple narrow resonance coupled to both channels, which we label $\pi\pi$ and $K\overline{K}$, where these channels are otherwise uncoupled. A simple Flatté form achieves this 999 in $S$-wave, $\mathcal{M}_{ij}(E^{\star})=16\pi\frac{g_{i}g_{j}}{m^{2}-E^{\star 2}-i\sum_{k}g_{k}^{2}(2q^{\star}_{k}/E^{\star})}$ , and gives us our canonical view of a multichannel narrow resonance: it appears as a bump in the channels to which it couples, and the height of the bump is simply related to the coupling to that decay channel. It has a nearby pole on the nearest unphysical sheet (since we are above both threshold in two-channel scattering, this is sheet III, $\mathrm{Im}\,q^{\star}_{\pi\pi}<0$, $\mathrm{Im}\,q^{\star}_{K\overline{K}}<0$) that is partnered by an ‘image’ pole Au et al. (1987) on a more distant Riemann sheet (in this case, sheet IV, $\mathrm{Im}\,q^{\star}_{\pi\pi}>0$, $\mathrm{Im}\,q^{\star}_{K\overline{K}}<0$) .

Focussing on the corresponding finite-volume spectrum we notice that, as was the case in elastic scattering, there is a typically an ‘extra’ state in the vicinity of the resonance mass. The characteristic ‘avoided level crossing’ behavior is again present. It is clear that computing in moving frames (such as $\mathbf{P}=[110]$ shown in the figure) provides a high density of levels which can be used to constrain the scattering matrix.

We end this section by pointing out that the physics of multichannel resonances in infinite and finite volumes is not always as simple as was implied by the previous example. Figure 7 shows a two-channel amplitude which features a resonant pole causing a rapid energy variation at the opening of the $K\overline{K}$ threshold 101010The scattering matrix here is constructed from a $K$-matrix parameterized via

\mathcal{K}^{-1}=\left(\begin{array}[]{cc}a&b+c(E^{\star})^{2}\\ b+c(E^{\star})^{2}&d+e(E^{\star})^{2}\end{array}\right),

(22)

where $a,b,c,$ and $d$ are constants, and we have replaced the phase space factor in Eq.(6) with the so-called Chew-Mandelstam phase space, first introduced in Chew and Mandelstam (1960) and described in detail in Wilson et al. (2015b).. It is clear that this resonance does not manifest itself as a bump, but rather as a dip in the magnitude of the $\pi\pi\to\pi\pi$ amplitude, and as a rapid turn-on of the $\pi\pi\to K\overline{K}$ and $K\overline{K}\to K\overline{K}$ amplitudes at threshold. There is a pole on sheet II whose presence is being felt at the $K\overline{K}$ threshold – a very distant pole on sheet III has little impact on the behavior at $K\overline{K}$ threshold. We will not discuss here what physics might cause such a pole distribution, only point out that this is a possible scattering amplitude, one which may be somewhat reminiscent of the experimental $f_{0}(980)$ resonance.

Now, examining the finite-volume spectra corresponding to this amplitude, we see that, while there can be large departures from the non-interacting spectrum in the region around the resonance position (which is very close to the $K\overline{K}$ threshold), there is no obvious ‘extra’ level present. This example illustrates the case that in multichannel scattering, a presence of an additional level might suggest the existence of a resonance in the theory, but the converse, that the absence of an ‘extra’ level implies the absence of a resonance, is certainly not true.

V Determining the finite-volume spectrum

It should be clear from the results of Section IV that one can learn about scattering amplitudes if one has determined the discrete spectrum of states in one or more volumes, in one or more frames. In this section we present a discussion of how the reliable extraction of a discrete excited state spectrum can be achieved in practical lattice QCD calculations.

V.1 Variational analysis of correlation matrices

From the form of Eq.(12),

[TABLE]

it would seem that by decomposing the time-dependence of a single correlation function in terms of a sum of exponentials one could get access to the entire spectrum of states. Unfortunately this approach is not practical for the determination of anything beyond the energy of the lightest state with the quantum numbers of $\mathcal{A}$, $\mathcal{B}^{\dagger}$. A technique that has proven extremely effective is to compute a matrix of correlation functions using a basis of operators,

[TABLE]

which can be analyzed variationally Michael (1985); Luscher and Wolff (1990); Blossier et al. (2009). The set of $N$ operators $\{\mathcal{A}_{i}\}_{i=1\ldots N}$ all have the same quantum numbers but will ideally have differing magnitudes of overlap $\langle 0|\mathcal{A}_{i}|E_{n}\rangle$, onto each state in the spectrum. It follows that there will be a certain linear superposition of operators, $v_{i}^{(n)}\mathcal{A}_{i}$, that will optimally (in the variational sense) interpolate each state $n$ in the spectrum, and these superposition weights can be shown to be the eigenvectors which solve a generalized eigenvalue problem,

[TABLE]

In this expression, $t_{0}$ is a suitably chosen reference time, which ideally is selected to be such that only $N$ states are required to saturate $C_{ij}(t_{0})$ to a good approximation Dudek et al. (2008). The eigenvalues feature the relevant state energies, ${\lambda_{n}(t,t_{0})\sim e^{-E_{n}(t-t_{0})}}$. Fitting the time dependence of each of these quantities, also known as principal correlators, therefore enables a determination of the discrete spectrum. There are various approaches to handling the time-dependence of the solution of Eq.(23), see e.g. Dudek et al. (2008); Bulava et al. (2009); Mahbub et al. (2013); Kiratidis et al. (2015, 2017).

V.2 Operator construction

There still remains the question of what form the operators, which must have the quantum numbers of the hadronic system of interest, but which need to be constructed from the basic quark and gluon fields of QCD, should take. A longstanding approach for the case of systems with meson quantum numbers is to make use of fermion bilinear operators, $\overline{\psi}\bm{\Gamma}\psi$, where the quark fields may be smeared over space, and where the $\bm{\Gamma}$ object controls the spin-structure and may, if desired, have spatial dependence. A recent application of this approach Dudek et al. (2009, 2010, 2011a); Thomas et al. (2012); Dudek et al. (2013a) considers a large basis of operators constructed out of the following basic form,

[TABLE]

where the $D_{i}$ are gauge-covariant derivatives in the spatial directions. Using up to three derivatives, a very large basis of operators can be constructed, and with these operators, large matrices of correlation functions computed. When these are analyzed variationally, spectra like those shown in Figure 8 can be obtained Dudek et al. (2013a). The set of states observed closely matches the expectations of a model where mesons are quark-antiquark constructions, with the addition of some states that also have gluonic content Dudek (2011). But this spectrum is surely incomplete — in particular when the volume is changed, the extracted spectrum varies relatively little, in contrast to the expectations of the previous section (where for example at least some levels approximately track the volume dependence of non-interacting levels). Subsequent calculations using an augmented operator basis have shown that this spectrum is indeed incomplete, and it is to the form of the required additional operators that we now turn.

V.2.1 The importance of “multi-hadron” operators

In order to resolve the complete low-lying spectrum of states through variational analysis of correlation matrices, it proves to be necessary to include in the operator basis some operators which ‘resemble’ the expected multi-hadron states in that energy region. This was illustrated in Wilson et al. (2015a) for the case of the $\rho$ resonance appearing in $\pi\pi$ scattering. The operator basis used (in the $\mathbf{P}=\mathbf{0}$ frame) included a large number of operators of the ‘single-meson’ $\overline{\psi}\bf{\Gamma}\psi$ type described in the previous section, and in addition several operators of form $\sum_{\hat{\mathbf{k}}}c_{\hat{\mathbf{k}}}\,\pi(\mathbf{k})\,\pi(-\mathbf{k})$, where $\pi(\mathbf{k})$ is a shorthand notation for a superposition of $\overline{\psi}\bf{\Gamma}\psi$ operators which optimally 111111Optimal operator determined through variational analysis of a matrix of correlators for the irrep in which the moving pion sits.. interpolates the lightest pseudoscalar with momentum $\mathbf{k}$ The product of two “pion” operators, in a fictitious world where pions do not interact, would closely resemble an eigenstate of the system with energy $2\sqrt{m_{\pi}^{2}+\mathbf{k}^{2}}$. These ‘meson-meson’ operators are qualitatively different from the ‘single-meson’ operators — each meson operator is projected into a definite momentum and this means that they effectively sample the entire lattice volume.

Figure 9 shows the result of variational analysis of a large matrix of correlation functions performed on an $L^{3}\sim(4\,\mathrm{fm})^{3}$ lattice. In the lattice units presented, the pion mass is $a_{t}m_{\pi}=0.0393$ and the kaon mass is $a_{t}m_{K}=0.0834$. The first column shows the spectrum extracted using 30 operators including several ‘$\pi\pi$-like’ and ‘$K\overline{K}$-like’ operators. The histograms show the relative size of overlaps, \big{\langle}0\big{|}\mathcal{A}(0)\big{|}E_{n},\mathbf{P};L\big{\rangle}, of each state onto each operator in the basis, and it is clear that the lowest two levels have strong overlap with both ‘single-meson’ operators subduced from $J=1$ (orange) and $\pi\pi$ operators (red). The third column indicates the spectrum which is extracted if meson-meson operators are not included in the basis, and it is immediately clear that the spectrum is completely different, in particular in place of two low-lying states in the complete spectrum, only one state appears. It appears that without the meson-meson operators the variational system cannot find two orthogonal combinations to overlap with the two states.

A rather simple explanation for this was presented in Dudek et al. (2013b). Imagine one could turn off the $\rho\to\pi\pi$ coupling and all $\pi\pi$ interactions, but retain the $\rho$ as a stable state in the spectrum — in this limit, the finite-volume eigenstates will be a single $|\rho,L\rangle_{0}$ state of definite mass, and free $|\pi(\mathbf{k})\pi(\mathbf{k}),L\rangle_{0}$ states at energies $2\sqrt{m_{\pi}^{2}+\mathbf{k}^{2}}$. These states would be interpolated by ‘single-meson’ operators and ‘meson-meson’ operators respectively. With a small but not zero coupling $\rho\to\pi\pi$, these basis states will admix to form finite-volume eigenstates, and in the case that the $\rho$ lies near just one $\pi\pi$ state we can treat this as a two-state system,

[TABLE]

and in the expected case that a ‘single-meson’ operator $\mathcal{O}$ overlaps significantly only with $|\rho,L\rangle_{0}$, the corresponding correlator would behave like

[TABLE]

We note that as soon as any dynamical coupling $\rho\to\pi\pi$ is turned on, the mixing angle $\theta$ will be particularly sensitive to $L$ as the $\pi\pi$ (non-interacting) energy passes the $\rho$ “mass”. If the free particles’ energy is in the near vicinity (loosely, within the width) of a resonance, both $\cos^{2}\theta$ and $\sin^{2}\theta$ will be appreciable. It is therefore physically not possible for a local $\rho$-like operator to have significant preference to couple to either $E_{1}$ or $E_{2}$. If we work with an operator basis that only has appreciable couplings to the compact $\rho$-like state it would prove impractical with finite statistics and a limited time-extent to resolve the presence of two exponentials. As a consequence we would end up concluding the presence of a single state whose mass is somewhere between $E_{1}$ and $E_{2}$.

It can be justified that the ‘single-meson’ operator has suppressed overlap onto the $\pi\pi$ state by appealing to the difference in spatial distributions of the operator and the state. The operator is a local object, which samples only a small region of space, while the two-meson state occupies the entire lattice volume, and hence the overlap is suppressed by the volume. We emphasize here that operators interpreted as being of ‘single-hadron’ type need not be just fermion bilinears (or three quark constructions in the case of baryons), but may contain any number of fermion fields (for example, tetraquark meson constructions, $\sim\overline{\psi}\overline{\psi}\psi\psi$). What distinguishes them from multi-hadron operators is how they sample the spatial volume, in particular that they are local to a restricted region of space rather that the full spatial volume.

It is clear that we need to include ‘meson-meson’ operators into our basis, but how many do we need? As an example, do we need to include $K\overline{K}$-like operators when studying $\pi\pi$ scattering, since these two channels can, in principle, couple. This is illustrated in the second column of Figure 9 where we observe that excluding $K\overline{K}$-like operators does not change the spectrum below $K\overline{K}$-threshold. On the other hand, above $K\overline{K}$-threshold, where we expect there to be states which contain an admixture of $K\overline{K}$, we observe that one state disappears when the relevant operators are not included.

We conclude this discussion by providing a ‘rule of thumb’ for the required operator basis: *in addition to ‘single hadron’ operators capable of overlapping onto any resonances you might anticipate, be sure to include multihadron operators for all expected non-interacting levels in the energy region you wish to study *.

In light of this, how should we view the spectra presented in Fig. 8, which were extracted without including any multihadron operators, but using a large basis of ‘single-hadron’ operators? We might guess that we are overlapping only onto the ‘single-hadron-like’ parts (pieces analogous to $|\rho,L\rangle_{0}$, but for other resonances) and hence the spectrum is indicating the presence of a resonance (which is probably narrow) in some energy region, without precisely determining its properties. We will return to this question in Section VII.1, where we derive a relevant result for how narrow resonances might manifest themselves in finite-volume spectra.

VI Examples of resonance determination

Recent years have seen an increasing number of lattice QCD calculations making use of the Lüscher approach in order to determine hadron resonance properties. In this section we will summarize this progress. We begin by reminding the reader that the procedure laid out in the previous two sections is not only useful in the extraction of resonances, scattering amplitudes in which no resonance appears can also be determined. $\pi\pi$ scattering with isospin=2 is the classic example, and this case has been studied in a number of lattice calculations  Sharpe et al. (1992); Gupta et al. (1993); Kuramashi et al. (1993); Fukugita et al. (1995); Aoki et al. (2002); Du et al. (2004); Yamazaki et al. (2004); Beane et al. (2006, 2008); Li et al. (2007); Feng et al. (2010); Yagi et al. (2011); Sasaki et al. (2014); Fu (2013); Helmes et al. (2015); Dudek et al. (2011b); Beane et al. (2012c); Dudek et al. (2012); Bulava et al. (2016). Several recent calculations Dudek et al. (2011b); Beane et al. (2012c); Dudek et al. (2012); Bulava et al. (2016) determine multiple energy levels and use these to determine the energy-dependence of the elastic phase-shift.

One application of these calculations has been to study the pion mass dependence of the scattering length Yamazaki et al. (2004); Beane et al. (2006, 2008); Feng et al. (2010); Beane et al. (2012c); Yagi et al. (2011); Sasaki et al. (2014); Fu (2013); Helmes et al. (2015), which may be compared to the expectations of chiral perturbation theory  Weinberg (1966); Colangelo et al. (2001), and these studies find that chiral perturbation theory describes the scattering length up to surprisingly large values of the pion mass, $m_{\pi}\sim 400$ MeV. The energy dependence of the $S$-wave $\pi\pi$ interaction has been studied by the NPLQCD Collaboration by considering low-lying rest and boosted lattice spectra across multiple volumes. These results, at a pion mass of $m_{\pi}\sim 390\,{\rm MeV}$, sufficiently constrain the NLO chiral expansion such that the effective range parameters and be extrapolated to the physical pion mass. The resulting phase shift is shown in Figure 10, where it is seen to be in reasonable agreement with experimental data.

Using the same gauge field configurations as the NPLQCD study, the Hadron Spectrum Collaboration Dudek et al. (2012) considered multiple irrep spectra in various moving frames to determine the energy dependence of the $S$ and $D$-wave elastic scattering phase-shifts up to the $4\pi$ threshold, as shown in Figure 11. The effect of the angular momentum barrier at threshold is clearly observed, with the $D$-wave being significantly reduced in magnitude with respect to the $S$-wave.

VI.1 Elastic resonances in $\pi\pi$ scattering

The most widely studied example of resonance extraction in lattice QCD is the $\rho$ resonance in $P$-wave isospin=1 $\pi\pi$ scattering  Aoki et al. (2007); Feng et al. (2011); Lang et al. (2011); Aoki et al. (2011); Dudek et al. (2013b); Pelissier and Alexandru (2013); Wilson et al. (2015a); Bali et al. (2016); Bulava et al. (2016). At all but the smallest pion masses, the $\rho$ is an elastic resonance, kinematically able to decay only to $\pi\pi$, and even at the physical pion mass, although the $\pi\pi\pi\pi$ channel is kinematically open, it has negligible coupling to the $\rho$ Patrignani et al. (2016).

The pioneering application of the Lüscher approach to this case was reported on in Aoki et al. (2007). This early work considered a relatively small volume, determining two energy levels in the rest frame using a basis of operators featuring a $\pi\pi$-like construction and a $\overline{\psi}\Gamma\psi$ construction with vector quantum numbers. The two energy levels were utilized to give two points on the phase-shift curve using Eq.(19), which is the minimum required to determine the two free parameters appearing in a Breit-Wigner description of the scattering amplitude, while of course not providing a measure of goodness-of-fit.

Later works, in particular Feng et al. (2011), computing for several values of the pion mass, extended the approach by making use of some moving frames to obtain more energy levels and hence more points on the phase-shift curve. With one volume per pion mass, Feng et al. (2011) had typically three points in the energy region corresponding to the resonance (and more outside this region), giving slightly more constraint on the resonance parameters. At about the same time, Aoki et al. (2011) also considered moving frames.

Dudek et al. (2013b) and Wilson et al. (2015a) perhaps best illustrate the power of extracting multiple energy levels using multiple volumes and/or several moving frames. Figure 12, taken from Dudek et al. (2013b), shows the elastic phase-shift determined at $m_{\pi}\sim 391$ MeV using energy spectra in three volumes, $16^{3},20^{3},24^{3}$ ($L\sim 2.0,2.4,2.9$ fm), in all relevant irreps, in all frames up to $\mathbf{P}=[200]$. There can be no doubt from these data points that there is a resonance present with the energy dependence of the phase-shift being mapped out in detail across the entire elastic region, and the resonance parameters being tightly constrained. Figure 13, taken from Wilson et al. (2015a) illustrates that the use of multiple moving frames on a single larger volume (which has a denser energy spectrum) can lead to the same level of detail in the mapping out of the elastic phase-shift. These calculations also computed spectra in irreps whose lowest angular momentum is not $\ell=1$, but rather $\ell=3$, and in this way they were able to place some constraint on the size of the $F$-wave phase-shift, confirming that it is indeed negligibly small in the energy region where the $\rho$ appears, justifying the neglect of higher partial-waves in the finite-volume analysis.

Also illustrated in Figure 13 are the resonance pole positions found for a variety of amplitude parameterizations constrained to describe the finite volume spectra. It is clear from the tiny degree of scatter in the pole position, that the presence of a pole at that particular complex energy value is required, regardless of the other details of the parameterization, to describe this scattering system.

A calculation on the same lattice configurations as Wilson et al. (2015a), but using a different correlator construction technique, operator basis and variational analysis method is presented in Bulava et al. (2016). As shown in Figure 14, the phase-shift extracted is compatible with that presented in Figure 13 within the larger statistical errors.

The $\rho$ is an isospin=1 resonance in $\pi\pi$ scattering in $P$-wave, and with this isospin, this is lowest allowed partial wave — $\pi\pi$ scattering in $S$-wave cannot occur with isospin=1, but can with isospin=0. Experimentally the $S$-wave is very different to the $P$-wave, featuring not a narrow resonance “bump”, associated with a rapid rise in phase-shift through $90^{\circ}$, but rather a gradual increase in phase-shift from $\pi\pi$ threshold up to the $K\overline{K}$ threshold, where some more rapid variation (associated with the $f_{0}(980)$ resonance) occurs. It has long been suspected that the slow increase in phase-shift is associated with a very broad resonance, the $\sigma$, whose pole is located far into the complex energy plane. Recently, a highly constrained analysis of experimental $\pi\pi$ scattering data has confirmed that such a pole is present and determined its position with some precision (see Pelaez (2016) for a review of the situation).

Very recently we have seen the first serious lattice QCD determination of elastic $\pi\pi$ scattering in the isospin=0 $S$-wave Briceno et al. (2017a). Lattice QCD calculation of this channel had long been considered extremely challenging owing to the need to compute diagrams in which all the quarks and antiquarks annihilate, leading to something which is completely disconnected. By computing a large number of propagation objects in the distillation framework Peardon et al. (2009), the Hadron Spectrum Collaboration were able to compute the required correlation functions and obtain finite-volume spectra at two pion masses, $m_{\pi}\sim 236,391$ MeV. The lattices are the same ones used in the $\rho$ extractions described above, with three volumes at the heavier mass and a single larger volume at the lighter mass.

Figure 15 shows the elastic scattering phase-shift determined from spectra on these lattices for the two pion masses, and a clear change is observed between the two. At the heavier quark mass, the behavior is that of a bound-state lying just below threshold, while at the lighter mass we observe something much closer to the experimental situation, with a slow increase in phase-shift over the elastic region.

At the heavier quark mass, all analytic parameterizations of the scattering amplitude capable of describing the finite-volume spectra feature a pole located on the real energy axis, on the physical sheet, at ${E^{\star}=758(4)}$ MeV, which is interpreted as a bound-state $\sigma$ (lying below the $\pi\pi$ threshold at $2\,m_{\pi}=782$ MeV). At the lighter quark mass, the situation is somewhat less clear — many different parameterizations are capable of describing the spectra, and while they do feature a pole far into the complex plane on the unphysical sheet, the position of that pole is not precisely determined, with considerable scatter observed as the amplitude parameterization is varied. This observation is not unique to the finite-volume situation — the same scatter in pole position is observed when a variety of amplitudes forms are constrained using only the experimental elastic phase-shift data. It is only when amplitude forms which build in the required constraints of analyticity and crossing-symmetry are utilized that the pole position can be pinned down with precision Pelaez (2016).

This brings us back to a point previously raised in the context of analyzing the spectra that lie above multiple two-body open channels, discussed in Sec. IV.3.2. We found that the quantization condition describing such systems, Eq.(20), relates a single energy level to multiple components of the scattering amplitude. We argued this issue can be circumvented by using flexible parametrizations of the energy dependence of the scattering amplitude. If a large density of states are determined, the systematic error in the resulting scattering amplitude on the real axis is expected to be small, because the quantization conditions constrains different linear combinations of the scattering amplitude in close proximity. That being said, it is quite possible that consistent results on the real energy axis might lead to significantly different values for the position of any resonant pole if it lies far into the complex plane. This motivates the need for further constraints on the scattering amplitude – as an example, for light quark masses, it might be important to implement low-energy constraints imposed by chiral perturbation theory, such as the Adler zero present in the isoscalar $\pi\pi$ channel, which plays an important phenomenological role.

A system rather similar to the $\rho$ in $\pi\pi$ is the narrow $K^{\star}$ meson in $P$-wave $\pi K$ scattering with isospin=1/2. There is an important difference, however, when this system is considered in a finite volume — the $\pi K$ $S$-wave system does not decouple from the $P$-wave in moving frames, which significantly complicates the analysis since the $S$-wave, like the $\pi\pi$ case, appears to house a broad resonance. A few lattice calculations have considered elastic scattering in the $\pi K$ channel Bali et al. (2016); Lang et al. (2012); Fu and Fu (2012); Prelovsek et al. (2013).

A number of papers have reported on the determination of elastic scattering amplitudes involving mesons containing a heavy quark Mohler et al. (2013); Prelovsek and Leskovec (2013a); Lang et al. (2014b, 2015, 2016) – these studies have typically considered only the rest frame in a single volume, giving a limited constraint on the energy dependence of the phase-shift. An example is shown in Figure 16 which displays the $DK$ $S$-wave elastic phase-shift at two pion masses computed in Lang et al. (2014b). The use of a single volume and only the rest frame limits the number of relevant energy levels to two, with an effective range formula used to interpolate between the two in order to locate the position of a bound-state $D_{s0}^{\star}$ meson.

The need to determine multiple energy levels in a limited energy region in order to map out the scattering amplitude has lead to suggestions for approaches other than the use of moving frames. One example is to consider lattices in which one spatial direction is longer than the others, i.e. $L\times L\times(L^{\prime}>L)$. In this way the unit of momentum in the $z$-direction is smaller than that in the $x,y$ directions, leading to a changed spectrum of non-interacting levels. In the interacting case, the spectrum is also changed, and the breaking of the cubic symmetry causes there to be more levels in the rest-frame spectrum. The complications are that, similar to the moving-frames case, operators in the irreps of a smaller symmetry group must be constructed, and unlike the moving frames case, new gauge fields are required for each new value of $L^{\prime}$. This approach was applied in explicit calculation of the $\rho$ resonance in  Pelissier and Alexandru (2013). In principle one could also consider moving frames in these asymmetric boxes Lee and Alexandru (2017).

Another proposed approach is to change the boundary conditions felt by the quark fields, by introducing a “twist” that can effectively induce a non-zero value of momentum Bedaque (2004); Sachrajda and Villadoro (2005). Typically the proposal, called partial twisting, is to only apply the twist to “valence” quarks, i.e. to the propagators that go into correlator construction, so that new gauge-field generation is not required for each new twist value. Complications include the need to show that the application of a different boundary conditions to the valence versus the sea quarks does not lead to an unacceptable breakdown in the unitarity of the theory121212This concern may be somewhat less if only heavy quarks are twisted, since their annihilation and ‘sea’ behavior is likely to be less relevant. Ozaki and Sasaki (2013) used partial twisting to obtain more points on the phase-shift curve in a calculation of $J/\psi\,\phi$ elastic scattering and similarly Chen et al. (2014) in the elastic scattering of isospin=1 $D\overline{D}$.

VI.2 Resonances in coupled-channel meson-meson scattering

To date there have been four lattice QCD calculations which have used the Lüscher approach to determine coupled-channel scattering amplitudes, all in the meson sector. The first was a study of coupled $\pi K$, $\eta K$ scattering, reported on in Dudek et al. (2014); Wilson et al. (2015b), in which spectra were computed on three lattice volumes, including a large number of moving frames. The resulting energy levels were used to tightly constrain the two-channel scattering matrix, and this led to the conclusion that there is very little coupling between the channels, with the $\eta K$ channel being weakly repulsive. The $\pi K$ channel was found to feature a broad resonance and a virtual bound-state in $S$-wave, a bound-state $K^{\star}$ in $P$-wave, and hints of a narrow resonance in $D$-wave.

The second reported coupled-channel study appeared in Wilson et al. (2015a), where the extension of the spectra used to determine elastic $\pi\pi$ scattering in $P$-wave, to an energy region above $K\overline{K}$ threshold, was used to constrain coupled-channel $\pi\pi,K\overline{K}$ scattering. As in the previous case, relatively little coupling between the two channels was observed.

The third case is the only one to-date involving heavy quarks, and the first to consider three-channel scattering. Using the same distillation technology, operator constructions and analysis approach as the studies above, Moir et al. (2016) studied the coupled $D\pi,D\eta,D_{s}\overline{K}$ system, finding a near-threshold bound-state in $S$-wave, a deeply bound $D^{\star}$ in $P$-wave, and evidence for a narrow $D$-wave resonance coupled dominantly to $D\pi$.

The fourth and most recent calculation, Dudek et al. (2016), considered the coupled $\pi\eta,K\overline{K}$ system (with additional limited consideration of the coupling with the $\pi\eta^{\prime}$ channel). Experimentally, the $S$-wave in this system features a strong enhancement at $K\overline{K}$ threshold that is usually ascribed to the presence of a resonance, the $a_{0}(980)$. In Dudek et al. (2016), the Hadron Spectrum Collaboration presented spectra from three volumes in a range of frames with $m_{\pi}\sim 391$ MeV. In total, below $\pi\eta^{\prime}$ threshold, they found 47 energy levels which were used to constrain the two-channel $\pi\eta,K\overline{K}$ scattering matrix. A large number of $K$-matrix parameterizations were found to be capable of describing the finite volume spectra – Figure 17 shows all such amplitudes determined in this way which have an acceptable goodness-of-fit.

Examining Figure 17, we see that the $\pi\eta\to\pi\eta$ amplitude shows a ‘cusp-like’ behavior at the $K\overline{K}$ threshold, and indeed a cusp is always allowed when a new kinematic threshold opens. However in this case, the strength of the effect, and the rapid turn-on of amplitudes leading to $K\overline{K}$ suggest there might be a resonance nearby. This can be examined by analytically continuing the amplitude parameterizations into the complex energy plane, and on doing so, a clear outcome emerges – all successful parameterizations feature a pole on sheet IV ($\mathrm{Im}\,q^{\star}_{\pi\eta}>0,\,\mathrm{Im}\,q^{\star}_{K\bar{K}}<0$), and as can be seen in Figure 18, there is very little scatter with variation in parameterization form. The residues of $\mathcal{M}$ at the pole can be factorized leading to couplings of the resonance to the $\pi\eta$ and $K\overline{K}$ channels which prove to be comparable in magnitude 131313A description of these energy levels using amplitudes motivated by unitarized chiral perturbation theory Guo et al. (2017) gives a quite similar result for the pole position and couplings, although the validity of using the effective theory at this pion mass is not clear..

Poles are also found on sheet III, which is closest to physical scattering above the $K\overline{K}$ threshold, but they are either far into the complex plane, or lie outside the energy region in which there are energy levels to constrain the amplitude. The position of these poles, which we emphasize are irrelevant for the behavior near the $K\overline{K}$ threshold, shows significant scatter with parameterization change, suggesting that they might be artifacts.

The same study found that the $P$-wave amplitude below $\pi\eta^{\prime}$ threshold is compatible with zero, and the $D$-wave scattering matrix features a narrow resonance with comparable couplings to both $\pi\eta$ and $K\overline{K}$ final states. The extraction of the amplitudes in $D$-wave is not at the same level of rigor as the $S$-wave: the three-body $\pi\pi\pi$ channel can couple to the $J^{P}=2^{+}$ partial wave, but not to $0^{+}$, and the $D$-wave resonance lies above the $\pi\pi\pi$ threshold. The calculation reported in Dudek et al. (2016) did not include operators resembling $\pi\pi\pi$, nor was any attempt made to account for such a scattering channel, indeed the formalism to do so does not yet exist in a complete form (this will be discussed further in Section IX).

This study illustrates some important points that are likely to be generally applicable to future studies of coupled-channel resonances using the Lüscher approach. The extracted amplitudes do not have a simple ‘canonical’ resonance behavior (as illustrated in Figure 6), and it is only by obtaining a dense spectrum of energy levels across the relevant energy region (illustrated by the dots at the bottom of Figure 17) that it is possible to constrain amplitude parameterizations sufficiently to determine the energy-dependence. The rest-frame spectrum in one volume alone (unless it were extremely large) would not be enough to understand this system. The use of a range of analytic parameterizations to describe the spectra, and the ability to continue these into the complex plane to search for pole singularities, leads to a description of the resonant physics in terms of a pole position and couplings, which one can argue is the least model-dependent approach possible.

In summary, Dudek et al. (2014); Wilson et al. (2015b, a); Moir et al. (2016); Dudek et al. (2016) present an approach to studying coupled-channel scattering in lattice QCD, based upon calculating the spectrum in many moving frames using a large basis of ‘single-hadron’ and ‘hadron-hadron’ operators, that can be extended to other scattering channels. The relevant finite-volume formalism is in place for systems featuring any number of coupled two-hadron channels, including those where the scattering hadrons have non-zero spin Briceno (2014). Nucleon-meson scattering such as $\pi N,\,\eta N\ldots$, in which excited nucleon resonances, $N^{\star}$s, are expected to appear, is one obvious example. Here the challenge is largely computational – with the increase in the number of quarks comes an increase in cost – this has restricted the scope of calculations Lang et al. (2017). A larger challenge, which becomes more acute at the mass of the pion is reduced, is the presence of three-hadron and higher multiplicity channels, for which a complete finite-volume formalism is not yet in place. We will return to this issue in Section IX.

VII Other approaches to resonance determination

We have presented an approach for resonance determination in lattice QCD which follows from the Lüscher quantization condition, Eqn. 10, and illustrated its successful application in a number of explicit lattice QCD calculations. Alternative approaches have been proposed and in some cases applied, and in this section we will briefly review them, beginning with a consideration of the simplest possible approach: directly associating a particular finite-volume energy level with a resonance.

VII.1 Resonances in the Lüscher formalism in the narrow-width approximation

We might wonder if, in the case of a narrow resonance, we can use the width of the resonance as a small parameter and obtain a simple result for the appearance of the resonance in a finite volume. We can explore this in the simplest case of elastic scattering with a single partial-wave dominating. In this case the Lüscher quantization condition, Eq.(10), takes the form,

[TABLE]

where the removal of the imaginary part of $F$ and $\mathcal{M}^{-1}$ which is completely specified by unitarity has lead to the introduction of the principal-value integration.

For an amplitude dominated by a narrow elastic resonance, we may use the Breit-Wigner form, ${\cot\delta_{\ell}(E^{\star})=\frac{m_{0}^{2}-E^{\star 2}}{E^{\star}\,\Gamma(E^{\star})}}$, with an energy-dependent width ${\Gamma(E^{\star})=\frac{g^{2}}{6\pi}\frac{(q^{\star})^{2\ell+1}}{E^{\star 2}}}$ that ensures the correct threshold behavior. The corresponding quantization condition will then be,

[TABLE]

and in the limit of the width being small (either due to a small coupling, or a small phase-space for decay), we expect a solution very close to $E^{\star}=m_{0}$, and further solutions very close to the non-interacting two-body energies, $E=\omega_{k}+\omega_{Pk}$. Specifically these “nearly non-interacting” levels are the energies which solve

[TABLE]

where $\nu$ is the degeneracy of the non-interacting state. An explicit equation for the energy of the level near $m_{0}$ can be found at lowest non-trivial order in the small width of the resonance,

[TABLE]

where the pole-position mass ($m_{R}$) and width ($\Gamma_{R}$) have been introduced.

Clearly, for a small width, this is a small volume-dependent shift away from $m_{R}$, unless it happens that $E_{R}=\sqrt{m_{R}^{2}+\mathbf{P}^{2}}$ is approximately equal to a non-interacting energy level, in which case the shift can be enhanced due to the effect of the pole. This is precisely the ‘avoided level crossing’ behavior observed in Figure 5 for a Breit-Wigner resonance.

Because the sum over discrete momenta is not invariant under boosts, the volume-dependent correction term will differ in different frames, and we would not expect the energy of the resonance to perfectly satisfy the dispersion relation expected for a single particle, i.e. E_{R}(\mathbf{P},L)\neq\sqrt{\big{(}E_{R}^{\star}(L)\big{)}^{2}+\mathbf{P}^{2}}. This was explicitly observed in the results presented in Wilson et al. (2015a), which we show in Figure 13 — we see that when we boost the moving frame energy levels that lie near the resonance position back to the cm-frame, they do not all lie at the same energy, and it is this effect which allows us to map out the energy-dependence of the phase-shift within the narrow width of a resonance.

It is worth pointing out that the corresponding analysis for a bound-state leads to the finite-volume corrections being not power law as above, but rather exponentially suppressed like $e^{-\kappa L}$ with $\kappa$ the binding momentum. Because $E$ for a bound-state is always below threshold, the function in the sum above is smooth and the difference with the integral is exponentially small Luscher (1986a); Koenig and Lee (2017); Hansen and Sharpe (2017); Korber and Luu (2016); Meissner et al. (2015); Briceno et al. (2014, 2013a); Koenig et al. (2011, 2012); Bour et al. (2011); Davoudi and Savage (2011); Beane et al. (2004)

It appears that we can conclude (at least in the elastic case) that if there is a narrow resonance of mass $m_{R}$, there will be a finite-volume energy level close to $m_{R}$, with a volume-dependent shift that is likely to be small unless a non-interacting level happens to lie very close by. In the limit of zero width due to there being no coupling to the decay channel, we would interpret this state as being like the $|\rho,L\rangle_{0}$ basis state we discussed in Section V.2.1, that will be interpolated well by ‘single-hadron-like’ operators. We can thus combine our observations in this section and those in Section V.2.1 to propose an explanation of the spectra extracted in calculations like Dudek et al. (2013a), which use only a large basis of ‘single-hadron-like’ operators and no multihadron operators, one example of which we previously presented in Figure 8. It is likely that these calculations are resolving the presence of relatively narrow resonances in the relevant energy regions, and the masses plotted should be viewed as being imprecise guides to the $m_{R}$ values of these resonances. What the decay properties of these states are, or whether there are broader resonances present, cannot be determined from such calculations.

VII.2 Resonances and ‘naive’ level counting

The previous section suggested that the presence of a narrow resonance might be inferred by there being an ‘extra’ energy level beyond those expected in the non-interacting spectrum. A number of recent calculations  Prelovsek and Leskovec (2013b); Prelovsek et al. (2015); Lee et al. (2014); Padmanath et al. (2015); Lang et al. (2015, 2016, 2017) make use of bases of operators that include multihadron operators and extract spectra that, in some cases, are likely to be reasonably close to the complete spectrum in a given energy region. Typically only the rest frame is considered, and the use of relatively small volumes limits the density of states extracted. Since these spectra commonly extend into energy regions where more than one scattering channel is kinematically open, a coupled-channel Lüscher analysis of the type described in Section VI.2 would be desired, yet there are insufficient numbers of levels to fully constrain multichannel parameterizations. Instead an approach is followed that appeals to narrow resonance arguments of the type we presented in the previous section, and as we have argued, at best this approach can suggest the presence of a narrow resonance, but cannot determine its properties with any precision.

The calculated finite-volume spectrum is compared to the known multichannel non-interacting spectrum, with each calculated level being associated with a non-interacting level, until all non-interacting levels in the relevant energy region have been exhausted. Any calculated level that appears in excess of this counting is thus a candidate to be due to the presence of a narrow resonance. Further confirmation of these assignments are sought from the values of the overlaps $\langle E_{n}|\mathcal{O}|0\rangle$ extracted from the variational solution, with large overlaps onto ‘single-hadron’ operators used to suggest a resonance assignment.

An example, albeit not an ideal one, owing to the complicating presence of scattering hadrons of non-zero spin, and three-body channels which are not considered, is the analysis presented in Padmanath et al. (2015). Using a basis of meson-meson operators and single-hadron operators of both quark-bilinear and tetraquark constructions, they extract a rest-frame spectrum in the $T_{1}^{+}$ irrep (where the $X(3872)$ experimental state would be expected to appear if it is a $J^{P}=1^{+}$ resonance) in a single volume. Observing an ‘extra’ level beyond the $\eta_{c}$‘$\sigma$’, $D\overline{D}^{*}$, $J/\psi\,\omega$ and $\chi_{c1}$‘$\sigma$’ levels expected in a narrow energy window, they choose to associate this state with the $X(3872)$.

The converse argument has also been made, that if no level appears in excess of those expected in the non-interacting spectrum, there is likely to be no resonance. An example is the calculation presented in Prelovsek et al. (2015) of the hidden-charm $I=1$ sector, where a basis of meson-meson operators is supplemented with tetraquark operators, and the rest-frame spectrum of the $T_{1}^{+}$ irrep is determined in a single volume. Each calculated lattice level is matched with an expected nearby non-interacting level, with no lattice levels left over, and it is concluded that there is no signal for a $Z_{c}$ resonance.

We conclude this section by raising a note of caution regarding the use of the absence of an ‘extra’ level to conclude the absence of resonant behavior. In Section IV.3.3 we presented an example of a two-channel amplitude which features a relatively narrow resonance corresponding to a pole singularity on the second sheet, but whose finite-volume spectrum does not obviously feature an isolated energy level that we can associate with the resonance (Fig. 7). In Figure 19 we present an example of how the approach described in the current section could be misleading when applied to this amplitude — a calculation of the rest-frame spectrum in a quite reasonable $(2.4\,\mathrm{fm})^{3}$ volume would give the spectrum shown, which would certainly not lead one to conclude the presence of a relatively narrow resonance close to the $K\overline{K}$ threshold and might, unfortunately, lead one to incorrectly conclude the absence of any resonance in this energy region.

In summary, the level counting approach is based largely on expectations from elastic scattering, and in the case of coupled-channels, while we might get simple narrow resonances which do behave like elastic resonances (e.g. Fig. 6), and which come with an ‘extra’ level, we can also have more subtle situations where a resonance is present, but no ‘extra’ level appears. The approach described in Sec. VI.2, where the complete spectrum in many moving frames and/or volumes is used to tightly constrain the energy dependence of the scattering matrix, does not suffer from the weaknesses outlined above.

VII.3 Finite volume EFT Hamiltonian approach

Following the diagrammatic approach presented in Section IV.2, given an explicit quantum field theory Lagrangian describing the interactions between hadron fields, one can proceed to determine and parametrize the propagators of multi-particle states in a finite-volume, and the location of the poles of these objects would give the prediction of the finite-volume spectrum for that given QFT.

The procedure we have been describing thus far is one where a particular hadron-level Lagrangian is never specified, but instead we use relations that hold for any unitary quantum field theory, which relate the finite-volume spectrum to the infinite-volume $S$-matrix.

Alternatively, one can specify a particular choice of hadron-level lagrangian, consider the finite-volume propagator, say of two-particle states, and tune the parameters of the theory in order to reproduce the observed finite-volume spectrum Bernard et al. (2008b); Beane et al. (2012a); Hall et al. (2013); Wu et al. (2014). Because this procedure assumes a suitable low-energy effective field theory that describes the interacting system in the region of interest being evaluated in a finite-volume, it is commonly referred as the finite volume EFT Hamiltonian approach.

Although the Lüscher quantization condition is never explicitly used in this approach, we can show that the finite-volume spectra of the EFT Hamiltonian have the same relationship with the infinite volume scattering amplitudes as given by the Lüscher quantization condition (up to exponentially suppressed corrections) — provided one works within the validity of the EFT and incorporates enough non-negligible partial waves 141414 A similar technique to those presented in the references above was used in Briceno et al. (2013b) to generalize the Lüscher condition to two-nucleon systems . We can illustrate the equivalence of the two approaches in the simple case of elastic scattering of equal mass spinless bosons in a single dominant partial-wave, in the overall rest-frame Wu et al. (2014). The development of the Lüscher approach presented in Section VII.1 leads to

[TABLE]

while the equivalent expression in the finite-volume EFT approach Hall et al. (2013) follows from the relation determining the relative momentum of the two-body system,

[TABLE]

via

[TABLE]

In these expressions, $\Delta_{0}$ and $g$ are bare parameters or functions in the EFT being considered and the self-energy, $\Sigma(q^{\star})=2\pi^{2}\int\!\frac{d\mathbf{k}}{(2\pi)^{3}}\,\frac{4\pi}{k^{2}-q^{\star 2}}\frac{g^{2}(k)}{\omega_{q}-\omega_{k}+i\epsilon}$. Eliminating the imaginary part that is fixed by unitarity, we are left with

[TABLE]

and if $g(k)$ is a smooth function of the momentum, then the difference between the sum and the integral will be an exponentially small quantity such that

[TABLE]

which is clearly equivalent to Eq.(26) for reasonably large volumes.

Given that the finite volume spectra of the Hamiltonian and the Lüscher formalism have the same relationship to on-shell scattering, neither method offers any more or less information about the physical scattering amplitude. The EFT Hamiltonian can be viewed as a particular amplitude parameterization. When constrained by the same lattice QCD results, all sufficiently flexible parameterizations should agree on the on-shell scattering amplitude. The alternative parameterizations will in general give different analytic continuation into the complex plane and this may manifest itself in variation in the position and couplings of any resonance poles.

The Hamiltonian formulation has been used in the analysis of recent numerical lattice QCD studies, with the study of the $\Lambda$(1405) generating significant attention Hall et al. (2015); Liu et al. (2017a). Owing to the computational cost, the numerical studies of this state have been limited to the use of only local interpolating operators Menadue et al. (2012); Hall et al. (2015). At heavy quark masses ($m_{\pi}\gtrsim 400\,{\rm MeV}$) the $\Lambda$ is clearly isolated as a deeply-bound state, far below the $\pi\Sigma$ threshold. As the quark masses are lowered, this state becomes resonant in $\pi\Sigma$ and lies somewhere near the vicinity of the $\overline{K}N$ threshold. Without having access to a dense finite-volume spectrum, it has not been possible to constrain a general parameterization of the coupled-channel scattering as described in Section VI.2. A more phenomenological approach is taken by introducing a Hamiltonian that describes the meson-baryon interactions through a “bare” state. The parameters of this Hamiltonian are then constrained by the quark-mass evolution from the bound-state region at larger quark masses, supplemented with a coupling strength informed by the physically-observed width.

This work has gone further to investigate the electromagnetic coupling to this finite-volume state. The claimed result is that the strangeness magnetic moment of the $\Lambda$(1405) vanishes in the vicinity of the physical point Hall et al. (2015). Such a result has highlighted the tremendous discriminating power of external probes to resolve the structure of resonant states. In this case, a vanishing strangeness moment would suggest that the strange quark is bound to a $\overline{K}$ meson, which cannot carry any magnetization in an $S$-wave. As discussed above, even in the narrow-width approximation, the quantization condition leads to finite-volume eigenstates which boost nontrivially. Finite-volume eigenstates in the $\Lambda$ system could then differ significantly in composition in boosted frames. The influence of this modified structure on three-point correlation functions remains to be quantified. A systematic approach to matrix elements will be discussed further in Section VIII.

The Hamiltonian formulation has also been used in further studies of the $N^{*}$(1535) Liu et al. (2016) and Roper Liu et al. (2017b); Kiratidis et al. (2017); Wu et al. (2017) resonances. In these exploratory investigations, two-body operators have not been used in the extraction of the spectrum. As discussed in Sec. V.2.1, this implies that the spectrum obtained is incomplete and one cannot trust that the spectra obtained are indeed accurate, which would subsequently the resulting scattering amplitudes and resonant poles.

VII.4 Unitarized chiral perturbation theory and chiral extrapolations

An approach closely related to that described in the previous section considers the particular EFT of chiral perturbation theory, unitarized using some prescription Oller et al. (1998); Dobado and Pelaez (1997); Oller et al. (1999); Gomez Nicola and Pelaez (2002); Pelaez and Rios (2006), in a finite-volume Agadjanov et al. (2016b); Chen and Oset (2013); Doring et al. (2012); Doring and Meissner (2012); Molina and Doering (2016); Doring et al. (2011); Garzon et al. (2014); Bernard et al. (2011); Nebreda et al. (2011); Rios et al. (2008); Guo et al. (2009); Hanhart et al. (2008). In practice, these works rederive or generalize the Lüscher spectral condition using $\chi$PT as a starting point. The final result always has the same functional form as Eq.(10) (up to exponentially suppressed corrections), with the scattering amplitude being that of $\chi$PT or some suitable unitarization (U$\chi$PT ). The advantage of this procedure is that it typically builds into the amplitudes some low-energy properties of QCD (relevant for light quarks), and it allows for the relation of observables obtained using different values of the quark masses, a “chiral extrapolation”.

This methodology has been applied in the study of non-resonant scattering lengths in the heavy-light meson sector Liu et al. (2013), as well as in chiral extrapolations of the $\rho$ Bolton et al. (2016); Hu et al. (2016). First exploratory studies of the chiral extrapolations of the $a_{0}(980)$ Guo et al. (2017) and $\sigma$ Doering et al. (2016) resonances have appeared, analyzing the finite-volume energy levels calculated by the Hadron Spectrum Collaboration Dudek et al. (2016); Briceno et al. (2017a).

From a technical standpoint there is no new insight coming from considering U$\chi$PT in a finite volume as compared to simply using Eq.(10) with amplitude forms given by infinite-volume U$\chi$PT. In the case of multichannel scattering there is not a unique procedure for unitarizing $\chi$PT, so one should be careful not to assume that the corresponding amplitudes are a direct consequence of QCD. Pioneering lattice QCD calculations of the excited state spectrum are currently being performed with relatively heavy quark masses, $m_{\pi}\gtrsim 300$ MeV, where the reliability of $SU(2)$ flavor $\chi$PT would be questionable, and furthermore, the application of $\chi$PT for channels with explicit strange-quark degrees of freedom requires making an expansion about the $SU(3)$ flavor symmetric point of QCD, and obtaining reliable estimates associated with this poor approximation of QCD is notoriously challenging.

The use of chiral effective field theories to extrapolate hadron results down in quark mass is not a new approach. In the past, when states which appear as resonances at the physical quark mass were only realized as bound states at higher quark masses, this approach was used to extrapolated lattice results for the $\rho$ Leinweber et al. (2001); Bruns and Meissner (2005); Allton et al. (2005); Armour et al. (2006) and $\Delta$ Leinweber et al. (2000); Young et al. (2002); Bernard et al. (2005); Pascalutsa and Vanderhaeghen (2006) resonances, for instance. In those works, it was necessary to put in the decay thresholds artificially, using the experimental width to constrain the relevant coupling parameters. The significant new feature moving into the true resonance regime is the ability to directly determine the width (and hence coupling) within the lattice calculation and thereby directly study the quark mass dependence of the complex pole location — such as highlighted in Figure 1.

VII.5 Other approaches

There have been several other approaches proposed to determine resonance properties in lattice QCD calculations, here we will briefly summarize them.

In what we might call the ‘tuned threshold’ or UKQCD method (see e.g. McNeile and Michael (2006)), the approach requires us to find a value of quark mass such that the resonance to be studied happens to lie rather close to the kinematic threshold of its two-body decay channel. That is, if the resonance decay is $R\to ab$, $m_{R}\approx m_{a}+m_{b}$. Under these conditions it is argued that the resonance to decay-channel coupling, $g_{R\to ab}$, can be determined by considering the correlation function with a ‘single-hadron’ operator with the quantum numbers of $R$ at the source, and an $ab$-like ‘meson-meson’ operator at the sink. The finite volume of the lattice appears to play no role in this approach. It was applied to the decays $b_{1}\to\pi\omega$ and $\pi_{1}\to\pi b_{1}$ in McNeile and Michael (2006), and has recently been applied to baryon decays in Alexandrou et al. (2013, 2016).

Underlying this framework are two necessary assumptions which might lead to important systematic errors. First, it is assumed that states interpolated using single-particle and multi-particle operators correspond solely to a resonance, $|R\rangle$, and a scattering state, $|ab\rangle$, respectively. As we saw in Sec. V.2.1, this in general not true. Second, one must assume that the resonance state is an asymptotic state, and the overlap between the multi-particle and resonance state, $\langle R|ab\rangle$, is equal to the scattering amplitude coupling these two “channels”. This, of course, is not true for actual resonances, but the process of tuning the threshold such that resonance cannot actually decay, making it a bound-state at threshold, increases the possible validity. Of course in practice it may not be possible to find a quark mass where the required near-degeneracy occurs, or that quark mass may be very far from the physical quark mass, making the extrapolation questionable.

The histogram method is based on the construction of a probability distribution which mimics the scattering cross section Bernard et al. (2008a). This method utilizes the low-lying eigenstate spectrum across a uniform distribution of finite volumes, binned as a function of energy, one then generates a histogram which counts the number of states. In the vicinity of a resonance there will be a pronounced enhancement in this “density of states” distribution. While offering a clear visualization of the resonance structure, a direct comparison has demonstrated that it is not as straightforward to implement as the conventional Lüscher technique Giudice et al. (2012), and it is not clear that the approach would be useful for the study of resonances which do not appear as a clear ‘bump’ in the scattering amplitude.

A recently proposed approach to coupled-channel scattering in a finite-volume makes use of an optical potential to convert the problem into one featuring just a single channel with the effect of the other channels absorbed into the optical potential Agadjanov et al. (2016b). The intention is to extend the power of the Lüscher method to high energies where many-particle states can go on-shell without having to explicitly take into account such effects (by hiding them in the optical potential), and although this appears to be a formally sound idea, putting this into practice even for the simplest of systems introduces some model-dependence. At this early stage it is not clear how much information must be supplied by the finite-volume spectrum to make this approach practical – in the coupled $\pi\eta$, $K\overline{K}$ case considered for illustration in Agadjanov et al. (2016b), the density of states required was extremely high, at the level where, if one actually had this set of energy levels, one could bin in energy, and in each bin solve coupled versions of Eq.(10), to find all elements of the scattering matrix without needing to introduce any auxiliary optical potential.

Lastly we consider the potential method, which has been championed by the HAL QCD collaboration Ishii et al. (2007); Aoki et al. (2010); Ishii et al. (2012); Aoki et al. (2012a, 2013). This is the only approach which does not make direct use of the discrete spectrum in a finite-volume. Instead, from correlation functions constructed using operators featuring spatially displaced hadron operators, one accesses a quantity which the authors associate with a non-relativistic inter-hadron potential. By solving a Schrödinger equation including this potential, one obtains scattering amplitudes and/or bound-state energies.

This formalism has largely been implemented in the study of baryon-baryon systems, where unlike the results obtained by other groups using fairly standard determinations of sub-threshold energy levels Berkowitz et al. (2017); Beane et al. (2012b, 2013); Yamazaki et al. (2012, 2015), the HAL QCD collaboration sees no evidence for a bound deuteron or di-neutron for heavier than physical quark masses. At this stage it is not clear what is the source of this discrepancy, whether it is that the long-distance part of the potential, to which shallow bound-states are most sensitive is not being accurately determined by HAL QCD, or whether, as HAL QCD have suggested Iritani et al. (2016, 2017), that the spectra obtained in Berkowitz et al. (2017); Beane et al. (2012b, 2013); Yamazaki et al. (2012, 2015), which do not make full use of the variational technique presented in Section V.1, are not reliable. 151515We point the reader to a response by the NPLQCD collaboration Beane et al. (2017) in regards to the remarks made by the HAL QCD collaboration Iritani et al. (2017). Another application of the potential method by the HAL QCD collaboration has been in a multi-channel study of the $Z_{c}(3900)$ Ikeda et al. (2016).

An important observation is that while systems potentially featuring shallow bound-states (baryon-baryon), and complicated coupled-channel meson systems ($Z_{c}$) have been studied in the potential approach, the simplest resonant systems have not. In particular, there has been no presentation of $\pi\pi$ $I=1$ scattering and the $\rho$ resonance, which would seem to be an ideal benchmark calculation. The sensitivity of the tightly-bound $\rho$ to the very short-distance part of a $\pi\pi$ ‘potential’ might be a problem for the method, as it is not clear that the short-distance region can be precisely determined.

VIII Coupling resonances to external currents

As we have seen already in this review, the field has witnessed a tremendous burst of activity providing progress in the ability to determine masses, widths and couplings to two-body hadronic decay channels of hadronic resonances. In particular, a general quantization condition for two-particle systems in a finite volume has been derived, Eq.(10), which relates the discrete spectrum in a finite volume to scattering amplitudes in infinite volume, and approaches based upon Eq.(10) have been successfully implemented in a number of numerical calculations.

Going beyond this, we can consider if there are further rigorous tools we can bring to bear to aid development of an understanding of hadron resonances. The response of resonances to external probes, such as electromagnetic or weak currents, has the potential to reveal substantial new information about the internal structure of resonant systems. For instance, as discussed above, it has been demonstrated that a vanishing strange-quark magnetization of the $\Lambda(1405)$ would suggest that the strange quark is most-probably bound in a kaon and hence a clear signature that the resonance has the structure of a $\bar{K}N$ molecule Hall et al. (2015). Similarly, the nature of the Roper resonance has challenged theorists for many years, and resolving the relevant electromagnetic transition form factors Wilson et al. (2012); Segovia et al. (2015), which can also be measured Aznauryan et al. (2008); Dugger et al. (2009); Aznauryan et al. (2009); Aznauryan and Burkert (2012), could help to distinguish between different proposed scenarios, such as a quark orbital excitation, explicit gluonic excitation, or hadronic molecule. The gluonic content of resonances could also be studied directly by evaluation of gluonic operators which have been previously considered for stable hadrons Meyer and Negele (2008); Horsley et al. (2012); Detmold and Shanahan (2016). Such matrix elements would provide a more rigorous framework with which to build on existing analyses where gluonic content has been inferred from operator overlaps Dudek et al. (2013a); Liu et al. (2012); Dudek and Edwards (2012); Dudek (2011); Edwards et al. (2011); Dudek et al. (2011a, 2010, 2009).

To isolate matrix elements of resonant states in lattice QCD, we once again must apply a formalism which maps between the quantities directly obtained from finite-volume correlation functions and the desired infinite-volume observables.

Firstly we will present a diagrammatic representation for an infinite-volume scattering amplitude featuring external currents, which we will refer to as a transition amplitude. In close analogy to scattering amplitudes, the diagrammatic representation of the transition amplitude describing $n$ incoming and $n^{\prime}$ outgoing hadrons, with $j$ insertions of external currents, $\mathcal{T}_{j}(n\to n^{\prime})$, can be expressed as

the sum over all diagrams with $n$ incoming and $n^{\prime}$ outgoing hadronic legs (that have been amputated and put on-shell) and $j$ insertions of the external current. The momenta flowing through the external currents are constrained to satisfy conservation of momentum, but are otherwise free. All intermediate hadron propagators are evaluated with the $i\epsilon$-prescription and all intermediate loop momenta are integrated.

Figure 20 illustrates this in the case of (a) a $0\to 2$ transition, relevant to e.g. the vector decay constant of the $\rho$ resonance decaying to $\pi\pi$ and (b) a $1\to 2$ transition, relevant to e.g. the transition form factor of the $\Delta$ baryon resonance in $\gamma^{\star}N\to\pi N$. The case of $2\to 2$ transitions, relevant to e.g. deuteron photo-disintegration $\gamma d\to np$ where the deuteron is a bound-state in $np$ scattering, can be found in Briceno and Hansen (2016).

Focussing on the $1\to 2$ transition, it is clear that $i\mathcal{H}$ should depend on the virtuality of the external current, $Q^{2}$. For a fixed value of $Q^{2}$, with the three external hadrons on-shell, the only remaining freedom is in the direction of the relative momentum of the two outgoing hadrons, and it follows that in analogy to the partial-wave expansion of the scattering amplitude, we can partial-wave expand the transition amplitude,

[TABLE]

where for simplicity we have assumed the external hadrons are all spinless, although the generalization is straightforward. The energy dependence of the partial-wave transition amplitude, $\mathcal{H}_{\ell m}(Q^{2},E^{\star})$, must be closely related to that of the corresponding scattering amplitude, $\mathcal{M}_{\ell}(E^{\star})$. For example if the transition process is $\gamma^{\star}\pi\to\pi\pi$ in $P$-wave, we would expect a close relationship with the scattering process $\pi\pi\to\pi\pi$ in $P$-wave. Indeed this relation can be seen by comparing the diagrams in Figure 20 with those in Figure 2, where we observe that the iterated ‘rescattering’ is the same in both cases, with the only difference being the presence of an additional kernel characterizing the coupling to the external current. In particular we would expect that if the strong rescattering in $\mathcal{M}_{\ell}$ gives rise to a pole in $E^{\star}$ corresponding to a composite particle (e.g. a bound-state or a resonance), that same pole singularity must also be present in $\mathcal{H}_{\ell m}$.

In fact we can rigorously define what we mean by a resonance transition form-factor by considering the behavior of $\mathcal{H}$ near the pole. Recall that for the scattering amplitude, we could factorize at the pole to define couplings to the incoming and outgoing channels, e.g.

[TABLE]

similarly we can factorize $\mathcal{H}$,

[TABLE]

where $F_{\gamma e\to R}(Q^{2})$ is the transition form-factor for $\gamma e\to R$ as a function of the virtuality of $\gamma$. In these expressions $R$ may be a bound-state, in which case $s_{0}$ is real and below kinematic threshold, and $g,F$ are (in a suitable convention) purely real, or $R$ may be a resonance, in which case $s_{0}$ is complex, and $g,F$ may also be complex.

VIII.1 Determining matrix elements in lattice QCD

Retaining our focus on amplitudes like $\mathcal{H}$ which feature hadrons in both incoming and outgoing states, with insertion of a single current, we are led to consider three-point functions of the following form,

[TABLE]

where $\mathcal{B}^{\dagger}$ has the quantum numbers of the desired incoming state, $\mathcal{A}$ has the quantum numbers of the desired outgoing state, and $\mathcal{J}$ is the current operator. As in the two-point function case, there is a dispersive representation in terms of the discrete eigenstates of QCD in a finite-volume,

[TABLE]

where the finite-volume transition matrix elements, \big{\langle}E_{n^{\prime}},\mathbf{P}^{\prime},L\big{|}\mathcal{J}(0)\big{|}E_{n},\mathbf{P},L\big{\rangle} are what we would like to determine in explicit calculations.

In the simplest case where the initial and final states are both QCD-stable hadrons, the finite-volume transition matrix elements have an immediate interpretation – they can be related to the infinite volume transition amplitudes with only exponentially suppressed finite-volume corrections. For example if $\mathcal{A}$ and $\mathcal{B}$ are both isospin-1 pseudoscalar operators and the current has vector quantum numbers 161616strictly speaking, operators in cubic or little-group irreps that contain subductions of the relevant helicities

, the lowest energy contributions to the correlator comes from the pion form-factor, with the matrix-element having a Lorentz covariant decomposition (in infinite volume 171717in finite volume there is an extra factor of \big{[}(2\omega L^{3})(2\omega^{\prime}L^{3})\big{]}^{-1/2} due to the different normalization of finite-volume states. ),

[TABLE]

where the virtuality $Q^{2}=-(P^{\prime}-P)^{2}$. Current conservation constrains the value of $F_{\pi}(0)$, but the form-factor is otherwise unconstrained, and non-perturbative calculations are required to determine the $Q^{2}$ dependence. Provided volumes such that $m_{\pi}L\gg 1$ are used, the exponentially suppressed finite-volume corrections can be neglected, and such stable-hadron form-factors can be computed in a relatively straightforward manner in lattice QCD. Well studied cases include the $\pi$ and nucleon form factors Green et al. (2015); Shultz et al. (2015); Chambers et al. (2017); Collins et al. (2011); Alexandrou et al. (2011); Capitani et al. (2015); Lin et al. (2010); Yoon et al. (2016).

Our interest here though is in composite particles, and in this case the relation between the finite-volume transition matrix element and the infinite-volume transition amplitude is not so simple, and the formalism to relate the two in general cases has only very recently been laid out.

VIII.2 Lellouch-Lüscher formalism and its generalizations

The first suggestion for how to relate transitions in infinite volume to the matrix elements calculated in a finite volume came from Lellouch and Lüscher, who described the case of the weak decay $K\to\pi\pi$ Lellouch and Luscher (2001) in the rest frame of the kaon. This idea has subsequently been generalized to systems in moving frames  Kim et al. (2005); Christ et al. (2005), coupled channels Hansen and Sharpe (2012), $\gamma^{\star}\to\pi\pi$ Meyer (2011); Feng et al. (2015), bound state photodisintegration Meyer (2012), $\gamma^{\star}N\to\Delta$ Agadjanov et al. (2014), and elastic resonance form factors Bernard et al. (2012) among others. In this subsection we will sketch the framework as presented in Briceno and Hansen (2015); Briceno et al. (2015b); Briceno and Hansen (2016) which holds for generic two-body systems. 181818 Agadjanov et al. (2016a) presented a rederivation of the general result presented in Briceno et al. (2015b) applied to the specific $B\to K^{\star}\to K\pi$ weak decays.

Briceno et al. (2015b) showed that the relationship between finite-volume matrix elements and infinite-volume transition amplitudes can be obtained using a diagrammatic representation of three-point functions, following a procedure similar to the one described in Sec. II.3 for two-point functions191919an earlier study using similar methods, Bernard et al. (2012), was restricted to the non-relativistic case.

As an example of this formalism, let us consider the three-point correlation functions with a single QCD-stable hadron interpolated by $\mathcal{A}$ and a two-hadron state interpolated by $\mathcal{B}^{\dagger}$ (relevant to e.g. $\gamma\pi\to\pi\pi$). In momentum space, this is diagrammatically depicted in Fig. 21. As was the case for two-point functions, we can express the correlation function in terms of on-shell infinite-volume quantities (in particular, the transition amplitudes $\mathcal{H}^{\rm in}(Q^{2},P)$), and finite-volume quantities,

[TABLE]

where $\Delta$ is the single-hadron propagator. One can then Fourier transform to Euclidean space-time, and match the resulting expression to Eq.(30). After a substantial amount of algebra, one arrives at the following expression relating the finite-volume matrix element of the external current to infinite-volume transition amplitudes,

[TABLE]

where the matrix $\mathcal{R}$ was defined in Eq.(IV.2). $\mathcal{H}^{\mathrm{in}}$ and $\mathcal{H}^{\mathrm{out}}$ are row and column vectors in the same space in which $\mathcal{R}$ resides, namely the space of angular momenta and open channels, and each element corresponds to the angular momentum component of the amplitude associated with the given channel. For example, when only one channel is kinematically open, we can explicitly write the product appearing in the equation above as,

[TABLE]

where prior to partial wave projection, $\mathcal{H}$ is defined diagrammatically in Fig. 20, and its partial wave projection was defined in Eq.(27).

This result demands some commentary: firstly, we observe that there is not a one-to-one mapping between the finite-volume matrix element and the infinite-volume transition amplitudes. In general, the reduction of rotational symmetry tells us that multiple angular momentum transition amplitudes will appear in a matrix element computed for a given lattice irrep (see e.g. Table 1 for the $\mathbf{P}=\mathbf{0}$ case), and this is implicit in Eq.(VIII.2). Secondly, we see from Eq.(IV.2) that $\mathcal{R}$ depends on the volume, on two-body kinematics, and importantly, also the two-body scattering dynamics in $\mathcal{M}$. As a result, in order to correctly extract transition amplitudes from finite-volume matrix elements, one must first determine the two-body spectrum and then use Eq.(10) to constrain the relevant scattering amplitudes. Furthermore, since $\mathcal{R}$ depends not just on $\mathcal{M}$, but also its derivative, we require a parameterization of the energy-dependence of the scattering amplitude, which can in principle introduce a source of systematic uncertainty.

Although we have sketched a case in which the asymptotic hadrons are spinless, the extension to hadrons with spin is straightforward Briceno and Hansen (2015). The related formalism for $0\to 2$ transitions was first introduced by Meyer in the context of extracting the timelike $\pi$ form factor via $\gamma^{\star}\to\pi\pi$ Meyer (2011), and the generalization to make it applicable to generic $0\to 2$ processes is presented in Briceno and Hansen (2015).

VIII.3 Applications

The primary application of the finite-volume formalism to date has been to the process $K\to\pi\pi$ where the weak interaction induces a hadronic decay of the otherwise QCD-stable kaon. The weak interaction does not conserve isospin, so both isospin=0 and isospin=2 $\pi\pi$ states can appear in this decay. Experimentally the $I=0$ channel is observed to be an order of magnitude larger than the $I=2$ channel, and this can be qualitatively explained by the presence of the resonant $\sigma$ in $I=0$ and the absence of resonant behavior in $I=2$. The role these amplitudes play in constraining sources of CP violation in the Standard Model has ignited a significant research program to determine the $K\to\pi\pi$ amplitudes precisely using lattice QCD Bai et al. (2015); Blum et al. (2012a); Boyle et al. (2013); Blum et al. (2011, 2012b).

Kinematically, the decay $K\to\pi\pi$ decay is rather simple – the cm-energy of the two-pion state must coincide with the mass of the kaon, which is possible since the pion momenta lie in a continuum in infinite volume. However in lattice QCD calculations one is not at liberty to specify the energies of $\pi\pi$ states in a given volume, since these are governed by the dynamics of the system, as we have previously discussed. What one can determine is the transition amplitude at nonzero virtualities, here defined as $Q^{2}=-(P_{K}-P_{\pi\pi})^{2}$, and interpolate to the physically relevant point, $Q^{2}=0$ and $E_{\pi\pi}=E_{K}$. Technical challenges have to-date prevented this approach from being applied fully, instead what has been done mostly is to implement carefully chosen boundary conditions (as described at the end of Section VI.1) to tune a $\pi\pi$ energy level such that it coincides with the kaon mass. Furthermore, in all of these studies a single operator is used to extract the finite-volume $\pi\pi$ state, which we have extensively discussed can lead to large uncontrolled systematic errors in the spectrum, matrix elements, and amplitudes.

For $0\to 2$ transitions, the case of $\gamma^{\star}\to\pi\pi$ has been recently explored Feng et al. (2015), partially motivated by the important role this amplitude plays in constraining the dominant hadronic contribution to anomalous magnetic moment of the muon. The infinite-volume transition amplitude is a Lorentz vector whose dynamics can be parametrized by a single dimensionless function which corresponds to the pion form-factor, Eq.(32), for time-like virtualities, $Q^{2}=-(E^{\star})^{2}<0$. By restricting attention to energies below the $K\overline{K}$ threshold, Feng et al. Feng et al. (2015) determined $F_{\pi}(E^{\star})$ at two values of the quark masses corresponding to $m_{\pi}=290,\,380$ MeV, shown in Figure 22. One clearly observes the presence of the $\rho$ resonance in the energy dependence of the form-factor, as we would expect since $\gamma^{\star}\to\pi\pi$ should have the same pole singularities as $\pi\pi\to\pi\pi$ in $P$-wave.

The prototype amplitude to illustrate $1\to 2$ scattering is $\pi\gamma^{\star}\to\pi\pi$, which plays a role in a wide range of phenomenology, principally through the isospin=1 $P$-wave amplitude, which features the $\rho$ resonance (see, for example, Colangelo et al. (2014b, a); Wess and Zumino (1971); Witten (1983); Huston et al. (1986); Capraro et al. (1987)). This amplitude was recently determined for the first time using lattice QCD in Briceno et al. (2015a, 2016) using the correlator construction technology laid out in Shultz et al. (2015). With a pion mass of $391$ MeV, a large number of three-point correlations functions were computed corresponding to a range of momenta for the $\pi$ operator, the $\pi\pi$ operator, and the current insertion. The $\pi\pi$ operators used were those found to optimally interpolate finite-volume energy levels through variational analysis in Dudek et al. (2013b) – they contain linear combinations of $\pi\pi$-like constructions and $\rho$-like ‘single-meson’ operators. From these correlation functions, finite-volume matrix elements of the type appearing on the left-hand-side of Eq.(VIII.2) were determined, and the corresponding infinite-volume amplitude $\mathcal{H}$ could then be obtained using parameterizations of the scattering amplitude capable of describing the finite-volume spectrum. A global analysis, making use of parameterization of the $E^{\star}$ and $Q^{2}$ dependence of the transition amplitude was performed.

The results for the transition amplitude as a function of $E^{\star}$ at two different values of $Q^{2}$ are shown in the left panel of Figure 23, where the enhancement due to the $\rho$ resonance is clearly visible. In the right panel of Figure 23 we show the $\rho\to\pi\gamma$ transition form-factor rigorously extracted from the residue of the $\rho$-pole in our parameterized transition amplitude as indicated in Eq.(29). The smallness of the imaginary part of $F$ is correlated with the narrow width of the $\rho$ at $m_{\pi}\sim 391$ MeV.

It is worth remarking the rapid progress there has been in this area of the field. In the last couple of years, we have seen the ideas first set in place by Lellouch and Lüscher Lellouch and Luscher (2001) generalized to increasingly complex reactions Briceno and Hansen (2015); Briceno et al. (2015b); Briceno and Hansen (2016). Furthermore, we are entering an era where the sophistication of these ideas are matched by that of the numerical tools. Two illustrated examples are the two discussed in this section. First, the idea of studying the $\gamma^{\star}\to\pi\pi$ amplitude was proposed by Meyer by considering a system of $\pi\pi$ at rest Meyer (2011). This was shortly after generalized by Feng et al. Feng et al. (2015) for $\pi\pi$ systems with nonzero momenta, and in the same study the authors performed the first numerical implementation. Similarly, the idea of studying $1\to 2$ reactions was first proposed in 2014 by two groups Agadjanov et al. (2014); Briceno et al. (2015b), and within a year these ideas were first put into practice in the context of the $\gamma^{\star}\pi\to\pi\pi$ amplitude Briceno et al. (2015a, 2016).

With the formalism now laid out and tested in the simplest of systems, one can expect to see application to a variety of more interesting systems. A contemporary example comes in the charmonium sector, where $XYZ$ resonances defy explanation in simple models. The vector $Y$ resonances can have their decay constants studied rigorously using the approach presented above, and the relative size of these can be compared to the production rates in $e^{+}e^{-}$ and may also inform suggestions as to the internal structure of these states. Radiative transitions between various $XYZ$ resonances and conventional charmonium states can also be studied rigorously, and the corresponding form-factors can similarly be used to infer details of internal structure.

IX Contemporary Extensions

In this document we have focused our attention on ideas that have attained a certain maturity and that have been tested in at least one numerical study. In this section we will discuss a selection of extensions that have not yet quite reached this level.

IX.1 Particles with nonzero intrinsic spin

The bulk of serious calculations to-date have considered meson-meson scattering systems where the scattering particles are spinless. The formalism to deal with stable scattering particles of non-zero spin is in place Briceno (2014), and with it one could consider systems featuring mesons like the $\omega$ or the $D^{\star}$, which are stable for pion masses larger than physical, or the $J/\psi$ which only has a small hadronic width and which can be rendered stable by neglecting $c\bar{c}$ annihilation. The formalism would also have an application in meson-baryon and baryon-baryon scattering.

Scattering systems where the hadrons are not spinless are typically more complicated, featuring coupled amplitudes even if only one kinematic channel is open. For example, in the case of pseudoscalar-vector scattering, the $S$-wave has $J^{P}=1^{+}$ quantum numbers (${}^{3}\!S_{1}$ in the ${}^{2S+1}L_{J}$ notation), but so does the $D$-wave (${}^{3}\!D_{1}$) and these channels will mix in a $2\times 2$ scattering matrix. The limited set of calculations that have been performed here so far have typically ignored this fact, for example in the study of $\pi\rho$ and $\pi\omega$ scattering in $J^{P}=1^{+}$ in Lang et al. (2014a), the ‘spinless’ Lüscher formalism was applied, effectively setting the $D$-wave amplitude to zero by fiat. In addition, not all possible meson-meson operator constructions corresponding to non-interacting levels in the energy region of interest were included, so it is not clear that the complete finite-volume spectrum was determined.

Similarly, to date, calculations of two-baryon systems Berkowitz et al. (2017); Beane et al. (2012b, 2013); Yamazaki et al. (2012, 2015) typically ignore the fact that different partial waves mix, either dynamically due to the tensor force, or due to the mismatch between the cubic boundary and rotational symmetry (in finite volume). The first attempt to constrain the contribution of the tensor force in the deuteron channel was carried out in Orginos et al. (2015), using ideas proposed in Briceno et al. (2013b, a), but the signal obtained, with heavier than physical light quarks, was consistent with zero.

The approach outlined in Sections V, in which the finite-volume spectrum is determined in a range of volumes and/or moving frames, providing a large number of energy levels to constrain the various channels which are coupled by spin, is likely to be the way forward to study such systems.

IX.2 Three-particle systems

A larger challenge is presented by the fact that most resonances in QCD with physical mass quarks are actually able to decay to three-body final states, and the formalism we have been primarily discussing in this document, which is expressed in Eq.(10), applies only to systems below three-body thresholds. This has motivated extensions of this formalism for energies where three-particle states can go on-shell Polejaeva and Rusetsky (2012); Briceno and Davoudi (2013b); Hansen and Sharpe (2015, 2014); Briceno et al. (2017b) 202020Interesting progress has also been made for 1+1 dimensional systems Guo (2016); Guo and Gasparian (2017)..

The current state of the art formalism is outlined in Briceno et al. (2017b), where a quantization condition was derived relating the finite-volume spectrum to scattering amplitudes for systems where two- and three-particle channels are coupled, for energies below the four-particle threshold. We do not present the result here, but it can be compactly written in form analogous to that of coupled two-particle systems, expressed in Eq.(20). In arriving at this result, two restrictions have been imposed, causing this to be not the most general equation for three-particle systems. First, all individual particles involved are assumed to be identical and carry no intrinsic spin, and as a result, there are only two channels that are kinematically open. This is a fairly mild assumption that is put in place only to simplify the kinematics. Second, the energies of the two-body subsystem (inside the three-particle system) must be below any pole in the two-particle $K$-matrix. In Hansen and Sharpe (2015, 2014) it was shown that poles in the $K$-matrix can lead to power-law finite-volume effects which have thus far not been accounted for. The removal of this restriction is not so simple, but it is not expected to be insurmountable.

Although this formalism has yet to be implemented in numerical studies, there have been two important formal checks performed for the simpler system considered in Hansen and Sharpe (2015, 2014), where two- and three-particle states decouple due to an additional symmetry in the system (e.g. $G$-parity in multi-pion systems). First, in Hansen and Sharpe (2016a, b) the authors demonstrated that this quantization condition reduces to previously determined perturbative expressions for the finite-volume corrections to the energy of the lowest lying state in a system composed of three weakly repulsive scalar bosons Beane et al. (2007); Tan (2008). Another important check was performed in Hansen and Sharpe (2017), where the authors showed that this quantization condition also reproduces the finite-volume corrections of Efimov bound states considered in Meissner et al. (2015).

IX.3 Elastic form factors of resonances

There are rigorously defined properties of resonances that can be computed in lattice QCD, which are not necessarily accessible in experiments. A good example are the form-factors of unstable particles, which may be extracted from the residue at the resonance pole of a $2\to 2$ transition amplitude featuring a current insertion of whatever construction is of interest. The required finite-volume formalism, whose construction is analogous to those presented in Section VIII is given in Briceno and Hansen (2016); Bernard et al. (2012). From the form-factors, one can proceed to rigorously determine radii or other moments for these states, and this could give a concrete tool for distinguishing compact and molecular states, for example. Furthermore, it was recently demonstrated by Ji Ji (2013) that quark distribution functions of QCD-stable states can be accessed non-perturbatively via lattice QCD. This idea, in conjunction with the formalism laid out in Briceno and Hansen (2016), could allow for the determination of the quark distribution functions of unstable states, allowing a view of the bound structure of resonance states.

X Outlook

The study of scattering processes and resonance properties within lattice QCD is entering an exciting period. In the last few years we have witnessed tremendous progress both in the development of the relevant formalism, and its application in explicit calculation. In this review we have presented some of the highlights of this program.

At the core of the formalism is the observation that infinite-volume scattering amplitudes control the discrete spectrum of eigenstates in a finite periodic box. We have discussed techniques for solving the inverse problem, where the discrete spectrum is extracted from a lattice QCD calculation and the scattering amplitudes are initially unknown. In order to get to the point where these techniques can be applied, it is necessary to reliably determine the excited state spectrum within a lattice QCD computation. We have illustrated the power of the variational method applied to matrices of two-point correlation functions computed using a large basis of operators, and shown the vital importance of including relevant multi-hadron operators in this basis.

Recent calculations have successfully determined elastic scattering amplitudes in several channels, and the extension into the coupled-channel sector has now been demonstrated in a number of cases.

Going beyond the simplest hadron-hadron scattering processes, we have seen progress in the development of formalism that handles reactions in which an electroweak current couples to a scattering system. As well as allowing rigorous first-principles QCD study of a range of phenomenologically interesting reactions, this formalism can also be used to compute processes which are experimentally inaccessible, such as the form-factors of hadron resonances.

Exploratory studies are being performed using unphysically large values of the quark masses, with reasoning beyond just the usual argument of decrease in cost: as the pion mass decreases, the threshold for three-hadron channels opening also decreases. Eq.(10) as presented is restricted in validity to energies below the lowest coupled three-hadron threshold. This unwanted constraint has motivated studies that are attempting to generalize Eq.(10) for systems where three or more particles can go on-shell. These works, which we have briefly reviewed, are advancing rapidly, but remain in the developmental stage. There is hope that the most general quantization condition involving two- and three-particle states will be derived in the next few years and its implementation will follow.

Acknowledgements.

We thank the authors of papers whose figures we reproduce in this article for their permission to do so, and S. Sharpe, M. Hansen and C. Thomas for their comments on a early draft of the manuscript. RAB and JJD acknowledge support from U.S. Department of Energy contract DE-AC05-06OR23177, under which Jefferson Science Associates, LLC, manages and operates Jefferson Lab. JJD acknowledges support from the U.S. Department of Energy Early Career award contract DE-SC0006765. RDY acknowledges support from Australian Research Council grants CE110001004, FT120100821 and DP14010306.

XI Numerical evaluation of $F$

The matrix $F(E,\mathbf{P};L)$, whose inverse appears in Eq.(10), is, in the case of single-channel scattering, an object whose matrix elements are given by Eq.(15). The appropriate generalization to multiple scattering channels, including the case where the scattering particles are not spinless is presented in Briceno (2014) – here we will focus on the numerical implementation of the simplest case of elastic scattering of spinless particles. In this case we can show that, for a frame momentum $\mathbf{P}=\frac{2\pi}{L}\mathbf{d}$,

[TABLE]

where $\gamma=E/E^{\star}$, b(\ell m;\bar{\ell}\overline{m}\,|\ell^{\prime}m^{\prime})=\sqrt{\frac{(2\ell+1)(2\bar{\ell}+1)}{2\ell^{\prime}+1}}\big{\langle}\ell m;\bar{\ell}\overline{m}\big{|}\ell^{\prime}m^{\prime}\big{\rangle}\big{\langle}\ell 0;\bar{\ell}0\big{|}\ell^{\prime}0\big{\rangle}, and where the Lüscher zeta functions are defined by

[TABLE]

In this expression the sum is performed over \mathcal{P}_{\mathbf{d}}=\big{\{}\mathbf{r}\in\mathbf{R}^{3}\big{|}\,\mathbf{r}=\hat{\gamma}^{-1}(\mathbf{n}-\alpha\mathbf{d})\big{\}} where $\mathbf{n}$ is a triplet of integers and $\alpha=\frac{1}{2}\left(1+\frac{m_{1}^{2}-m_{2}^{2}}{E^{\star 2}}\right)$ for scattering particle masses $m_{1},m_{2}$. The operation $\hat{\gamma}^{-1}$ is defined as ${\hat{\gamma}^{-1}\mathbf{x}=\frac{1}{\gamma}\mathbf{x}_{||}+\mathbf{x}_{\perp}}$ with $\mathbf{x}_{||}$, $\mathbf{x}_{\perp}$ being the components of $\mathbf{x}$ parallel to and perpendicular to the direction of $\mathbf{d}$. A longer discussion can be found in Rummukainen and Gottlieb (1995) and Leskovec and Prelovsek (2012).

The zeta function in Eq.(36) contains the factor $\left(|\mathbf{r}|^{2}-\left(\frac{q^{\star}L}{2\pi}\right)^{2}\right)^{-1}$, which is singular if ${q^{\star}=\frac{2\pi}{L}|\mathbf{r}|}$, which corresponds to a two-particle state taking an allowed non-interacting energy in the $L\times L\times L$ volume. We can easily illustrate this in the rest-frame ($\mathbf{d}=\mathbf{0}$) where $q^{\star}=\frac{2\pi}{L}|\mathbf{n}|$ and hence ${E^{\star}=\sqrt{m_{1}^{2}+\left(\frac{2\pi}{L}\right)^{2}|\mathbf{n}|^{2}}+\sqrt{m_{2}^{2}+\left(\frac{2\pi}{L}\right)^{2}|\mathbf{n}|^{2}}}$, which are precisely the non-interacting energy levels. Strategies for numerically evaluating the zeta function can be found in Rummukainen and Gottlieb (1995); Kim et al. (2005); Fu (2012); Leskovec and Prelovsek (2012); Yamazaki et al. (2004).

The matrix $F$ becomes block-diagonal if we consider the remaining symmetries of the cube (or the boosted cube), using the operation known as subduction into irreducible representations (irreps) of the relevant symmetry group. Subduction from $\ell,m$ into the $n^{\mathrm{th}}$ embedding of $\ell$ into the irrep $\Lambda$ (row $\rho$) is achieved by applying

[TABLE]

where $R$ is the rotation which carries a vector along $[001]$ into the direction of $\mathbf{d}$, and where we are summing over a quantity $\lambda$ that can be associated with helicity. The helicity-based subductions, $\mathcal{S}^{\Lambda\rho n}_{\ell\lambda}$ can be found in Thomas et al. (2012); Dudek et al. (2012). It follows that we can subduce $F$ as follows:

[TABLE]

where the orthogonality of different irreps and different rows within an irrep is expressed explictly on the left-hand-side. We can write the subduced $F$ in the form,

[TABLE]

and above threshold ($q^{\star}>0$), $f$ can be defined to be real, so that the imaginary part of $F$ is $\frac{1}{16\pi}\frac{2\,q^{\star}}{E^{\star}}\delta_{\ell,\ell^{\prime}}\delta_{n,n^{\prime}}$. Since elastic unitarity ensures that $\mathrm{Im}\,\frac{1}{\mathcal{M}_{\ell}}=-\frac{1}{16\pi}\frac{2\,q^{\star}}{E^{\star}}$, the Lüscher quantization condition, Eq.(10), which can be recast into the form \det\big{[}\mathcal{M}^{-1}+F\big{]}=0, is effectively then a real equation. Importantly if $\mathcal{M}$ is parameterized in a manner which does not respect unitarity, we can see that this equation will have an imaginary part which is volume-independent and which cannot be solved. These arguments generalize naturally to the coupled-channel case, and it is there that the condition that $\mathcal{M}$ satisfy unitarity is a practical constraint.

In Figure 24 we illustrate a typical behavior for $f^{\Lambda}$ in the case of the four one dimensional irreps, $A_{1},B_{1},B_{2},A_{2}$ of $\mathbf{d}=[110]$. Note that there are many divergences in these functions, and these divergences correspond to the possible non-interacting particle pairs in these irreps. Since the non-interacting spectrum is different in different irreps, the set of divergences differs depending upon $\Lambda$. Note also that there can be very tight regions in $(q^{\star}L/2\pi)^{2}$ between neighboring divergences – when we come to find solutions to Eq.(10) using numerical rootfinding techniques it will be important to ensure that we pay attention to this.

In Figure 25 we illustrate the behavior of some elements of a typical $f^{\Lambda}$ below threshold. We observe that the diagonal elements $f^{\Lambda}_{\ell n;\ell n}$, as we go further below threshold, we tend to value of $i$, and the off diagonal elements tend to [math]. This ensures that $F^{\Lambda}\to 0$ as we go some way below threshold. It follows that bound-states lying far below threshold appear in finite-volume at an energy that is very close to their infinite-volume energy. Another consequence is that in the case of multi-channel scattering, a kinematically closed channel does not have a significant impact on the quantization condition at energies far below its threshold, so it is not necessary to include in Eq.(10) the continuation of distant closed channels. Note that there is a narrow region below threshold where the finite-volume functions are not zero, and this indicates that the effect of a closed channel can ‘leak’ down a little way below its threshold – the effect is exponentially suppressed with increasing volume.