The $b_1$ resonance in coupled $\pi\omega$, $\pi\phi$ scattering from lattice QCD

We present the first lattice QCD calculation of coupled $\pi\omega$ and $\pi\phi$ scattering, incorporating coupled $S$ and $D$-wave $\pi\omega$ in $J^P=1^+$. Finite-volume spectra in three volumes are determined via a variational analysis of matrices of two-point correlation functions, computed using large bases of operators resembling single-meson, two-meson and three-meson structures, with the light-quark mass corresponding to a pion mass of $m_\pi \approx 391$ MeV. Utilizing the relationship between the discrete spectrum of finite-volume energies and infinite-volume scattering amplitudes, we find a narrow axial-vector resonance ($J^{PC}=1^{+-}$), the analogue of the $b_1$ meson, with mass $m_{R}\approx1380$ MeV and width $\Gamma_{R}\approx 91$ MeV. The resonance is found to couple dominantly to $S$-wave $\pi\omega$, with a much-suppressed coupling to $D$-wave $\pi\omega$, and a negligible coupling to $\pi\phi$ consistent with the `OZI rule'. No resonant behavior is observed in $\pi\phi$, indicating the absence of a putative low-mass $Z_s$ analogue of the $Z_c$ claimed in $\pi J/\psi$. In order to minimally present the contents of a unitary three-channel scattering matrix, we introduce an $n$-channel generalization of the traditional two-channel Stapp parameterization.


I. INTRODUCTION
Contemporary studies of hadron spectroscopy seek to relate the spectrum of hadron resonances, including their decay properties, to the fundamental theory of quarks and gluons, quantum chromodynamics. The most successful theoretical technique to achieve this has proven to be lattice QCD which considers the theory on a discretized space-time grid of finite size, allowing numerical calculation of correlation functions through averaging over Monte-Carlo generated field configurations. The discrete spectrum in a finite volume corresponding to a particular choice of quantum numbers can be extracted from a matrix of correlation functions, constructed using a basis of operators which resemble the hadronic system being studied. The fact that lattice QCD studies the theory in a finite volume can be turned to our advantage -an approach introduced by Lüscher relates the discrete spectrum in a finite volume to hadron-hadron scattering amplitudes. Initially this was only for elastic scattering of spinless particles with the system overall at rest with respect to the lattice [1][2][3][4], but subsequent extensions generalize the formalism to describe coupled-channels, particles with intrinsic spin, and moving frames [5][6][7][8][9][10][11][12].
This approach has been applied to a number of cases in which several coupled pseudoscalar-pseudoscalar channels are present, for example πη, KK in which the scalar a 0 appears as a resonance [13], or ππ, KK, ηη where scalar f 0 and tensor f 2 resonances appear [14]. Pseudoscalar-pseudoscalar scattering with relative orbital angular momentum defines the 'natural parity' sequence, J P = 0 + , 1 − , 2 + , . . ., where J is the angular momentum and P is the parity. To observe resonances with two-body decays in the 'unnatural parity' sequence, J P = 0 − , 1 + , 2 − , . . ., we must consider the scattering of mesons with non-zero spin. An experimentally-observed example [15] is the b 1 (1235) resonance which is dominantly seen through its decay to the πω final state, where the ω is the lightest isoscalar vector meson which has a very small decay width to three pions.
Once we move into the pseudoscalar-vector scattering sector, there can often be more than one partial-wave construction having a particular J P . For example, in the 1 + case relevant for the b 1 , we can have the π and ω in a relative S-wave or a relative D-wave -indeed, by studying the angular distribution in the decay of the b 1 , experiments have estimated the amplitudes of these two partial-waves [16].
The finite-volume formalism to handle pseudoscalarvector scattering is in place [5], and has been tested previously in a channel which did not feature a resonance, namely πρ scattering in isospin-2 with quark masses sufficiently heavy such that the ρ resonance becomes a boundstate, kinematically stable against decay to two pions [17]. That first calculation determined the S-and D-wave J P = 1 + amplitudes and their dynamical mixing, finding relatively weak effects as expected in this exotic isospin channel.
In this paper we will report on a study of the J P = 1 + I G = 1 + channel, where I is the isospin and G is Gparity, in which we expect to see a b 1 resonance decaying to πω. We make use of N f = 2 + 1 lattice configurations generated with a light-quark mass such that the pion has a mass around 391 MeV. With this light-quark mass, the ω meson is found to have a mass around 881 MeV [18,19], and hence is stable against decay to three pions.
To study the b 1 we have computed matrices of correlation functions in three lattice volumes, in several moving frames (i.e. where systems have overall non-zero momentum with respect to the lattice). To robustly determine the finite-volume spectrum, a wide range of operators resembling both single-hadron and multi-hadron structures were included in the basis. These correlation functions provide information which constrains the energy dependence of the I G = 1 + J P = 1 + scattering matrix whose channels are πω in S and D-wave, and, in addition, πφ which is kinematically open in the considered energy region 1 .
A previous lattice QCD study [20] of the b 1 limited itself to the rest-frame in one rather small volume. By considering only two degenerate flavors of light quarks and no strange quarks, any physics associated with the πφ channel was disallowed. A very small operator basis was used, such that only one usable energy level was obtained and this had a statistical uncertainty at the percent level. Enforcing elastic S-wave scattering only, ignoring any effect from the D-wave, and fixing the decay coupling of an assumed b 1 resonance at a value equal to that extracted from experimental measurements, a crude estimate of the b 1 mass was made in the case that the pion mass is 266 MeV. An earlier study [21] used a different approach in which the light-quark mass was tuned such that the b 1 decay to πω is exactly at kinematic threshold. From the time-dependence of a single correlation function, an estimate of the decay coupling was inferred.
In this calculation, we determine a large number of finite-volume energy levels in multiple volumes and moving frames. We use up to 36 of these levels, each typically having statistical uncertainty at the tenth of a percent level, to constrain the coupled-channel scattering matrix.
As well as the πω and πφ channels, we pay attention to the fact that three-body channels, ππη and πKK, which have relatively low thresholds even for m π ≈ 391 MeV, can in principal play a role. Experimentally, three-body decays of resonances are found to be dominated by twobody isobar resonances. For example, in a ππη final state at relatively small total energy, the Dalitz plot will be expected to have the bulk of the events in narrow horizontal and vertical bands around m ππ ∼ m ρ and m πη ∼ m a0 . 2 We will explore the role of these three-body channels by including operators in our bases whose construction resembles a meson coupled to a two-body resonance, in a way which respects the symmetries of the finite cubic lattice. No finite-volume formalism capable of rigorously 1 the φ is stable against decay to KK and πππ at the light-quark mass considered 2 There will also be a diagonal 'reflection' of the a 0 band. incorporating three-body scattering channels is yet sufficiently mature to be applied in the current case, but there has been significant recent developments [22][23][24][25][26][27][28][29]. Our explorations will yield evidence that suggests that the three-body channels have a negligible effect in this particular case of a low-lying b 1 resonance.
To convert the finite-volume spectra calculated in lattice QCD into scattering amplitudes, we consider parameterizations of the energy dependence of the scattering t-matrix and the parameters are found which best describe the finite-volume spectra. This approach allows us to explore the resonance content of each J P in a rigorous way by searching for the presence of pole singularities in t(s) at complex values of s = E 2 . Poles lying relatively close to the real energy axis typically have the real and imaginary parts of their pole position interpreted in terms of the mass and width of the resonance, and from the residue of t(s) at the pole we can determine the relative couplings of the resonance to its decay channels.
A relatively light b 1 resonance is expected based upon an earlier set of calculations, performed on the same lattice configurations used in this paper, in which the operator basis was restricted to a set of fermion bilinears [18,19,30]. The resulting spectrum, which we expect to be incomplete owing to the lack of multi-meson operators, nevertheless featured a J P C = 1 +− state near 1400 MeV, which had strong overlap with, in particular, those operators which resemble the qq spin-singlet, P -wave structure expected for the b 1 in the quark model. Such a calculation can do no more than indicate to us the likely presence of a narrow resonance -in the current calculation we will rigorously determine its presence and properties.
It has been suggested [31] that the πφ channel, coupled to πω, may feature a Z s resonance analogous to the Z c enhancement that has been claimed in the πJ/ψ final state [32,33]. We will find no evidence of a Z s resonance in this work.
The remainder of this paper is structured as follows. In Section II we briefly review the calculation of finite-volume spectra from correlation functions and describe our single-, two-and three-meson operator constructions. The lattice setup used and relevant hadron masses and thresholds are presented in Section III, and in Section IV we discuss the partial waves which are present and our choice of operator bases. The finite-volume spectra are presented and commented on in Section V. In Section VI we discuss the techniques used to relate these spectra to scattering amplitudes and apply them to determine πω 3 S 1 , πω 3 D 1 and πφ 3 S 1 amplitudes, and in Section VII we examine the pole singularities of these amplitudes. Systematic tests of our analysis are given in Section VIII where we examine the effects of additional partial-waves, including those that mix due to the reduced symmetry of the finite-volume, and additional channels that resemble ππη and πKK. An interpretation of the results is provided in Section IX and we conclude with a summary in Section X.

II. SPECTRAL DETERMINATION AND OPERATOR CONSTRUCTION
Working in a cubic volume of size L × L × L with spatially periodic boundary condition discretises momenta, restricting to values P = (2π/L)(n x , n y , n z ) where n i ∈ Z. For particles at rest with respect to the lattice, the infinitevolume O(3) spatial symmetry is broken to that of the double cover of the octahedral group with parity, O D h , and total angular momentum J and parity P labelling the irreducible representations (irreps) of O(3) are replaced by Λ P , the irreps of O D h . In this work we only encounter integer spin and therefore irreps of the single cover O h . For particles "in-flight", i.e. moving with respect to the lattice, parity is no longer a good quantum number and the irreps Λ are those of the little group of symmetries, LG( P ), as discussed in Ref. [34]. We write lattice irreps P Λ with shorthand P = [n x n y n z ], omitting units of (2π/L) for brevity.
In order to robustly determine the discrete finitevolume energy eigenstates in each irrep, P Λ, we first compute a large matrix of two-point correlation functions, C(t) ij = 0|O i (t + t src )O † j (t src )|0 , by employing a diverse basis of operators O i . These operators are constructed with the desired flavour structure and subduced into the irrep P Λ [30,35]. A variationally optimal determination of the spectrum [3,36] follows from solving the generalized eigenvalue problem for each irrep, The energy levels E n are determined by fitting principal correlators λ n (t) to the form, where the second term soaks up any residual excited state contamination. The eigenvector v n can be used to construct a variationally optimised operator, Ω † n = i v n i O † i , efficient at interpolating the n th eigenstate in the spectrum. We refer the reader to Refs. [30,37] for further details of our implementation and techniques for selecting a reasonable value of t 0 .
Previous calculations [17,38,39] have demonstrated the importance of having sufficiently 'complete' operator bases in order to reliably determine the complete spectra in a given energy region. The region we study includes the opening of several multi-hadron thresholds: πω, πφ, ππη and πKK, and we find that this necessitates the inclusion of two-meson-like and three-meson-like operators in our basis, as well as single-meson operators of fermion-bilinear form which we expect to have good overlap with any bound state or relatively-narrow resonance present. Fourmeson thresholds lie beyond the energy region we consider, and previous calculations suggest that local tetraquarklike operators have little effect on the spectra [39,40], so neither of these types of operators are included in the basis. The construction of interpolating operators resembling single-meson, two-meson and three-meson structures is discussed in the subsections which follow.

A. Single-meson operators
The construction of 'single-meson-like' operators follows the procedure detailed in Refs. [30,35]. To summarise, fermion bilinearsψΓ ← → D ... ← → D ψ are constructed with definite J P and z-component of angular momentum M by appropriately coupling products of gauge-covariant derivatives ← → D and Dirac γ-matrices Γ. These are then projected onto definite momentum P and appropriate linear combinations yield continuum single-meson operators O †JM M ( P , t) of definite flavor, labelled by M. Schematically, where for P = 0 we use helicity operators, labelled by helicity λ rather than M , as discussed in Ref. [35]. Singlemeson operators, transforming irreducibly under the symmetry of the lattice grid and boundary, O †Λµ M ( P ), are obtained by subducing, where S JM Λµ are subduction coefficients tabulated in Refs. [30,35].
A large basis of operators can be constructed by combining γ-matrices with various numbers of derivativeshere we use up to three derivatives for operators with zero momentum and up to two otherwise. Single-meson operators are written asψΓψ for the remainder of this article.
Optimised operators for the stable ω (Ω † ω ) and φ (Ω † φ ) in each relevant irrep follow from variational analysis of a matrix of correlation functions constructed using a basis of quark bilinears with both hidden-light (ūΓu +dΓd) and hidden-strange (sΓs) flavor structure. The required 'annihilation' diagrams are computed but, as shown in Figures 4 and 5 of Ref. [19], they prove to be small in the vector channel in line with the experimentally-motivated 'OZI rule'. In each irrep, the ω appears as the ground state, dominated by overlap withūΓu +dΓd, and the φ as the first excited state, dominated bysΓs.
The same flavor basis is used to determine the optimum η operator (Ω † η ) in each irrep, but here significant mixing between light and strange is observed through the annihilation diagrams (see Figures 2 and 3 in Ref. [19]), indicating, as is well known, that the OZI rule does not apply in the pseudoscalar channel.
The need to account for 'in-hadron annihilation' when considering isoscalar mesons will reappear when the optimized operators are used in two-meson and three-meson constructions as discussed below.
We construct two-meson operators with definite flavor and momentum in irrep Λ (row µ) by taking appropriate linear combinations of the products of optimised singlemeson operators Ω † M , each independently constructed to transform irreducibly in some lattice irrep. Schematically, where the sum is over the rows µ i of the irreps Λ i and the sets of momenta { p i } * , containing all momenta related to p i by an allowed lattice rotation with the total momentum p 12 = p 1 + p 2 fixed -see Eq. 3.3 of Ref. [17]. For | p i | 2 < 9(2π/L) 2 , the set { p i } * is equivalently labelled by the magnitude of the momentum | p i |. The sum is weighted by lattice Clebsch-Gordon coefficients, For energies below three-meson thresholds, previous calculations suggest that a sufficient set of operators for a reliable calculation of the spectra consists of singlemeson and two-meson operators. Two-meson operators O †Λµ M1M2 ( p 12 ) are efficient at interpolating the finite-volume energy levels near to the associated non-interacting energies, and truncating the two-meson operator bases when the corresponding non-interacting energies are beyond the energy region of interest has been demonstrated to be sufficient for a robust determination of the spectra [13,14,17,38,39,[41][42][43][44][45][46]. Two-meson operators are written in all tables and figures for the remainder of this work. The fact that a vector meson in flight is subduced into multiple irreps means that there can be multiple MM constructions for a single non-interacting energy. For example, π 001 ω 001 subduced into the [000] T + 1 irrep (which contains J P = 1 + ) can be constructed independently from π(A 2 )⊗ω(A 1 ) or from π(A 2 )⊗ω(E 2 ). Cases such as these where the multiplicity of operators is greater than one are discussed in detail in Ref. [17], and we indicate them with a notation {n}.
Correlation functions with MM operators at the source and/or sink feature Wick contractions in which quarks annihilate either within an isoscalar meson or between two mesons. Considering a basis with overall I = 1 as relevant here, with M =ūΓd and MM = {πω, πφ}, we need to evaluate diagrams whose structure is similar to those shown in Figure 1 of [46].

C. Three-meson operators
Three-meson operators 3 can be constructed by iteratively applying the two-meson operator construction outlined above. Schematically, where O †Λµ M1M2 is a two-meson operator constructed from a product of optimised single-meson operators as in Section II B. Note that it does not matter with which optimised single-mesons we formed the intermediate , as the tensor product is associative. An argument for determining a sufficient set of three-meson operators, analogous to that presented previously, would suggest calculating the corresponding non-interacting energies and enforcing a similar truncation on the basis. While this approach has the advantage of being straightforward, it pays no attention to the fact that we expect certain two-meson pairs to feature resonating behavior, the finitevolume analogue of the Dalitz-plot enhancements mentioned in the introduction. Consider the example of πππ in isospin-2. Following the construction above, we would be attempting to describe energy eigenstates of the ππ isospin-1 subsystem using O †Λµ M1M2 constructed using only 'ππ'-like operators. To reliably determine the isovector ππ spectra, i.e. the ρ spectra, an operator basis including bothψΓψ and ππlike operators is needed as shown in Figure 1 of Ref. [42]. An alternative approach, based upon this observation and used in Ref. [39], utilizes an optimised two-meson operator which will be a linear combination ofψΓψ and ππ-like operators. We denote such an optimised operator Ω † R , where R indicates the meson with the corresponding quantum numbers, i.e. Ω † ρ for the example above. 4 In general, multiple optimised operators may be relevant -Ω † R n denotes the optimal interpolating operator for the n th excited state in the relevant meson-meson subsystem. 3 and operators with a structure resembling more than three mesons 4 Lattice irreps contain more than one spin but for convenience we choose the label R corresponding to the lightest such meson, e.g. in [000]T − 1 we choose ρ.
Combining these operators with an optimized singlemeson operator yields an alternative set of three-meson operators, given schematically by By design, we anticipate that these three-meson operators will efficiently interpolate finite-volume levels in the region of an energy value where E Λ12 R n

12
( p 12 ) are finite-volume energies calculated in the two-meson subsystem in irrep [ p 12 ]Λ 12 , i.e. they will efficiently capture interaction in the two-meson subsystem assuming weak residual interaction with the third meson. Calculating E (2+1) n.i. energies, for all possible combinations of two-meson subsystems that together with the third meson give the desired quantum numbers, and truncating at a desired energy, provides a procedure for selecting which of these three-meson operators to include in the basis.
To illustrate the construction presented above, consider the example of a three-meson operator resembling ππη in the irrep [000] T + 1 with I G = 1 + . We begin with the construction shown in Eq. 4. For p 1 = p 2 = p 3 = 0, there is only one possible irrep, so no non-interacting ππη level, or corresponding operator, appears in [000] T + 1 at threshold. If the pions are both given one unit of momentum, p 1 = p 2 = [001] and p 3 = 0 (recalling that directions of momenta p i are summed over as detailed in Section II B), the product n.i. = 2 m 2 π + 2π L 2 + m η . Now we consider bound-states and resonances in the ππ and πη two-meson subsystems and construct operators according to Eq. 5. Unlike in the previous construction, the order in which we combine the single-meson operators does matter as the intermediate Ω † R depends on the flavor structure of the two-meson subsystem. As before, for p 1 = p 2 = p 3 = 0 there is no [000] T + 1 , while for p 1 = p 2 = [001] and p 3 = 0, there are two possible distinct two-meson subsystems.
First, for the ππ subsystem, there are three possible flavor combinations, I G = 0 + , 1 + , 2 + , and three possible irreps with momentum p 12 = 0, namely [000] A + 1 , [000] T − 1 and [000] E + . When combined with the η, only the ππ subsystem with I G = 1 + transforming in [000] T − 1 gives the desired overall flavor and irrep. This ππ subsystem contains quantum numbers corresponding to the ρ and the construction is, schematically, Calculating the E (2+1) n.i. energies amounts to determining the ρ-like energy eigenstates in [000] T − 1 with I G = 1 + and adding these to the η energy according to Eq. 6, where we recall that ρ n denotes the n th energy eigenstate within the irrep. In many cases, including here, only the lowest energy two-meson state (n = 0) yields an operator below the energy cut-off.
The second possible construction considers the πη subsystem where there is only one flavor combination, I G = 1 − , and one possible irrep, [001] A 1 . These quantum numbers correspond to the a 0 meson. Schematically, and, as before, we determine the E (2+1) n.i. energies by calculating the a 0 -like energy eigenstates in [001]A 1 with I G = 1 − and add these to the π energy according to Eq. 6, For each E (2+1) n.i. below some energy cut-off we can construct operators of the form O † a0π via Eq. 5. The E A1 a0 n ([001]) energies are an example of a case where it may be prudent to consider multiple states (n ≥ 0) in the two-body sector. Figure 4 of Ref. [13] shows the [001]A 1 spectra corresponding to the E A1 a0 n ([001]) energies -there are many nearby low-lying energy levels on each volume. Following the construction given in Eq. 5 leads to multiple operators of the form O † a0π corresponding to similar E (2+1) n.i. .
The use of RM operators to efficiently interpolate finitevolume states above three-meson thresholds requires the calculation of a large number of diagrams. As an example, consider the case of an a 0 π operator at the sink, where the optimized a 0 operators are linear superpositions ofūΓd, πη and KK constructions (see Table VII). This leads to the diagram components shown in Figure 1, which need to be connected to the quark lines from the π and the source operator to form complete Wick contractions. It follows that even in the simple case of b 1 − a 0 π correlators we would have diagrams with the structures shown in Figure 2.

III. LATTICE SETUP
Correlation functions were computed on anisotropic lattices of spatial volumes (L/a s ) 3 = 16 3 , 20 3 and 24 3 each having temporal extent T /a t = 128, where the temporal lattice spacing, a t , is finer than the spatial lattice spacing, a s ∼ 0.12 fm, with an anisotropy ξ = a s /a t ∼ 3.5. Gauge fields were generated from a tree-level Symanzik-improved gauge action and a Clover fermion action with N f = 2 + 1 flavors of dynamical quarks where the strange quark is tuned to approximately its physical mass and the degenerate light quarks are such that m π ∼ 391 MeV [47,48]. We utilize the distillation framework [49] to compute correlation functions as successfully demonstrated in many previous works. All relevant Wick contractions were calculated within this framework without requiring additional propagator inversions beyond the basic set of t src − t and t − t 'perambulators' for light and strange quarks which were computed for use in previously reported calculations. The very large number of diagrams incurs only a combinatoric cost associated with the contraction of the perambulators with the operator constructions.
Correlation functions were computed using the number of distillation vectors, gauge configurations and timesources shown in Table I. Typically, we calculated all the elements of the matrix of correlation functions, including the transposes, C ij and C ji , which are related by hermiticity. In a few cases where there are a particularly large number of diagrams contributing, we made use of hermiticity to infer C ji from the computed C ij .
Masses of relevant stable hadrons are shown in Table II, where π, K, η ( ) and σ masses are taken from Refs. [41], [46], [13] and [14] respectively. Using energy levels on three lattice volumes, we determine the masses and anisotropies of the ω and φ mesons from fits to the dependence of the energy of a stable hadron of momentum, p = (2π/L) n, up to discretisation effects, as shown in Figure 3. We observe the same characteristic splitting between the |λ| = 0, 1 components as was found for the stable ρ meson in Ref. [17] at larger quark mass, and we attribute this splitting to discretisation effects given that the finitevolume effects here are small. The values of a t m ω , a t m φ and ξ we use are obtained by taking the largest variations within one standard deviation of the means across the different helicities. This yields the masses given in Table II and an anisotropy ξ = 3.443 (48) which is consistent with the anisotropies previously determined for π, K and η [13,41,46].

IV. PARTIAL WAVES AND OPERATOR BASES
In this study we are principally interested in irreps that contain J P = 1 + . For irreps at rest, J P = 1 + subduces only into T + 1 . However, for in-flight irreps, different helicity components of J P = 1 + are subduced across multiple irreps as shown in Table II   λ = 0 and ±1 subduce into A 2 and E 2 respectively for overall momentum P = [001]. Furthermore, at non-zero momentum parity is no longer a good quantum number and so many irreps contain both J + and J − , e.g. 1 + and 1 − .
We will restrict our attention to P A 2 in-flight irrepsthese contain subductions of the λ = 0 part of J P = 1 + but, because reflection parityη = P (−1) J is a good quantum number for λ = 0, they do not contain J P = 1 − . In contrast, [001]E 2 contains J P = 1 − as well as J P = 1 +the latter gives comparatively lower-lying J P = 1 − levels, TABLE III. Partial-wave J P ( 3 J ) content for pseudoscalarvector scattering in irreps P Λ containing J P = 1 + , transcribed from Ref. [17]. A subscript [N ] indicates that this J P has N embeddings in that irrep.
as seen in Ref. [42], and so will lead to a dense spectrum of mixed J P = 1 + and 1 − energy eigenstates. Considering only P A 2 allows us to avoid the complication of disentangling the J P = 1 + and 1 − scattering amplitudes.
The partial-wave content of a pseudoscalar-vector system for irreps [000] T + 1 and P A 2 with | P | 2 ≤ 4(2π/L) 2 is given in Table III. There we make use of the 2S+1 J notation for meson-meson scattering, where 2S + 1 = 3 reflects the unique spin-coupling in pseudoscalar-vector scattering, and is the relative orbital angular momentum. The use of the − S basis, over say the helicity basis, is for convenience; in particular, the threshold behavior of a partial-wave of definite is known. Table III includes cases where two 2S+1 J constructions appear with the same J P -in these cases the scattering matrix is 2 × 2 in the case of single meson-meson channel, e.g. for J P = 1 + scattering of πω, the t-matrix is where the symmetric nature of the matrix follows from time-reversal invariance.
K * K operators are constructed with definite G-parity analogous to the KK operators in Ref. [38]. For the ππσ-threshold, three-meson operators resembling ρσ and a 1 π were considered for inclusion. These appear in a relative P -wave in the [000]T + 1 and P A 2 irreps at values of E (2+1) n.i. that lie far above ππππ-threshold. Similarly, relevant ππσ non-interacting energies, E n.i. , are far beyond ππππ-threshold. Although the construction of operators resembling four-mesons could be done analogously to the three-meson operator construction described above, we do not include these in our basis and choose to restrict to energies below the ππππ-threshold.
The operator basis used for the [000] T + 1 irrep on each lattice volume is presented in Table IV. Included are all two-meson and three-meson operators corresponding to E (2) n.i. and E (2+1) n.i. below ππππ-threshold. 6 The operator lists for P A 2 irreps with P = 0 are presented in Appendix B -we include, as well as all low-lying two-meson operators, also the lowest three-meson (RM) operator in each irrep, with the intention of robustly determining the spectra up to the lowest E n.i. energy. As well as providing many more energy levels with which to constrain the scattering matrix, moving frames are required to determine the sign of the off-diagonal element, t πω 3 S 1 |πω 3 D 1 , as previously explored for πρ scattering in Ref. [17].
In order to estimate the strength of partial-waves with J ≥ 2 that appear alongside our desired J P = 1 + in [000] T + 1 and P A 2 , on the largest volume we also computed spectra in irreps [ partial-wave. The operator bases used for these irreps are presented in Appendix B.
In summary, because we are considering the G-parity positive isovector sector, the neutral channels have chargeconjugation C = −. The contributing J P C includes our target 1 +− where we expect a low-lying b 1 resonance, which in the quark model would be a qq spin-singlet in a P -wave. 2 −− and 3 −− are expected to resonate at a somewhat higher energy, corresponding to ρ 2 , ρ 3 resonances which would be spin-triplet D-waves in the quark model. Still higher we might have a 3 +− resonance, b 3 , as a spin-singlet F -wave qq. 0 −− and 2 +− are exoticthey do not appear in the qq quark model and previous lattice calculations [30] suggest that they may resonate in the form of hybrid mesons at much higher energy. Because they do not resonate, and feature at least a Pwave threshold suppression, it follows that we expect all partial waves except J P = 1 + to be small at low energies, and indeed we will find this to be the case below.
V. FINITE-VOLUME SPECTRA The spectra determined from variational analysis of [000] T + 1 correlation matrices on three volumes using the operator bases in Table IV are presented in Figure 4. For the largest lattice volume (L/a s = 24), the principal correlators and operator-state overlaps, Z n i = n|O † i (0)|0 , are also provided for illustration. The typical magnitude of statistical uncertainty on the energy levels, even relatively high in the spectrum, is at the level of a few tenths of a percent. It should be clear from the operatorstate overlaps that our operator basis is rather efficiently 'latching on' to the finite-volume eigenstates. In some cases an eigenstate has a dominant overlap with only one operator, suggesting that the state closely resembles that particular operator structure.
Consider first the number of energy levels expected below a t E cm ≈ 0.27 on each volume. In the absence of residual meson-meson interactions we would expect four on each lattice volume: one at each of the two E  Table IV. Solid curves are two-meson non-interacting energies, atE (2) n.i. , short dashed horizontal lines are atE (2+1) n.i.
, and long dashed horizontal lines show the two-, three-, and four-meson thresholds. Multiplicities (if greater than one) are shown as {n}. For each energy level on the largest volume, we show the principal correlators, plotted as λn(t, t0) e En(t−t 0 ) for t0 = 10 at, and histograms showing the operator-state overlap factors, Z n i = n|O † i (0)|0 , for the MM = πω (dark blue), πφ (green) and RM = ρη (blue-green), K * K (purple) operators along with a sample set of single-meson operators subduced from J P = 1 + (red) and J P = 3 + (orange). The overlaps are normalized such that the largest value for any given operator across all energy levels is equal to one. Right: The spectrum extracted when ρη and K * K operators are excluded from the basis (black) compared with the complete spectrum (gray).
as short dotted horizontal lines. Counting the number of energy levels actually extracted, we find five, with an 'additional' level appearing near πφ threshold. This may suggest the existence of a narrow resonance, as seen in calculations of the ρ resonance [38,42], with a mass close to πφ threshold. 7 On the largest volume, the effect of there being two ways to construct π 001 ω 001 can be seen: two energy levels are found, one very close to the noninteracting energy and one somewhat higher in energy.
In Figure 4 we also present an investigation of the importance of including RM operators in the basis. The rightmost panel shows the spectrum extracted when ρη and K * K operators are excluded, compared to the spectrum extracted with the full basis -with the smaller basis we see that typically the levels close to the ρη and K * K 'noninteracting' energies are no longer found. The spectrum at lower energies shows only modest discrepancies, except on the smallest lattice volume (L/a s = 16) where we might indeed expect the finite-volume effects associated with ρη and K * K to be largest. Finding 'incorrect' spectra due to 'incomplete' operator bases has been demonstrated in previous works. One example can be seen in Figure 1 of Ref. [42] where including bothψΓψ and ππ operators is shown to be essential in order to robustly determine the ρ spectrum. Figure 4 demonstrates an analogue of this for the case of three-meson operators.
Some qualitative observations about the spectrum can be gleaned from the operator-state overlap factors shown in Figure 4. The energy level just below πω threshold on all volumes has significant overlap onto both π 000 ω 000 and ψΓψ operators, as one might expect if a qq-like resonance lies nearby. For the two levels in close proximity to πφ threshold, one appears dominated byψΓψ operators with some overlap onto πω, ρη and K * K operators, while the other is completely dominated by π 000 φ 000 . Furthermore, we observe that all other levels have very small overlaps with the π 000 φ 000 operator, reflecting the fact that the matrix of correlation functions is approximately block diagonal with respect to π 000 φ 000 . This suggests that πφ is essentially 'decoupled', as might be expected from the 'OZI rule' which postulates that qq pairs in isoscalar mesons prefer not to annihilate. The states close to the ρη and K * K 'non-interacting' energies are observed to have large overlap with ρη and K * K operators respectively. The highest two states shown, near to the π 001 ω 001 twofold degenerate non-interacting energy, differ somewhat in their overlaps. The level shifted up has overlap with both the π 001 ω 001 andψΓψ operators, while the other, which lies on the non-interacting energy, has significant overlap only with the π 001 ω 001 operators.
In Figure 5, we present the cm-frame finite-volume spectrum for irreps [000] T + 1 and P A 2 on the three volumes, with only those levels found below the lowest E n.i. shown. 8 Points in grey are levels that prove to be sensitive to the presence of ρη, K * K and a 0 π operators in the basis, or which are very close to the energy cut-off, and these levels are excluded from the main scattering analysis in Section VI. Although we take a conservative approach and exclude these levels, we will find in Section VI that they are mainly well described by the scattering amplitudes, and we re-examine these levels in Section VIII.
For irreps P A 2 with P = 0, the density of energy levels is much higher than in irreps at rest -more momentum combinations for two-and three-mesons with associated E (2) n.i. , E (2+1) n.i. and E n.i. lying below the ππππ-threshold are possible. This can make identifying an 'additional' level more challenging in these irreps. However, in the [111] A 2 irrep we can clearly see an additional energy level on each volume relative to the number expected from counting the non-interacting two-meson energies. We also observe an 'avoided level crossing' where the π 000 ω 111 non-interacting energy crosses a t E cm ∼ 0.25, another hint that we may have a narrow resonance in this energy region.
In Figure 6, we present finite-volume spectra on the largest lattice volume for irreps [ where the lowest contributing spin is J = 2. We observe very little deviation of the extracted energy levels from non-interacting πω energies, suggesting that the πω scattering amplitudes in J ≥ 2 partial-waves are very small in this energy region. We also find levels in [001]B 1 and [001]B 2 consistent with non-interacting ππ energies and with dominant overlaps onto ππ operators. This is in line with the results of Ref. [42] where the ππ{ 1 F 3 } amplitude (J P = 3 − ) was found to be consistent with zero in this energy region. We also find a level in [001]B 1 consistent with the non-interacting KK energy and with 8 Errorbars on the energy levels include estimates of systematic uncertainly coming from varying t 0 and fitting time ranges, and reasonable variations of the operator basis. Also included is the effect of the uncertainty on the anisotropy which appears when we boost back from the 'lab' energy to the cm frame.
dominant overlap onto KK operators, suggesting that the opening of the KK threshold does not enhance the scattering in J P = 3 − .

VI. SCATTERING ANALYSIS
Finite-volume energy levels and infinite-volume scattering amplitudes are related through a quantisation condition derived by Lüscher [1,2,4] and extended by many others [5][6][7][8][9][10][11][12]50] to accommodate the most general case of two particle scattering. The quantisation condition, subduced into lattice irrep P Λ, can be expressed as the determinant of a matrix in the space of intrinsic spin S = S 1 ⊕ S 2 , orbital angular momenta , total angular momenta J, the embedding number n of a particular partial-wave in lattice irrep P Λ and hadron-hadron channel a. Written compactly, following the notation of Ref. [17], where the determinant over intrinsic spin is trivial in the current case as S takes only the value 1 for vectorpseudoscalar scattering. Here t(E cm ) is the scattering t-matrix, 9 diagonal in J with components t Ja, Jb . The diagonal matrix of phase-space factors ρ(E cm ) has components where k (a) is the cm-frame momentum for hadron-hadron channel a, .
Both t and ρ, being infinite-volume quantities, are diagonal in embedding number n and we have dropped this index for brevity. Lastly, M(E cm , L) is a matrix of known functions, diagonal in hadron-hadron channel, describing the kinematics of the system in a finite cubic volume, with components M n Ja,n J b . The subduced quantisation condition in Equation 11 reflects the little-group symmetry. The finite-volume spectrum calculated in irrep P Λ depends upon the various partial-wave amplitudes present in that irrep (see, for example, Tables III and V). In this way computing spectra in multiple irreps offers additional constraints on scattering. Further details can be found in Appendix C of Ref. [17]. Equation 11 is limited to describing two-body scattering -developments in the pursuit of a corresponding threebody formalism [22][23][24][25][26][27][28][29][51][52][53][54][55] n.i. . Black points are used in the scattering analysis in Section VI while gray points are excluded from the main analysis as discussed in the text. Solid curves are two-meson non-interacting energies, atE progress towards a general quantisation condition, but they are not yet mature at the level where we could apply them in the case considered in this paper. We therefore mainly restrict our attention to the two-body channels, πω and πφ, and in Section VIII we estimate the systematic effects of neglecting the three-body channels ππη and πKK in the energy region considered, finding them to be small.
In the case of elastic scattering, where a single mesonmeson channel appears in a single partial-wave, the t-matrix can be expressed in terms of a single real energy- In this case we can invert Eq. 11 to obtain a one-to-one relation between E cm and δ(E cm ). Given a set of discrete energy levels below the inelastic threshold, we can hence obtain a set of phase-shift points. In the case considered in this paper, there is never rigorously elastic scatteringas soon as the πω threshold opens, in J P = 1 + there are always two coupled partial-waves, 3 S 1 and 3 D 1 . However, at low energies the angular momentum suppression of the D-wave may make the system effectively elastic in S-wave.
For scattering with more than one partial-wave or hadron-hadron channel, there is no longer a one-to-one relation between E cm and elements of the t-matrix and we choose to make progress by parameterizing the energy dependence of t(E cm ). In order to calculate the scattering amplitudes using Eq. 11 we follow the successful approach detailed in [41]. In brief, taking an appropriate parameterization of the scattering t-matrix, we calculate the energy spectrum in each irrep using Eq. 11. By varying the free parameters in the parameterization, we find the best description of the finite volume spectra by minimizing a χ 2 , described in Eq. 9 of Ref. [42], measuring the agreement between finite-volume spectra obtained in the lattice calculations and those found by solving Eq. 11 with the parameterized t-matrix. To ensure that we have not introduced bias by any particular choice of t-matrix parameterization, we repeat the analysis for a range of parameterization forms, establishing which features of the resulting amplitudes are robust.
A very convenient approach to building parameterizations of the t-matrix is to work in terms of a real symmetric and I ab (s) = I a (s) δ ab is a matrix diagonal in hadronhadron channel. Unitarity of the S-matrix is guaranteed if Im I a (s) = −ρ a (s) above threshold in channel a and zero below. A simple choice is I a (s) = −iρ a (s). Alternatively, the Chew-Mandelstam prescription [56] defines Re I a (s) through a dispersive integral featuring ρ a (s)this has improved analytic structure as we transition across thresholds and move away from the real energy axis. A detailed discussion of our implementation can be found in Ref. [46]. One parameterization we utilize expresses the components of K −1 (s) as polynomials in s, where c (n) is a real symmetric matrix. Flexibility in this form comes from varying N and allowing parameter freedom in different combinations of c (n) Ja, Jb coefficients. An alternative approach is to parameterize the components of K(s) directly, using a parameterization of the form where m is a real parameter, g Ja (s) is some real polynomial in s, and γ (n) is a symmetric matrix of real parameters. These forms assume nothing about a nearby resonance or bound state but the pole can efficiently describe such behavior where it is present. These and similar K-matrix parameterizations have been successfully used in previous lattice QCD calculations of resonant and non-resonant scattering [13,14,17,38,43,45,46]. As an explicit example, one that we will make use of later, consider a K-matrix parameterization suitable for describing the dynamically-coupled J P = 1 + channels πω 3 S 1 , πω 3 D 1 and πφ 3 S 1 . 10 One possible choice, with 7 free parameters, is g πω{ 3 S1} g πω{ 3 D1} g πω{ 3 S1} g πφ{ 3 S1} g πω{ 3 S1} g πω{ 3 D1} g 2 πω{ 3 D1} g πω{ 3 D1} g πφ{ 3 S1} g πω{ 3 S1} g πφ{ 3 S1} g πω{ 3 D1} g πφ{ 3 S1} g 2 where this form allows mixing between πω and πφ channels only through g πφ{ 3 S1} . To include additional partial-waves that contribute as 10 In principle, we should also consider πφ in the 3 D 1 partial-wave; however, suppression due to the centrifugal barrier factor, compounded with strong OZI suppression of πφ, suggests it will be negligibly small and we find later in Section VIII that the amplitude is consistent with zero in the energy region we consider.
a consequence of the finite-volume but which do not mix in an infinite-volume, i.e. those with distinct J P as seen in Table III for irreps [000] T + 1 and P A 2 , we write the t-matrix in block-diagonal form with each block corresponding to a J P . We refer the reader to Ref. [17] for more details.
Statistical uncertainties on the scattering parameters and parameter correlations are determined by calculating the second derivatives of the correlated χ 2 at its minimum.
We make a conservative estimate of systematic uncertainties on each scattering parameter due to the uncertainties on stable hadron masses and the anisotropy by repeating the χ 2 minimization fitting procedure at all the various combinations of ξ ± δξ and m i ± δm i . 11 For each of these minimizations, we keep the finite-volume energies, E cm , their corresponding uncertainties, δE cm , and correlations between energy levels fixed, where and E lat is the energy in the lattice frame. For each scattering parameter, the largest change in the central value is quoted as its systematic uncertainty. Utilizing the approach outlined above, we now determine scattering amplitudes starting with a single partial wave using energy levels below πφ threshold, and progressing to a larger set of partial waves using the full set of energy levels.
A. "Elastic" πω 3 S1 scattering Below πφ threshold, the kinematically-open hadron channels are the two-body πω and three-body ππη. We expect ππη to become an important channel near the lowest E (2+1) n.i.
where the ρ resonance enhances ππ as discussed in Section II C. Below this energy, we expect the need to have a P -wave to get overall J P = 1 + will strongly suppress the amplitude. The lowest E (2+1) n.i.
in each of the irreps we consider is typically much higher in energy than πφ threshold, and so we will initially propose that we can ignore ππη.
In this energy region only slightly above πω threshold, the centrifugal barrier suppresses contributions of higherpartial waves, t J, J ∼ k + cm , such that we expect the 3 D 1 contributions to the coupled 3 S 1 , 3 D 1 partial-waves to be rather small. Similarly, πω scattering amplitudes in other partial-waves that appear in these irreps due to the finite-volume, as shown in Table III, are expected to be suppressed relative to the 3 S 1 amplitude and to have no significant resonant enhancement below πφ threshold. It follows that we can attempt an "elastic" analysis in terms of pure πω 3 S 1 → πω 3 S 1 scattering at low energy. We use 20 levels, all at least 1σ below the πφ threshold. Specifically, for each irrep these correspond to the lowest level on each of the (L/a s ) = 16 and 20 volumes and the 11 Values of the anisotropy, masses and uncertainties are given in Section III.  Table XIV of Appendix C with only the statistical uncertainties shown. The point size (small to large) of the discrete phase-shift point encodes the lattice volume (small to large).
lowest two levels on the (L/a s ) = 24 volume 12 , shown as the black points below πφ threshold in Figure 5. The resulting discrete phase-shift points are plotted in Figure 7, where we see that the trend is for them to increase toward a value close to 90 • as they approach the energy cutoff at πφ threshold. This is certainly consistent with a resonance located somewhere near to that energy. Instead of extracting discrete phase-shift points, we can also fit the spectrum using energy-dependent parameterizations of elastic scattering; a selection of choices which describe the finite-volume spectra well are included as gray curves in Figure 7 with the details of the parameterizations presented in Appendix C. One description, chosen as a reference amplitude and plotted as the blue curve in Figure 7, is given by, using the Chew-Mandelstam prescription for I(s) with Re I(s = m 2 ) = 0 -see Appendix B of Ref. [46]. The best fit description of the finite-volume spectrum is where the first uncertainty is statistical and the second is systematic as discussed above, and where the matrix shows the correlations between the parameters.
B. Dynamically-coupled πω 3 S1 , πω 3 D1 scattering Now we relax the assumption of negligible πω 3 D 1 contributions and perform a coupled-channel analysis on the dynamically-coupled πω 3 S 1 and πω 3 D 1 system, restricted to the same low energy region below πφ threshold as in Section VI A. Motivated by the suggestion of resonant behavior in the πω 3 S 1 phase-shift in the pre-vious section, we should allow for a resonance to have a πω 3 D 1 coupling as this could significantly enhance the πω 3 D 1 contribution above what might be expected on the basis of angular momentum suppression at threshold.
An example of a two-channel parameterization capable of describing the finite-volume spectra is The parameters m and g πω{ 3 S1} are compatible with those of the reference amplitude in Eq. 19 and we find g πω{ 3 D1} to be consistent with zero within uncertainties. In Figure 8 we present the πω 3 S 1 and πω 3 D 1 phase-shifts and thē (πω 3 S 1 |πω 3 D 1 ) mixing-angle as defined in the Stappparameterization [57] and given in Eq. A7 of Appendix A. A number of different K-matrix parameterizations were explored and are plotted as the gray curves in Figure 8 and listed in Table XV of Appendix C. We observe that all descriptions exhibit a πω 3 S 1 phase-shift compatible with the behavior seen in Section VI A, a πω 3 D 1 phase-shift that is very small, and a mixing-angle that is consistent with zero within a modest uncertainty over this energy range.
C. Coupled πω 3 S1 , πω 3 D1 , πφ 3 S1 scattering We now consider scattering amplitudes in an energy region up to the ππππ-threshold. In this region, πω, ππη, πφ, πKK and ππσ are all kinematically open, but we expect the three-body channels to have only a small effect. By using only energy levels below the lowest E (2+1) n.i. or E n.i. in each irrep, and excluding any energy levels which show significant sensitivity to the presence of ρη, K * K and a 0 π operators, we propose that we can effectively neglect the effect of three-body channels. In Section VIII we will explore possible effects of relaxing this assumption.
We proceed with a total of 36 energy levels -all the black points shown in Figure 5.
Both πω and πφ are vector-pseudoscalar channels dynamically-coupled in 3 S 1 and 3 D 1 partial-waves. However, considering the centrifugal barrier for the heavier threshold and the lack of mixing observed in the histograms presented in Figure 4, we assume that πφ 3 D 1 will have negligible impact at low energies. Subsequently, we are left with a system of three coupled channels: πω 3 S 1 , πω 3 D 1 and πφ 3 S 1 . Many other partial-waves can contribute to the finite-volume spectra as can be seen from Tables III and V, but, as discussed in Section IV, we expect these to be negligibly small and we will explicitly show this in Section VIII.
To parameterize the energy dependence of the threechannel t-matrix, we use K-matrices of the form in Eq. 15 restricted to linear expansions in g Ja (s) and γ -the parameterizations used are presented in full in Table XVI of Appendix C. It should be noted that, while use of the K-matrix guarantees unitarity, it does not guarantee good analytic properties. Indeed, we found that some parameterizations, which successfully describe the finite-volume spectra, have t-matrix pole singularities at complex energies on the physical sheet. Such poles are forbidden by causality, and these parameterizations must be rejected as giving rise to unphysical solutions. A list of such parameterizations is provided in Table 0 A somewhat minimal parameterization, used with the Chew-Mandelstam prescription with Re I a (s = m 2 ) = 0, proves to be capable of the describing the finite-volume spectra. The best-fit parameters are πφ{ 3 S1},πφ{ 3 S1} = 0.90 ± 0.24 ± 0.27 We found no improvement in the description of the finitevolume spectra by including freedom in g πφ{ 3 S1} and subsequently fixed this parameter to be zero in the reference amplitude.
There is no established method to minimally display the S-matrix in three-channel scattering. Plotting the real and imaginary parts of the elements of the S-matrix contains redundancy as it does not account for the constraints provided by unitarity. Plotting the magnitudes via ρ a ρ b t ab 2 has the advantage of being closely related to a differential cross-section, but discards important phase information. In the two channel case, the Stapp parameterization is minimal with regard to unitarity and reduces to single-channel phase-shifts when the channels decouple, but to our knowledge there is not a generalization to more channels that reduces to the two-channel Stapp parameterization. In Appendix A we provide such a generalization to n-channels where, if k are decoupled, the scattering S-matrix naturally block diagonalises into an (n − k) coupled-channel block and a diagonal block containing k decoupled phase-shifts. The phase-shifts and mixing-angles are plotted in Figure 9 for the amplitude in Eqs. 21 and 22 (colored curves) and the many other parameterizations listed in Table XVI of Appendix C (gray curves). We observe that the behavior of the πω 3 S 1 phase-shift is in close agreement with the results of Section VI B, and the πω 3 D 1 phase-shift is once again very small. The πφ 3 S 1 phase-shift shows a small positive tendency indicative of a weak attraction. The mixing-angle¯ (πω 3 S 1 |πω 3 D 1 ) is small but likely non-zero, while the mixing angles¯ (πω 3 S 1 |πφ 3 S 1 ) and¯ (πω 3 D 1 |πφ 3 S 1 ) are around two orders of magnitude smaller and statistically consistent with zero everywhere.
The same amplitudes are plotted as ρ a ρ b |t Ja, Jb | 2 in Figure 10. We observe a significant bump-like enhancement in the πω 3 S 1 → πω 3 S 1 process which would be a canonical indication for a resonance in a scattering cross-section measurement.
In Figure 11 we present the energies calculated using Eq. 11 with the reference amplitude of Eqs. 21 and 22 which, as suggested by the small χ 2 , are seen to be in good agreement with the lattice finite-volume energy levels. Notably, for levels not included in the fits, shown in gray, the predicted spectra on the (L/a s ) = 20, 24 volumes appear to be mainly in reasonable agreement, while on the (L/a s ) = 16 volume there is a larger discrepancy. This may be attributed to more significant contributions from three-meson amplitudes on smaller volumes, further supported by the observation that there is a much larger variation in the spectrum in the [000] T + 1 irrep on the smaller volume when three-meson like operators are removed -see Figure 4.
A final comment concerns the effect on the scattering results of the uncertainty placed on the anisotropy due to the observed dependence on vector-meson helicity in Section III. Unlike in the ρπ isospin-2 case presented in Ref. [17], where the weak nature of the scattering led to the anisotropy uncertainty being the largest systematic effect, here the interactions are strong and the anisotropy uncertainty contributes relatively little as can be seen from the relative sizes of the inner and outer bands in Figures 9 and 10.
To summarise, the characteristic "bump" we found in the scattering magnitudes in Figure 10 and the clearly observed avoided level crossing in the [111]A 2 spectrum seen in Figure 5 strongly suggests a resonance. To demonstrate this rigorously, we proceed to determine the pole singularities of our scattering amplitudes.  Table XV of Appendix C with only statistical uncertainties shown. Middle: As upper but for the mixing-angle,¯ (πω 3 S1 |πω 3 D1 ). Lower: Black points are the finite-volume energy levels used to constrain the fit and orange points are the energy levels calculated using Eq. 11 for the reference amplitude in Eq. 20.

VII. POLE ANALYSIS FOR COUPLED-CHANNEL AMPLITUDES
At each threshold, unitarity necessitates a branch point singularity and the corresponding branch cut divides the complex s-plane into two Riemann sheets, so for n open thresholds there are 2 n sheets. Riemann sheets can be labelled by the sign of the imaginary component of the cmframe momentum k (a) in each hadron channel a. We identify the physical sheet, where physical scattering occurs just above the real energy axis, as having Im(k (a) ) > 0 for all a. Sheets with other sign combinations are referred to as unphysical, and it is on these sheets that pole singularities corresponding to resonances lie, in complex-conjugate pairs, off the real energy axis. Poles off the real axis on the physical sheet indicate causality violating amplitudes and signal an unacceptable description of the scattering process.
For poles off the real axis, we define the real and imaginary parts of the pole singularity at s = s 0 in terms of the mass m R and the width Γ R of a resonance respectively, by √ s 0 = m R ± i 2 Γ R . For narrow resonances, with a single dominant decay mode, these definitions of the resonance mass and width agree well with the location Upper: As in Figure 8 but for the πω 3 S1 (blue), πω 3 D1 (purple) and πφ 3 S1 (green) phase-shifts for the reference amplitude in Eqs. 21 and 22, and for other parameterizations presented in Table XVI of Appendix C (gray). Middle: As upper but for the mixing-angle¯ (πω 3 S1 |πω 3 D1 ). The other mixing-angles, (πω 3 S1 |πφ 3 S1 ) and¯ (πω 3 D1 |πφ 3 S1 ), are extremely small and consistent with zero for all parameterizations and are not plotted. Lower: The energy levels used to constrain the scattering amplitude (black) and their corresponding description by the amplitude in Eqs. 21 and 22 (orange). and full-width at half-maximum of the "bump" seen in scattering cross sections. The advantage of associating the pole singularity with the resonance is that this definition is still useful in complicated coupled-channel cases, such as those seen in the lattice calculations of the a 0 [13] and f 0 [14], where the resonance does not appear as a clear isolated bump for real energies.
In the current case, the hadron-hadron channels πω and πφ lead to four sheets, (sign(Im k πω ), sign(Im k πφ )) = I(+, +), II(−, +), III(−, −), IV(+, −) . Close to the πφ threshold, all of sheets II (lower half-plane), III (lower half-plane) and IV (upper half-plane) are close to physical scattering. A single resonance can appear as a pole in slightly different positions on multiple sheets -some discussion of this in the context of a simple coupled-channel amplitude model can be found in Ref. [13].
For complex energies close to a pole singularity at s 0 , the scattering t-matrix can be written in the factorised form  Table XVI of Appendix C are in gray with bands reflecting only the statistical uncertainties. ρaρ b |t Ja, Jb | 2 not plotted are significantly smaller than those shown and consistent with zero.
where the complex valued couplings c Ja reflect the strength of the resonance coupling to channel a{ 3 J }. For each coupled hadron-hadron channel, the coupling is determined only up to a sign which gives no change to the physics. In the current case this leads to a sign ambiguity between the πω and πφ couplings, but conversely the relative sign between the 3 S 1 and 3 D 1 partial-waves in πω can be unambiguously determined and physically would lead to different angular decay shapes depending on its value. In Ref. [17], it was shown that in a finite volume, moving-frame spectra are required to constrain this sign.
For each amplitude parameterization we considered, using the best-fit values of parameters, we perform a search across all Riemann sheets over a large range of complex s, finding any pole singularities present and determining the couplings by factorizing the residue of the pole. Uncertainties on the pole positions and couplings are estimated by appropriately propagating through the uncertainties and correlations on the fit parameters. For the case of the reference amplitude presented in Eqs. 21 and 22, poles were found in complex conjugate pairs on sheet II at a t √ s 0 II = 0.2435(13)(10) ± i 2 0.0175(20) (19), (24) where the first uncertainty is statistical and the second is systematic. A complex conjugate pair of poles was also found on sheet III in agreement with Eq. 24 up to the precision shown. The couplings for the pole in the lower half-plane are a t c(πω 3 S 1 ) II = 0.106 (6)(6) and c(πφ 3 S 1 ) II is exactly zero, a result of the choice of reference amplitude. Considered as a ratio we have (37)(20) arg c(πω 3 D 1 ) II /c(πω 3 S 1 ) II = −π 0.103(26) (24).
That the poles on sheets II and III are in essentially the same position is a consequence of the πφ channel being almost completely decoupled from the πω channel as discussed in Section V. For each three-channel parameterization presented in Table XVI of Appendix C, we found poles and couplings broadly consistent with those given above. We show these in Figure 12, observing that the scatter over different parameterizations is in this case not significantly larger than the uncertainty on the reference amplitude.

VIII. SYSTEMATIC TESTS
To test the robustness of the extracted scattering amplitudes and the determination of the resonant pole and couplings, we consider two sources of potential systematic uncertainties due to possibilities we have so far neglected. First, we examine the partial-waves that mix as a consequence of the finite-volume, which we neglected based on observations discussed in Section V, and the πφ{ 3 D 1 } amplitude which we asserted was negligible. Second, we examine the dependence of the energy levels on the πω{ 3 D 1 }, πω{ 3 P 0 } and πω{ 3 P 2 } parameters to demonstrate that we are able to constrain these amplitudes. Lastly, we make a crude estimate of the possible size of effects due to the neglected three-body channels.

A. Additional Partial-Waves
We first consider the πω{ 3 P 0 } and πω{ 3 P 2 } amplitudes that enter in the P A 2 irreps as shown in Table III. Since a P -wave has less threshold suppression than a D-wave, we might expect these waves to be at least as important as πω 3 D 1 , though they are not expected to be resonant at such low energies. Augmenting the reference amplitude as defined in Eq. 21, we allow a non-zero amplitude in the πω{ 3 P 0 } and πω{ 3 P 2 } channels by including a constant γ-term for each in the K-matrix and for these additional channels we set Re I a (s = (m π +m ω ) 2 ) = 0 in the Chew-Mandelstam phase-space. The resulting t-matrix is block diagonal in J P reflecting the fact that this mixing is a result of the reduced symmetry on the lattice. We fit to the same 36 energy levels as in Section VI C and, allowing all parameters to vary, find m = (0.2466 ± 0.0007) · a −1 t g πω{ 3 S1} = (0.105 ± 0.007) · a −1 t g πω{ 3 D1} = (1.12 ± 0.46) · a t γ (0) πω{ 3 S1},πω{ 3 S1} = −0.34 ± 0.19 γ (0) πω{ 3 S1},πφ{ 3 S1} = 0.79 ± 0.25 γ (0) where correlations between the πω 3 S 1 , πω 3 D 1 and πφ 3 S 1 parameters are compatible with those shown in Eq. 22, and correlations between these and πω{ 3 P 0 } and πω{ 3 P 2 } parameters are small. We observe that the amplitudes in both πω{ 3 P 0 } and πω{ 3 P 2 } are consistent with zero. A similar approach allowing for πω{ 3 D 2 } and πω{ 3 D 3 } parameter freedom finds no evidence for large amplitudes as one would expect given the larger angular momentum suppression and lack of low-energy resonances with J P C = 2 +− and 3 +− . In order to investigate the possible effect of the previously excluded πφ{ 3 D 1 }, we take the reference amplitude in Eq. 21 and extend it to include a constant diagonal γ-term in πφ{ 3 D 1 } in the K-matrix. Once again fitting to the 36 energy levels and allowing all parameters to vary, we find the πφ{ 3 D 1 } parameter to be consistent with zero, as expected, with all other parameters compatible with those presented in Eq. 22.
B. Spectrum dependence on πω{ 3 P0}, πω{ 3 P2} and πω 3 D1 It is worth illustrating at this stage how particular energy levels in the finite-volume spectra depend upon the strength in πω 3 D 1 and the πω P -waves. For πω 3 D 1 this is shown in Figure 13, where the curves present the finite-volume energy spectrum for the reference amplitude in Eqs. 21 and 22, varying the value of g πω{ 3 D1} while keeping all other parameters fixed. In each irrep, we see a level near the lowest πφ non-interacting energy which appears to be independent of the value of g πω{ 3 D1} , as expected given the near complete decoupling of πφ.  Table XVI of Appendix C. Black ellipses correspond to the reference amplitude in Eq. 22 reflecting the statistical (inner) plus systematic (outer) uncertainties. Bottom: As top but for the corresponding couplings, c(πω 3 S1 ) II (blue), c(πω 3 D1 ) II (purple) and c(πφ 3 S1 ) II (green). Black ellipses again correspond to the couplings of the reference amplitude in Eq. 22 where c(πφ 3 S1 ) II = 0.
are levels observed to be consistent with the two-fold degenerate non-interacting πω energies, which show no visible dependence on g πω{ 3 D1} .
Interestingly, the position of these same levels proves to be strongly dependent on the amplitude strength in the πω{ 3 P 0 } and πω{ 3 P 2 } partial-waves, so the lattice computed energies allow us to confidently limit the amplitude of these P -waves to be very small in this energy region. Figures 14 and 15 show the analogue of Figure 13 but for varying πω{ 3 P 0 } and πω{ 3 P 2 } channel parameters respectively. In these two cases, the reference amplitude in Eq. 21 is augmented, as described in Section VIII A, to include a constant γ-term in the K-matrix for channels πω{ 3 P 0 } and πω{ 3 P 2 }.

C. Three-body channels
For the light-quark masses used in this calculation, the resonant behavior is found to occur between the relatively low-lying ππη threshold and the somewhat higher-lying πKK threshold. As such, we might worry that these channels could have a significant impact on the physics in this region. We previously presented some arguments for why we do not expect this to be the case, but noted that in Figure 4 there appeared to be deviations in the finite-volume spectra depending on whether or not threemeson operators were included in the bases, most notably in the smallest, (L/a s ) = 16, volume. As a precaution, we ensured that we only made use of those energy levels which lie below the lowest E (2+1) n.i.
value and which show no significant dependence on the presence/absence of ρη, K * K or a 0 π operators.
In this section, we attempt to quantify the size of possible contributions from the three-body sector on our scattering amplitudes and resonance pole by treating the scattering system as though ππ in ππη can be completely replaced by a stable ρ (with a fixed mass a t m ρ = 0.1509) and πK in πKK can be completely replaced by a stable K * (a t m K * = 0.1648). In this way we augment our scattering matrix with two extra channels ρη 3 S 1 and K * K 3 S 1 .
This approach cannot be expected to completely describe the finite-volume spectra because, for example, whenever the ρ has non-zero momentum, we expect there to be more than one corresponding energy level, as indicated by Figure 1 in Ref. [42] -a stable ρ model cannot capture this and will not even give the right number of energy levels in the 'three-body' spectrum. However, for the ρ at rest the nearest non-interacting ππ energy is much higher, and there is effectively only one finite-volume level which lies very close to the ρ resonance mass. In this case, the stable ρ may be a reasonable first approximation to the true three-body physics.
As a final test of the effects of the ρη{ 3 S 1 } and K * K{ 3 S 1 } channels, we find the resonance pole and corresponding couplings. There are 16 Riemann sheets and several 'mirror poles', but the closest pole is located at where we exclude the meaningless phase on a t c(ρη{ 3 S 1 }) II as the magnitude is consistent with zero and where c(πφ{ 3 S 1 }) II = 0 by choice of amplitude. The coupling to ρη{ 3 S 1 } is small but has a large uncertainty, while the coupling to K * K{ 3 S 1 } is larger. We conclude that although we cannot currently rigorously handle three-body contributions due to ππη and πKK, we do not see any evidence to suggest that they significantly affect the results reported in this paper.

IX. INTERPRETATION
All J P = 1 + amplitude parameterizations used, that prove to be capable of describing the finite-volume spectra in the energy region we are considering, had the same characteristic resonant bump in the πω 3 S 1 to πω 3 S 1 amplitude squared, with little strength in the diagonal πω 3 D 1 and πφ 3 S 1 elements. The off-diagonal amplitudes were all found to be relatively small. In every case we found that the bump is associated with a complex conjugate pair of poles on sheets II and III, which we interpret as the effect of a single resonance. As in previous calculations, to quote results in physical units, we choose to set the scale using the Ω-baryon mass measured on these lattices, a t m Ω = 0.2951 [58], and the physical Ω-baryon mass, m phys Ω = 1672 MeV [15]. This gives a −1 t = m phys Ω /(a t m Ω ) = 5666 MeV and stable hadron masses m π ≈ 391 MeV, m K ≈ 549 MeV, m η ≈ 587 MeV, m ω ≈ 881 MeV and m φ ≈ 1017 MeV.
Using this scale setting, we summarise the scattering amplitudes resulting from this work in Figure 17 In Figure 18 we plot the position of the pole found in this calculation compared to the experimental b 1 resonance, with mass m b1 = 1230(3) MeV and width Γ b1 = 142(9) MeV [15], and a lattice calculation at the SU(3) F point with m π ≈ 700 MeV [30]. In the latter calculation, the b 1 forms part of an axial-vector octet with mass around 1525 MeV; the pseudoscalar-vector threshold corresponding to πω is at roughly 1695 MeV, and thus the b 1 is stable at this pion mass. We observe that the trajectory of the pole with varying pion mass appears to be similar to that of the ρ meson shown in Ref. [38], as may be expected for a reasonably narrow resonance.
Since we find the b 1 to be a narrow resonance a moderate distance above πω threshold, it is reasonable to compute theoretical 'branching fractions' for its decay to πω. For channels a{ 3 J } these are given by [15], As mentioned in Ref. [14], the sum of these partial branching fractions does not necessarily give unity. We obtain where using the definition in Eq. 30 the πφ{ 3 S 1 } branching fraction is zero as the channel is kinematically closed (m R < m π + m φ ).
A crude extrapolation of the couplings to the physical value of the light quark masses comes if we assume them to be independent of the light quark masses once the threshold behavior is removed. Such a behavior is not guaranteed, but has been observed in lattice calculations Top: The scattering amplitudes-squared, ρaρ b |t Ja, Jb | 2 , transcribed from Figure 10 with the energy axis converted to physical units. Below the amplitudes are the energy levels used to constrain the amplitudes (black points). Bottom: The best estimate of the resonant pole position, where uncertainties combine statistical and systematic uncertainties with variations across parameterizations. The histograms show the best estimate of the magnitude of each coupling with the lightly-shaded region reflecting the combined uncertainties. The πφ{ 3 S1} coupling is an estimate of the upper bound.
where the cm-frame momentum is evaluated at the resonance pole position, and where we use the values presented above on the right-hand side, and the experimental b 1 mass to compute k phys. πω , gives a prediction of c phys. πω{ 3 D1} = 146(101) MeV. Subsequently, we obtain an estimate for the ratio of couplings at the physical pion mass of, c phys.
c phys. πω{ 3 S1} = 0.27 (20). Blue shows the ground-state mass of the axial-vector octet from a lattice calculation with mπ ≈ 700 MeV [38], red shows the estimate from this work with mπ ≈ 391 MeV and black is the experimentally determined mass and width of the b1 resonance [15].
The PDG [15] reports a ratio of D-wave to S-wave amplitudes for the b 1 resonance of magnitude 0.277 (27), which is not computed at the complex pole position and therefore not precisely the same quantity as we quote.

X. SUMMARY
This paper has reported on the first lattice QCD calculation of coupled πω, πφ scattering, the first time coupled pseudoscalar-vector scattering amplitudes have been computed. This large-scale calculation made use of a significant number of operators resembling single, two and three-meson constructions to extract finite-volume spectra which were used to constrain the coupled-channel scattering amplitudes.
Analysis of the obtained finite-volume spectra required consideration of coupled 3 S 1 − 3 D 1 partial-waves in πω scattering. A clear b 1 resonance was observed, visible as a rapid increase in the πω 3 S 1 phase-shift through 90 • or correspondingly as a bump in the magnitude of the πω 3 S 1 → πω 3 S 1 t-matrix element. More rigorously, we found pole singularities on unphysical Riemann sheets relatively close to the real energy axis with couplings that are large for the πω 3 S 1 final state, significantly smaller for πω 3 D 1 and compatible with zero for πφ. The mass and width of the b 1 resonance found in this calculation, with light-quark masses such that m π ≈ 391 MeV, appear to be compatible with a smooth interpolation between a stable state for much larger quark mass, and the experimental resonance at lower quark mass.
We explored the role of three-body channels by including operators in our bases whose construction resembles a meson coupled to a two-body resonance, utilizing earlier calculations of mesonmeson scattering channels [13,42,43]. There is no sufficiently-mature finite-volume formalism capable of rigorously incorporating three-body scattering channels here. However, as a systematic test, the finite-volume formalism, which in principle can handle any number of coupled meson-meson channels, was applied in a limited study of five coupled-channels πω 3 S 1 , πω 3 D 1 , πφ 3 S 1 , ρη 3 S 1 , K * K 3 S 1 . Our investigations suggested that the three-body channels have a negligible effect in this particular case of a low-lying b 1 resonance. Furthermore, observations were made of how particular finite-volume energy levels depend upon the various partial-waves which 'mix' due to the cubic nature of the lattice boundary.
In order to provide a way to minimally present nchannel scattering on the real energy axis, a generalization of the two-channel Stapp parameterization was presented in which a unitary S-matrix is expressed in terms of n phase-shifts and n(n−1)/2 mixing-angles. This parameterization was used to present the three-channel πω 3 S 1 , πω 3 D 1 , πφ 3 S 1 J P = 1 + scattering matrix in which the b 1 resonance appears. The construction provided conveniently reduces to the Stapp form in the case that one channel decouples from the others (as approximately found here).
As expected, no I G = 1 + resonances are observed with a mass comparable to the b 1 in J P = 0 − , 2 − . Notably, no resonating behavior is observed in a largely decoupled πφ channel, suggesting the absence of a Z s state which might be proposed as an analogue of the Z c state claimed in πJ/ψ. This work has advanced lattice techniques for studying coupled-channel scattering involving hadrons with nonzero spin and operators which effectively interpolate three hadrons. Looking forward, once a three-hadron scattering formalism is practical to use, a future calculation would enable the rigorous determination of the ππη and πKK scattering amplitudes. Furthermore, utilizing such a formalism would allow the calculation of the G-paritynegative axial-vector, the a 1 , which has a dominant decay to the pseudoscalar-vector meson pair πρ, for which the ρ is unstable at this pion mass, and would make for an interesting comparison. Moving on from the simplest lowlying resonances, and as the light-quark mass approaches its physical value, it becomes more important to reliably determine such three-hadron scattering processes.
These conventions mean that δ 1 is equal to arg(S 11 ), which is in agreement with the conventions in Refs. [13,14,45] where the phase-shift is defined as δ i = arg(S ii ). However, we see for δ 2 and δ 3 there are corrections to the phase due to the imaginary components ∝ σ 12 s 13 s 23 c 23 in the expressions for S 22 and S 33 , given in Eq. A8. For a very weakly mixed channel these corrections are very small and δ i ≈ arg(S ii ) for i = 2, 3.

n = 5
For the limited five coupled-channel analysis of πω 3 S 1 , πω 3 D 1 , πφ 3 S 1 , ρη 3 S 1 , K * K 3 S 1 given in Section VIII C, we calculate the five phase-shifts and ten mixing-angles. We find that seven of the mixingangles, all featuring either πφ 3 S 1 and/or ρη{ 3 S 1 }, are extremely small and consistent with zero, in agreement with the observation that both are decoupled from the resonance as shown in Eq. 29. This illustrates the natural reduction from the five-channel parameterization to the three-channel parameterization in the case that two channels decouple. The five phase-shifts and the remaining three non-zero mixing-angles are presented in Figure 19.
TABLE VI. Single-meson and two-meson operators used to compute optimised ρ operators in the [000]T − 1 irrep and P A1 irreps at various overall momenta on the three volumes. Momentum labels on the π's that form the ππ operators are omitted for brevity. Upper: As in Figure 8 but for the πω 3 S1 (blue), πω 3 D1 (purple), πφ 3 S1 (green), ρη{ 3 S1} (orange) and K * K{ 3 S1} (red) phase-shifts for the reference amplitude in Eq. 27. The faded error bands reflect the statistical uncertainty on the scattering parameters. The ρη and K * K "thresholds" are calculated using the ρ and K * masses given in Section VIII C. Lower: As upper but for the mixing-angles (πω 3 S1 |πω 3 D1 ) (blue),¯ (πω 3 S1 |K * K{ 3 S1}) (gray) and¯ (πω 3 D1 |K * K{ 3 S1}) (brown). All other mixing-angles are extremely small and consistent with zero as discussed in the text.    2 ) 2 . Otherwise, we set Ia(s) = −iρa(s). Results of fits to these parameterizations can be found in the Supplemental Material.