# Toward $N$ to $N\pi$ matrix elements from lattice QCD

Lorenzo Barca<EMAIL_ADDRESS>Gunnar Bali<EMAIL_ADDRESS>Sara
Collins<EMAIL_ADDRESS>Fakultät für Physik, Universität Regensburg,
Universitätsstr. 31, 93053 Regensburg, Germany

(28th August 2024)

###### Abstract

QCD matrix elements of axial and vector currents between nucleons are required
for the Monte Carlo reconstruction of the energy of neutrinos that are
detected in long baseline oscillation experiments in the quasielastic regime.
The cleanest approach for determining the axial matrix elements is lattice
QCD. However, the extraction of these from the corresponding correlation
functions is complicated by very large excited state contributions, that are
related to transitions from the nucleon to a nucleon-pion pair. In this pilot
study with a pion mass $m_{\pi}=429\leavevmode\nobreak\ \mathrm{MeV}$, we
demonstrate for the first time that these contributions can be removed by
including five-(anti)quark operators into the basis of interpolators used to
create the nucleon. The same techniques will be needed to compute transition
matrix elements between the nucleon and nucleon-pion scattering states that
are relevant in the resonance production regime.

††preprint: APS/123-QED

##  I INTRODUCTION

The groundbreaking discovery of atmospheric and solar neutrino oscillations
more than two decades ago by the Super-Kamiokande [1] and SNO [2] experiments,
respectively, required an adjustment of the Standard Model to accommodate
massive neutrinos. The present generation of terrestrial long baseline
neutrino oscillation experiments, aimed at a more precise determination of the
neutrino masses and mixing parameters, NOvA [3] and T2K [4] as well as the
future DUNE [5] experiment and the upgrade of T2K to the Hyper-Kamiokande
detector [6] determine the fluxes of muon and antimuon neutrinos via their
interaction with nuclear targets in a near and a far detector. The neutrino
energies of the scattering events are reconstructed via Monte Carlo event
generators [7, 8], which require knowledge of the differential neutrino-
nucleon cross section. For neutrino energies below $1\leavevmode\nobreak\
\mathrm{GeV}$ this is dominated by (quasi-)elastic scattering, while from
about $400\leavevmode\nobreak\ \mathrm{MeV}$ onwards also resonance production
with nucleon-pion ($N\pi$) final states sets in [9]. Focusing on low energies,
the cross section is proportional to the square of a combination of
nonperturbative nucleon vector and axial matrix elements, which can be
parameterized in terms of form factors. The two vector form factors as
functions of the squared four-momentum transfer ($Q^{2}$) are sufficiently
well known from experiment. However, the two isovector axial form factors
$G_{A}(Q^{2})$ and $\widetilde{G}_{P}(Q^{2})$ are much less well constrained
experimentally, apart from $G_{A}$ in the forward limit ($Q^{2}=0$, i.e. the
axial charge $g_{A}$ [10]) and
$\widetilde{G}_{P}(0.88\,m_{\mu}^{2})=g_{P}^{*}$ at the muon capture point of
muonic hydrogen [11]. Fortunately, these form factors can be computed directly
from QCD via lattice simulation. However, there is a tension [12, 13] between
recent lattice results [14, 15, 16, 17, 18] and analyses of neutrino-deuteron
scattering experiments [19]. Therefore, it is important to establish the
reliability of the lattice determinations. This requires the investigation of
all systematics and, in particular, the one associated with extracting the
nucleon matrix elements from correlation functions at finite Euclidean times.
The latter receive contributions also from single- and multiparticle states
with the same quantum numbers as the nucleon (normally referred to as excited
states). At zero momentum, the lowest excitations with positive parity include
$N\pi$ P-wave and $N\pi\pi$ S-wave scattering states, whereas, at nonvanishing
momentum, parity is not a good quantum number and also $N\pi$ in an S-wave can
contribute. Towards small pion masses, the mass gap between the ground state
and the first excitation decreases and the spectrum becomes more dense.
Bearing in mind that the signal-to-noise ratio of correlation functions
decreases exponentially with the Euclidean time separations, it can be very
challenging to reliably extract nucleon matrix elements. In order to control
the leading excited state contributions to these nucleon to nucleon form
factors, we will, for the first time, explicitly calculate matrix elements
that are also related to nucleon to $N\pi$ transition form factors, which are
required for a firm understanding of the resonance production regime.

Reliable continuum limit results for the axial form factors should reproduce
the experimentally known values of $g_{A}$ and $g_{P}^{*}$ and also be
consistent with the partially conserved axial current (PCAC) relation (also
referred to as the axial Ward identity (AWI)), which relates the axial form
factors to the pseudoscalar form factor. In many previous simulations $g_{A}$
was reproduced, however, $g_{P}^{*}$ (defined at $Q^{2}>0$) was found to be
smaller than the experimental value and also the AWI between form factors was
significantly violated [20, 21, 22, 23, 16]. Since the AWI was found to be
satisfied on the level of the correlation functions in the continuum limit
[24], the inconsistency had to be related to the difficulty of isolating the
ground state contribution when extracting the form factors [24, 15, 14]. While
the interpolator that is used to create the nucleon was found to have little
overlap with excited states, as evidenced by analyses of two-point functions,
transition matrix elements between different states, contributing to the
spectral decomposition of the three-point function, appeared to be enhanced.
Indeed, in chiral perturbation theory (ChPT) the axial and pseudoscalar
currents directly couple to the pion. Regarding the pseudoscalar current or
the time-component of the axial current, $N$ to $N\pi$ transitions can
contribute substantially to the three-point functions [25, 26, 27] (see also
Refs. [28, 29]). At a small but nonvanishing momentum, the leading such
contribution increases in proportion to the ratio of the nucleon mass over the
pion energy, $m_{N}/E_{\pi}$ [25, 26]. These terms were taken into account in
recent analyses of the Euclidean time dependence of lattice correlation
functions, where form factors were obtained, that are consistent with the AWI
[14, 15, 17]. However, the size of the excited state contamination, found in
these analyses, is quite large in some channels for the Euclidean times that
are accessible at present. A more reliable approach would be to construct
optimized interpolators to minimize the dominant ($N\pi$) excited state
contributions.

In this work, we take into account directly the $N\pi$ contribution by
constructing nucleon-pion-like interpolators [$(qqq)(\bar{q}q)$ with the
quarks $q\in\\{u,d\\}$] $\mathrm{O}_{5q}$, and computing the associated two-
point and, for the first time, three-point correlation functions between the
standard three-quark nucleon interpolator $\mathrm{O}_{3q}$ and
$\mathrm{O}_{5q}$. Using this basis, that has good overlap both with the
nucleon ground state and the lowest lying $N\pi$ state, nucleon to nucleon
three-point functions can be constructed with minimized $N\pi$ contributions,
enhancing the reliability of the extraction of the nucleon matrix elements. As
mentioned above, this is the first step towards determining nucleon to
nucleon-pion matrix elements, associated with neutrino scattering in the
resonance production regime [13, 30, 9]. In this pilot study we carry out the
analysis for a single unphysical pion mass $m_{\pi}=429\leavevmode\nobreak\
\mathrm{MeV}$. It turns out that even at this relatively large value the
$N\pi$ contribution is very significant and that this can effectively be
removed with our approach. We expect this method to work even better at the
physical pion mass: ChPT becomes more reliable as the pion mass is reduced,
and the tree-level $N\pi$ contribution is even more dominant, which is
consistent with the observations made in Refs. [14, 15, 25, 27, 26].

##  II DEFINITION OF THE FORM FACTORS

We define local isovector pseudoscalar and axial currents,
$\mathcal{P}=\bar{d}\gamma_{5}u$ and
$\mathcal{A}_{\mu}=\bar{d}\gamma_{\mu}\gamma_{5}u$, respectively. The Lorentz
decompositions into form factors of the respective matrix elements read as

$\displaystyle\langle
n_{\mathbf{p}^{\prime}}|\mathcal{P}|p_{\mathbf{p}}\rangle$
$\displaystyle=\bar{u}_{\mathbf{p}^{\prime}}G_{P}(Q^{2})\gamma_{5}u_{\mathbf{p}},$
(1) $\displaystyle\\!\\!\langle
n_{\mathbf{p}^{\prime}}|\mathcal{A}_{\mu}|p_{\mathbf{p}}\rangle$
$\displaystyle=\bar{u}_{\mathbf{p}^{\prime}}\\!\\!\left[\gamma_{\mu}G_{A}(Q^{2})+\frac{q_{\mu}}{2m_{N}}\widetilde{G}_{P}(Q^{2})\right]\\!\\!\gamma_{5}u_{\mathbf{p}},$
(2)

where we assume isospin symmetry (i.e. $m_{N}=m_{p}=m_{n}$ and
$m_{\ell}=m_{u}=m_{d}$ for the quark masses), $u_{\mathbf{p}}$ is the spinor
of a nucleon with three-momentum $\mathbf{p}$,
$q_{\mu}=p^{\prime}_{\mu}-p_{\mu}$ is the four-momentum transfer and
$Q^{2}=-q_{\mu}q^{\mu}$. Note that the above decomposition of the axial matrix
element does not hold if the two states differ in their mass, e.g., if a
nucleon is on the right-hand side and a $N\pi$ on the left. The AWI
$\partial_{\mu}\mathcal{A}^{\mu}=2im_{\ell}\mathcal{P}$ implies the relation
between form factors,
$m_{N}G_{A}(Q^{2})=m_{\ell}G_{P}(Q^{2})+(Q^{2}/4m_{N})\widetilde{G}_{P}(Q^{2})$,
which is exact in the continuum limit but will be affected by moderate
discretization effects at our lattice spacing $a\approx 0.098$ fm [24, 14]. In
addition, the pion pole dominance (PPD) assumption gives the approximate
relation, $\widetilde{G}_{P}(Q^{2})\approx
4m_{N}^{2}G_{A}(Q^{2})/(m_{\pi}^{2}+Q^{2}).$ While this only holds exactly for
$m_{\pi}=0$, in Ref. [14] it was found to hold within uncertainties of 1%–2%
at the physical point in the continuum limit, with violations of less than 3%
up to $m_{\pi}\approx 420$ MeV. Deviations from these relations can be
quantified in terms of the differences from unity of the combinations

$\displaystyle r_{\rm{PCAC}}$
$\displaystyle=\frac{4m_{N}m_{\ell}G_{P}(Q^{2})+Q^{2}\widetilde{G}_{P}(Q^{2})}{4m_{N}^{2}G_{A}(Q^{2})},$
(3) $\displaystyle r_{\rm{PPD}}$
$\displaystyle=\frac{(m_{\pi}^{2}+Q^{2})\widetilde{G}_{P}(Q^{2})}{4m_{N}^{2}G_{A}(Q^{2})}.$
(4)

##  III ANALYSIS

We construct the matrices of two- and three-point correlation functions (see
the Supplemental Material),

$\displaystyle C_{2pt}(\mathbf{p},t)_{ij}$
$\displaystyle=\left\langle\mathrm{O}_{i}(\mathbf{p},t)\leavevmode\nobreak\
\bar{\mathrm{O}}_{j}(\mathbf{p},0)\right\rangle,$ (5) $\displaystyle
C^{\mathcal{J}}_{3pt}(\mathbf{p}^{\prime},t;\mathbf{q},\tau)_{ij}$
$\displaystyle=\left\langle\mathrm{O}_{i}(\mathbf{p}^{\prime},t)\leavevmode\nobreak\
\mathcal{J}(\mathbf{q},\tau)\leavevmode\nobreak\
\bar{\mathrm{O}}_{j}(\mathbf{p},0)\right\rangle,$ (6)

where we indicate the three-momentum transfer in the argument of the local
current $\mathcal{J}\in\\{\mathcal{P},\mathcal{A}_{\mu}\\}$. The interpolators
$\mathrm{O}_{i}\in\\{\mathrm{O}_{3q},\mathrm{O}_{5q}\\}$ are projected onto
the $G_{1}$ representation of the double cover of the cubic group
${}^{2}O_{h}$ (or, for nonvanishing momentum, the relevant little group) [31,
32, 33], corresponding to spin and helicity $1/2$ in the continuum, as well as
to definite momentum and isospin. For instance, $I_{3}=-1/2$ corresponds to
$\mathrm{O}_{3q}\sim n$ and
$\mathrm{O}_{5q}\sim\sqrt{1/3}\,n\pi^{0}-\sqrt{2/3}\,p\pi^{-}$. The Wick
contractions of the correlation functions are evaluated using the sequential
method [34] for quark-line connected topologies, while the stochastic “one-
end-trick” [35, 36, 37] is used for disconnected diagrams.

For the results shown here, about 200 propagators (see the Supplemental
Material) for each of the six source-sink separations have been computed on
800 gauge configurations. In view of the computational cost, we carry out the
analysis on a single coordinated lattice simulations (CLS [38]) ensemble
(A653, see Ref. [39]) with the spatial volume $L^{3}=(24a)^{3}$, employing
$N_{f}=3$ nonperturbatively improved Wilson fermions with a lattice spacing
$a\approx 0.098\leavevmode\nobreak\ \mathrm{fm}$ and the pion mass
$m_{\pi}=429\leavevmode\nobreak\ \mathrm{MeV}$. The best results were
obtained, using extended (smeared) quark fields in the nucleon and pion
interpolators. (For details on the smearing, see appendix C.1 and Table 15 of
Ref. [39]).

Figure 1: Top: the effective energies obtained solving the GEVP for
$\mathbf{p}=\mathbf{0}$ and $t_{0}=0.2\leavevmode\nobreak\ \mathrm{fm}$,
compared with the nucleon mass and the lowest energy levels of noninteracting
$N\pi$ P- and $N\pi\pi$ S-waves. Bottom: the moduli of the corresponding
generalized eigenvector components.

We extract the generalized eigenvalue and eigenvector matrices
$\Lambda(\mathbf{p},t;t_{0})=\mathrm{diag}\left(\lambda^{1}(\mathbf{p},t;t_{0}),\lambda^{2}(\mathbf{p},t;t_{0})\right)$
and
$V(\mathbf{p},t,t_{0})=\left(v^{1}(\mathbf{p},t;t_{0}),v^{2}(\mathbf{p},t;t_{0})\right)$,
respectively, by solving the generalized eigenvalue problem (GEVP) [40, 41,
42, 43] for the matrix of two-point functions,
$C_{2pt}(\mathbf{p},t)V(\mathbf{p},t,t_{0})=C_{2pt}(\mathbf{p},t_{0})V(\mathbf{p},t,t_{0})\Lambda(\mathbf{p},t,t_{0})$,
for fixed reference times $t_{0}$, where we employ the normalization
$v^{\alpha\intercal}C_{2pt}(t_{0})v^{\alpha}=1$. For large times $t$ the
eigenvalues will decay exponentially with the energy of the state:
$\lambda^{\alpha}(\mathbf{p},t;t_{0})\rightarrow
d^{\alpha}(\mathbf{p};t_{0})e^{-E^{\alpha}(\mathbf{p})(t-t_{0})}$, where
$d^{\alpha}\lesssim 1$.

The effective energies $E_{\rm
eff}^{\alpha}(t)=a^{-1}\ln[\lambda^{\alpha}(t)/\lambda^{\alpha}(t+a)]$ are
shown in Fig. 1 for $\mathbf{p}=\mathbf{0}$ and $t_{0}=0.2\leavevmode\nobreak\
\mathrm{fm}$. The lowest energy coincides with the nucleon mass on this
ensemble [39], while the second level is close to the sum of the nucleon and
pion energies for the lowest P-wave momentum combination. Therefore, we will
identify $N$ with $\alpha=1$ and $N\pi$ with $\alpha=2$. Note that the
eigenvectors are very stable in $t$ and that the contribution of
$\mathrm{O}_{5q}$ (subscript $i=2$) to the nucleon state is suppressed by more
than $1$ order of magnitude relative to $\mathrm{O}_{3q}$. Nevertheless, as we
will see, the impact on three-point functions can be significant. We also
solve the GEVP for moving frames, in particular for
$\mathbf{p}=\mathbf{e}_{z}\eqqcolon\frac{2\pi}{L}(0,0,1)$
($|\mathbf{e}_{z}|\approx 530\leavevmode\nobreak\ \mathrm{MeV}$). Regarding
$\mathrm{O}_{5q}$, we consider the combinations
$\mathrm{O}_{3q}(\mathbf{e}_{z})\mathrm{O}_{\bar{q}q}(\mathbf{0})$ and
$\mathrm{O}_{3q}(\mathbf{0})\mathrm{O}_{\bar{q}q}(\mathbf{e}_{z})$, with
$\mathrm{O}_{\bar{q}q}$ being a pion interpolator. Solving the GEVP, in both
cases we find the effective energy of the second eigenvalue for
$t>0.5\leavevmode\nobreak\ \mathrm{fm}$ to be consistent with the
$N(\mathbf{e}_{z})\pi(\mathbf{0})$ and $N(\mathbf{0})\pi(\mathbf{e}_{z})$
noninteracting energies, respectively.

Considering these results, we employ the generalized eigenvectors for
$t_{0}=0.2\leavevmode\nobreak\ \mathrm{fm}$ and $t=0.5\leavevmode\nobreak\
\mathrm{fm}$ to construct the GEVP-optimized correlation functions,

$\displaystyle\mathcal{C}_{2pt}(\mathbf{p},t)^{\alpha}$
$\displaystyle=v_{i}^{\alpha}(\mathbf{p})C_{2pt}(\mathbf{p},t)_{ij}v_{j}^{\alpha}(\mathbf{p}),$
(7)
$\displaystyle\mathcal{C}_{3pt}^{\mathcal{J}}(\mathbf{p}^{\prime},t;\mathbf{q},\tau)^{\alpha\beta}$
$\displaystyle=v_{i}^{\alpha}(\mathbf{p}^{\prime})C^{\mathcal{J}}_{3pt}(\mathbf{p}^{\prime},t;\mathbf{q},\tau)_{ij}v_{j}^{\beta}(\mathbf{p}),$
(8)

where $\alpha,\beta\in\\{N,N\pi\\}$. Note that here we only present results
for $\alpha=\beta=N$ and we neglect the $i=j=2$ element ($5q$ to $5q$) of
$C_{3pt}^{\mathcal{J}}$ that, in this case, is suppressed by the second power
of the small eigenvector component $v_{5q}^{N}$. In addition, from ChPT we
would only expect nondiagonal elements of the matrices of correlators to be
enhanced. The nucleon matrix elements of interest are then extracted by
forming the GEVP ratios [44, 45, 46]

$\displaystyle R_{\mathcal{J}}(\mathbf{p}^{\prime},t;\mathbf{q},\tau)$
$\displaystyle=\frac{\mathcal{C}_{3pt}^{\mathcal{J}}(\mathbf{p}^{\prime},t;\mathbf{q},\tau)^{NN}}{\mathcal{C}_{2pt}(\mathbf{p}^{\prime},t)^{N}}$
$\displaystyle\times$
$\displaystyle\sqrt{\frac{\mathcal{C}_{2pt}(\mathbf{p}^{\prime},\tau)^{N}\mathcal{C}_{2pt}(\mathbf{p}^{\prime},t)^{N}\mathcal{C}_{2pt}(\mathbf{p},t-\tau)^{N}}{\mathcal{C}_{2pt}(\mathbf{p},\tau)^{N}\mathcal{C}_{2pt}(\mathbf{p},t)^{N}\mathcal{C}_{2pt}(\mathbf{p}^{\prime},t-\tau)^{N}}}$
$\displaystyle\propto\langle
N_{\mathbf{p}^{\prime}}|\mathcal{J}|N_{\mathbf{p}}\rangle\qquad(t\gg\tau\gg
0),$ (9)

where excited state contributions of the type $N\rightarrow N\pi$ and
$N\pi\rightarrow N$ are explicitly removed. Forming the same ratio for the
usual two- and three-point functions, $C_{2pt,11}$ and
$C_{3pt,11}^{\mathcal{J}}$, will give the same result in the limit of large
$t$ and $\tau$. Any time dependence observed for these ratios (GEVP-improved
or not) is an indication of remaining excited state contamination.

##  IV RESULTS IN THE FORWARD LIMIT

For $\mathbf{p}^{\prime}=\mathbf{p}$ the combination under the square root in
Eq. (9) cancels. We consider two kinematic combinations:
$\mathbf{p}^{\prime}=\mathbf{p}=\mathbf{0}$ and
$\mathbf{p}^{\prime}=\mathbf{p}=\mathbf{e}_{z}$. Regarding the rest frame, the
three-point functions vanish due to parity for $\mathcal{J}=\mathcal{A}_{4}$
and $\mathcal{J}=\mathcal{P}$, while $\mathcal{J}=\mathcal{A}_{i}$ (with the
spin projected in the $i$ direction) at large Euclidean time separations gives
the axial charge $g_{A}$. Contamination from the coupling to $N\pi$ states
exists, however, only as a loop effect in ChPT. Indeed, even when using the
standard ratio, for $t>0.8\leavevmode\nobreak\ \mathrm{fm}$ the data near
$\tau=t/2$ show no time dependence within their errors. Fitting the
(unimproved) ratio for $1.15\leavevmode\nobreak\
\mathrm{fm}<t<1.4\leavevmode\nobreak\ \mathrm{fm}$, we find $g_{A}=1.156(7)$
at our unphysical pion mass.

Figure 2: Comparison of the ratios in Eq. (9) in the forward limit,
constructed from the standard and GEVP-optimized correlation functions, with
$\mathbf{p}=\mathbf{p}^{\prime}=\mathbf{e}_{z}$ as a function of the current
insertion time $\tau$ for a number of different source-sink separations $t$,
where $a\approx 0.098\leavevmode\nobreak\ \mathrm{fm}$. Top: renormalized
ratio for the time component of the axial current. Also shown are the standard
ratios for $\mathcal{A}_{z}$ at $\mathbf{p}=\mathbf{e}_{z}$ and $t=10a$ (black
stars). The dark and light green bands correspond to $g_{A}$, extracted from
this ratio at $\mathbf{p}=\mathbf{0}$ and at $\mathbf{p}=\mathbf{e}_{z}$,
respectively, using different $t$. Bottom: the same for the pseudoscalar
current. The green band highlights the expected result.

Also in the moving frame the $N\pi$ contributions to $\mathcal{A}_{i}$ only
appear as loop effects, and we find that the corresponding standard ratio is
almost constant (black stars in Fig. 2). However, regarding $\mathcal{A}_{4}$
and $\mathcal{P}$, $N$ to $N\pi$ transitions appear at tree-level and are
enhanced by one power of $m_{N}/E_{\pi}$, relative to the $N$ to $N$ matrix
elements of interest. The ratio $R_{\mathcal{A}_{4}}$ will be proportional to
$g_{A}$ at large times too; however, using the $\mathrm{O}_{3q}$
interpolators, we find substantial excited state contamination, which is
indicated by its strong dependence on the source-sink separation, see the blue
symbols in the upper panel of Fig. 2. A difference between the ratios for
$\mathcal{A}_{i}$ and $\mathcal{A}_{4}$ at
$\mathbf{p}^{\prime}=\mathbf{p}\neq\mathbf{0}$ was also observed, e.g., in
Ref. [47], using standard interpolators. In contrast, the GEVP-improved ratios
(red symbols) already agree for $t>0.6\leavevmode\nobreak\ \mathrm{fm}$ with
the value extracted from the standard ratios for $\mathcal{A}_{z}$ obtained at
$\mathbf{p}\neq\mathbf{0}$ (light green band) and at $\mathbf{p}=\mathbf{0}$
(dark green band). A similar, dramatic reduction of the excited state
contamination is observed for the pseudoscalar current for the GEVP-optimized
correlation functions. For this current, due to charge conjugation symmetry,
diagonal matrix elements vanish. Therefore, the deviation from zero of the
standard ratio in the lower panel of Fig. 2 is entirely an excited state
effect, which is removed within the present errors for the GEVP-optimized
ratio. This demonstrates that the dominant contribution is from $N$ to $N\pi$
transitions, which is consistent with the tree-level ChPT expectation.

##  V RESULTS FOR NONVANISHING MOMENTUM TRANSFER

When determining the axial and pseudoscalar form factors, the excited state
contamination is prominent for correlation functions involving the currents
$\mathcal{J}\in\\{\mathcal{A}_{4},\mathcal{P}\\}$ [25, 26, 14, 27], that can
transfer the momentum to a pion at tree-level in ChPT.

Figure 3: The same as in Fig. 2 (but as a function of the insertion time
$\tau$) for $\mathbf{p}^{\prime}=\mathbf{0}\neq\mathbf{p}=-\mathbf{e}_{z}$,
corresponding to $Q^{2}\approx 0.3\leavevmode\nobreak\ \mathrm{GeV}^{2}$. The
green bands represent results of a simultaneous fit to the shown GEVP-
optimized ratios and the standard ratio for $\mathcal{A}_{z}$.

This contribution, which is proportional to $m_{N}/E_{\pi}$ [25, 26], is
largest at small momentum transfer. Therefore, we consider these two currents
and set $\mathbf{q}=\mathbf{e}_{z}$ to a single unit of lattice momentum
($|\mathbf{q}|\approx 530\leavevmode\nobreak\ \mathrm{MeV}$). With
$\mathbf{p}^{\prime}=\mathbf{0}$ and $\mathbf{p}=-\mathbf{q}$, this
corresponds to $Q^{2}\approx 0.3\leavevmode\nobreak\ \mathrm{GeV}^{2}$. Our
results for the two ratios for the standard and optimized correlation
functions are shown in Fig. 3. Clearly, the time dependence is much reduced
for the GEVP-optimized results: at the source, excited states are effectively
removed; however, at the sink (that is at rest) there are clearly residual
effects from higher excitations. The ratios at large times (green bands) are
proportional to the respective matrix elements which, using the decompositions
(1) and (2), are related to a linear combination of $G_{A}$ and
$\widetilde{G}_{P}$ for $\mathcal{A}_{4}$ and $G_{P}$ for $\mathcal{P}$,
respectively.

The axial form factor $G_{A}=0.91\pm 0.01$ is extracted from the standard
correlation functions with $\mathcal{A}_{i}$ and
$\mathbf{e}_{i}\perp\mathbf{q}$. These show ground state dominance within our
range of $t$ and $\tau$. Indeed, the large tree-level $N\pi$ ChPT diagrams do
not contribute to this channel and only $N\pi$ loop diagrams appear [25]. The
ground state matrix element is proportional only to $G_{A}$ and we use the
value we extract as a prior in fits to the other channels. We fit the GEVP-
optimized ratios for the pseudoscalar and the temporal axial currents
simultaneously to constants plus exponentials $\propto e^{-\Delta E(t-\tau)}$.

We find $\Delta E\approx 2m_{\pi}$ for the gap between the nucleon ground
state and this first excitation. The resulting matrix elements then give the
pseudoscalar and induced pseudoscalar form factors at $Q^{2}\approx
0.3\leavevmode\nobreak\ \mathrm{GeV}^{2}$, the latter after subtracting the
$G_{A}$ contribution. The results for $G_{P}$ and $\widetilde{G}_{P}$ as well
as for the PCAC and PPD ratios of Eqs. (3) and (4) are shown in Fig. 4. We
also include results that are obtained using the ChPT guided methods

Figure 4: Results for the form factors and the PCAC and PPD ratios (3) and (4)
at $Q^{2}\approx 0.3\leavevmode\nobreak\ \mathrm{GeV}^{2}$ from the GEVP-
optimized correlators, in comparison to results obtained from the standard
correlation functions, using the ChPT guided multistate analysis techniques of
Ref. [14] ($\mathrm{M}_{1}$) and inspired by [15] ($\mathrm{M}_{2}$).

of Ref. [14] ($\mathrm{M}_{1}$) and a simultaneous fit to the channels
$\mathcal{A}_{4}$, $\mathcal{P}$ and $\mathcal{A}_{i}$ with
$\mathbf{q}=\mathbf{e}_{i}$, inspired by [15] ($\mathrm{M}_{2}$). In spite of
the large excited state contributions, the results using modern multistate
analysis techniques agree within errors with the GEVP results, at least at
$m_{\pi}=429\leavevmode\nobreak\ \mathrm{MeV}$. At this single lattice
spacing, the ratio $r_{\mathrm{PCAC}}$ somewhat differs from $1$.

##  VI CONCLUSIONS

Given the current tension [12, 13] between results for the axial form factor
obtained from lattice QCD and from reanalyses of historical neutrino-deuteron
scattering experiments [19], it is important to rigorously investigate the
systematics associated with the lattice approach. The PCAC relation between
form factors has only recently been verified in some studies [14, 15, 17] and
this provides an important cross check. We have shown that the very large
excited state contributions encountered can be removed by including
$N\pi$-type interpolators. This confirms ChPT expectations, even at our
relatively large pion mass, and supports assumptions made in recent
determinations of the axial form factor [14, 15, 16, 17, 18]. In the near
future, we will repeat the study at a smaller pion mass where excited state
contributions are even larger, with the aim of also determining $N$ to $N\pi$
matrix elements that are relevant for the scattering of neutrinos with
energies larger than $400\leavevmode\nobreak\ \mathrm{MeV}$.

###### Acknowledgements.

The authors thank M. Göckeler, M. Padmanath, S. Prelovsek, P. Wein and T. Wurm
for discussion. L.B. is grateful to G. Herdoiza for engaging discussions
during his secondment to IFT Madrid and to A. Kronfeld, M. Wagman, W. Detmold
and P. Shanahan for fruitful discussions and their hospitality at Fermilab and
MIT. This project has received funding from the European Union’s Horizon 2020
research and innovation program under the Marie Skłodowska Grant Agreement No.
813942 (ITN EuroPLEx) and the EU Grant Agreement No. 824093 (STRONG 2020).
Support from the MISTI Global Seed Fund MIT-Germany - University of Regensburg
project “Quark and gluon structure of the proton” is gratefully acknowledged.
The authors thank the Gauss Centre for Supercomputing (GCS) for providing
computing time through the John von Neumann Institute for Computing (NIC) on
the Booster partition of the supercomputer JURECA [48] at Jülich
Supercomputing Centre (JSC). G. C. S. is the alliance of the three national
supercomputing centres HLRS (Universität Stuttgart), JSC (Forschungszentrum
Jülich), and LRZ (Bayerische Akademie der Wissenschaften), funded by the BMBF
and the German State Ministries for Research of Baden-Württemberg (MWK),
Bayern (StMWFK) and Nordrhein-Westfalen (MIWF). Simulations were also
performed on the QPACE 3 computer of SFB/TR-55, using an adapted version of
the Chroma [49] software package.

Supplemental Material

## Construction of the two- and three-point functions

In order to employ the generalized eigenvalue (and eigenstate) approach, the
following two-point correlation functions (see Eq. (5))

$\displaystyle\langle\mathrm{O}_{3q}(\mathbf{p},t)\leavevmode\nobreak\
\bar{\mathrm{O}}_{3q}(\mathbf{p},0)\rangle,$ (S1)
$\displaystyle\langle\mathrm{O}_{5q}(\mathbf{p},t)\leavevmode\nobreak\
\bar{\mathrm{O}}_{3q}(\mathbf{p},0)\rangle,$ (S2)
$\displaystyle\langle\mathrm{O}_{5q}(\mathbf{p},t)\leavevmode\nobreak\
\bar{\mathrm{O}}_{5q}(\mathbf{p},0)\rangle$ (S3)

and three-point correlation functions (see Eq. (6))

$\displaystyle\langle\mathrm{O}_{3q}(\mathbf{p}^{\prime},t)\leavevmode\nobreak\
\mathcal{J}(\mathbf{q},\tau)\leavevmode\nobreak\
\bar{\mathrm{O}}_{3q}(\mathbf{p},0)\rangle,$ (S4)
$\displaystyle\langle\mathrm{O}_{5q}(\mathbf{p}^{\prime},t)\leavevmode\nobreak\
\mathcal{J}(\mathbf{q},\tau)\leavevmode\nobreak\
\bar{\mathrm{O}}_{3q}(\mathbf{p},0)\rangle$ (S5)

need to be evaluated. In these expressions, $\mathrm{O}_{3q}$ represents a
nucleon-like interpolating operator with a $3$-quark $qqq$ structure and
$\mathrm{O}_{5q}$ is nucleon-pion-like with a $(qqq)(\bar{q}q)$-structure.

We first discuss the construction of the three-point functions, where we
consider transitions from an $I=I_{3}=1/2$ state (e.g., the proton $p$) via a
charged current $\mathcal{J}=\bar{d}\Gamma u$ to an $I=-I_{3}=1/2$ state
(e.g., the neutron, $n$), i.e. $\bar{\mathrm{O}}_{3q}$ has the flavour
structure $\bar{u}\bar{u}\bar{d}\sim\bar{p}$, whereas $\mathrm{O}_{3q}$
corresponds to $udd\sim n$. To project $\mathrm{O}_{5q}$ onto $I=-I_{3}=1/2$,
the combination $\sqrt{1/3}\,n\pi^{0}-\sqrt{2/3}\,p\pi^{-}$ must be formed,
where $\pi^{-}\sim\bar{u}d$ and $\pi^{0}\sim 1/\sqrt{2}(\bar{u}u-\bar{d}{d})$.
Like $\mathrm{O}_{3q}$ also $\mathrm{O}_{5q}$ must be projected onto the
lattice irreducible representation $G_{1}$, see, e.g., Refs. [33, 32]. In the
rest frame, for the spin-up component, we form the combination

$\displaystyle\mathrm{O}_{5q}^{G_{1},\uparrow}(\mathbf{0})=$
$\displaystyle\leavevmode\nobreak\
\mathrm{O}_{3q}(-\mathbf{e}_{x})\mathrm{O}_{\bar{q}q}(\mathbf{e}_{x})-\mathrm{O}_{3q}(\mathbf{e}_{x})\mathrm{O}_{\bar{q}q}(-\mathbf{e}_{x})$
$\displaystyle-i\mathrm{O}_{3q}(-\mathbf{e}_{y})\mathrm{O}_{\bar{q}q}(\mathbf{e}_{y})+i\mathrm{O}_{3q}(\mathbf{e}_{y})\mathrm{O}_{\bar{q}q}(-\mathbf{e}_{y})$
$\displaystyle+\leavevmode\nobreak\
\mathrm{O}_{3q}(-\mathbf{e}_{z})\mathrm{O}_{\bar{q}q}(\mathbf{e}_{z})-\mathrm{O}_{3q}(\mathbf{e}_{z})\mathrm{O}_{\bar{q}q}(-\mathbf{e}_{z}),$
(S6)

where $\mathbf{e}_{i}$ corresponds to one unit of lattice momentum in the
$i$-direction and $\mathrm{O}_{\bar{q}q}$ is a pion interpolator. Regarding a
moving frame with $\mathbf{p}^{\prime}=\mathbf{e}_{z}$, we employ two
combinations that, in the continuum limit, will project on the helicity
$+1/2$:

$\displaystyle\mathrm{O}_{5q}^{G_{1},\uparrow}(\mathbf{e}_{z})^{1}$
$\displaystyle=\leavevmode\nobreak\
\mathrm{O}^{\uparrow}_{3q}(\mathbf{e}_{z})\mathrm{O}_{\bar{q}q}(\mathbf{0}),$
(S7) $\displaystyle\mathrm{O}_{5q}^{G_{1},\uparrow}(\mathbf{e}_{z})^{2}$
$\displaystyle=\leavevmode\nobreak\
\mathrm{O}^{\uparrow}_{3q}(\mathbf{0})\mathrm{O}_{\bar{q}q}(\mathbf{e}_{z}).$
(S8)

Following tree-level ChPT, in the forward limit ($\mathbf{q}=\mathbf{0}$), in
the moving frame ($\mathbf{p}=\mathbf{p}^{\prime}=\mathbf{e}_{z}$), Eq. (S7)
is the relevant interpolator, while for off-forward kinematics
($\mathbf{p}^{\prime}=\mathbf{e}_{z}$, $\mathbf{p}=\mathbf{0}$) Eq. (S8) is
used.

Figure S1: Topologies in position space that appear after performing the Wick
contractions of the three-point correlation functions
$\langle\mathrm{O}_{5q}(x,y)\leavevmode\nobreak\
\mathcal{J}(z)\leavevmode\nobreak\ \bar{\mathrm{O}}_{3q}(0)\rangle$. The three
black lines represent point-to-all propagators, the two red lines correspond
to all-to-all propagators and the red square indicates the current insertion
at $z=(\mathbf{z},\tau)$. The $\mathrm{O}_{5q}(x,y)$ interpolator contains a
$qqq$-structure ($\mathrm{O}_{3q}(x)$) at the spacetime position
$x=(\mathbf{x},t)$ and a $\bar{q}q$-structure ($\mathrm{O}_{\bar{q}q}(y)$) at
$y=(\mathbf{y},t)$.

Performing the Wick contractions, we encounter the four different quark-line
diagram topologies shown in Fig. S1. Each topology represents a number of
different Wick contractions. Note that for the charged current we are
interested in, the

disconnected diagram does not contribute to the transition $p\rightarrow
n\pi^{0}$, whereas all four diagrams contribute to $p\rightarrow p\pi^{-}$. In
each diagram there are two all-to-all propagators (red lines).

The quark-line disconnected diagram on the bottom-right is computed, combining
a point-to-all nucleon two-point function with a stochastically estimated
meson two-point function, using the “one-end-trick” [35, 36, 37]. For the
quark-line connected diagrams, we use the sequential source method [34] to
compute the all-to-all propagators. The Wick contractions for the two-point
functions Eq. (S3) have the same topologies and are computed in a similar way,
replacing the current by a smeared pion interpolating operator at $\tau=0$.

Most computer time is spent on the quark-line connected diagrams, where we
employ the sequential-source method twice for each combination of momentum
projections at $y$ and at $x$111 We use
$\mathcal{D}^{-1}(x,z)\gamma_{5}\mathcal{D}^{-1}(z,x)\gamma_{5}$ (and the same
for $\mathcal{D}^{-1}(x,y)$). Some building blocks are shared between
different diagrams. as well as for each spin polarization of
$\mathrm{O}_{3q}$. The disconnected diagram is much less expensive. The total
cost for each nucleon source position and fixed source-sink separation on a
given configuration is equivalent to the computation of about $200$
propagators. This can be compared to three propagators for the standard three-
point function, when setting $\mathbf{p}^{\prime}=\mathbf{0}$ and considering
the unpolarized case and one polarization. Including
$\mathbf{p}^{\prime}=\pm\mathbf{e}_{i}$ (as we do here) increases this cost
five-fold (15 propagators) and there would be an additional factor if further
polarizations were evaluated. Therefore, the addition of the $\mathrm{O}_{5q}$
interpolators increases the computational complexity by about one order of
magnitude.

## References

  * Fukuda  _et al._ [1998] Y. Fukuda _et al._ (Super-Kamiokande), Evidence for oscillation of atmospheric neutrinos, Phys. Rev. Lett. 81, 1562 (1998), arXiv:hep-ex/9807003 .
  * Ahmad  _et al._ [2001] Q. R. Ahmad _et al._ (SNO), Measurement of the rate of $\nu_{e}+d\to p+p+e^{-}$ interactions produced by 8B solar neutrinos at the Sudbury Neutrino Observatory, Phys. Rev. Lett. 87, 071301 (2001), arXiv:nucl-ex/0106015 .
  * Acero  _et al._ [2022] M. A. Acero _et al._ (NOvA), Improved measurement of neutrino oscillation parameters by the NOvA experiment, Phys. Rev. D 106, 032004 (2022), arXiv:2108.08219 [hep-ex] .
  * Abe  _et al._ [2020] K. Abe _et al._ (T2K), Constraint on the matter–antimatter symmetry-violating phase in neutrino oscillations, Nature 580, 339 (2020), [Erratum: Nature 583, E16 (2020)], arXiv:1910.03887 [hep-ex] .
  * Abi  _et al._ [2020] B. Abi _et al._ (DUNE), Deep Underground Neutrino Experiment (DUNE), far detector technical design report, volume I introduction to DUNE, JINST 15 (08), T08008, arXiv:2002.02967 [physics.ins-det] .
  * Bian  _et al._ [2022] J. Bian _et al._ (Hyper-Kamiokande), Hyper-Kamiokande experiment: a Snowmass white paper, in  _2022 Snowmass Summer Study_ (2022) arXiv:2203.02029 [hep-ex] .
  * Andreopoulos  _et al._ [2010] C. Andreopoulos _et al._ (GENIE), The GENIE neutrino Monte Carlo generator, Nucl. Instrum. Meth. A 614, 87 (2010), arXiv:0905.2517 [hep-ph] .
  * Tena-Vidal  _et al._ [2021] J. Tena-Vidal _et al._ (GENIE), Neutrino-nucleon cross-section model tuning in GENIE v3, Phys. Rev. D 104, 072009 (2021), arXiv:2104.09179 [hep-ph] .
  * Formaggio and Zeller [2012] J. A. Formaggio and G. P. Zeller, From eV to EeV: neutrino cross sections across energy scales, Rev. Mod. Phys. 84, 1307 (2012), arXiv:1305.7513 [hep-ex] .
  * Märkisch  _et al._ [2019] B. Märkisch _et al._ , Measurement of the weak axial-vector coupling constant in the decay of free neutrons using a pulsed cold neutron beam, Phys. Rev. Lett. 122, 242501 (2019), arXiv:1812.04666 [nucl-ex] .
  * Andreev  _et al._ [2013] V. A. Andreev _et al._ (MuCap), Measurement of muon capture on the proton to 1% precision and determination of the pseudoscalar coupling $g_{P}$, Phys. Rev. Lett. 110, 012504 (2013), arXiv:1210.6545 [nucl-ex] .
  * Meyer  _et al._ [2022] A. S. Meyer, A. Walker-Loud, and C. Wilkinson, Status of lattice QCD determination of nucleon form factors and their relevance for the few-GeV neutrino program, Annu. Rev. Nucl. Part. Sci. 72, 205 (2022), arXiv:2201.01839 [hep-lat] .
  * Simons  _et al._ [2022] D. Simons, N. Steinberg, A. Lovato, Y. Meurice, N. Rocco, and M. Wagman, Form factor and model dependence in neutrino-nucleus cross section predictions, (2022), arXiv:2210.02455 [hep-ph] .
  * Bali  _et al._ [2020] G. S. Bali, L. Barca, S. Collins, M. Gruber, M. Löffler, A. Schäfer, W. Söldner, P. Wein, S. Weishäupl, and T. Wurm (RQCD), Nucleon axial structure from lattice QCD, J. High Energy Phys. 05, 126, arXiv:1911.13150 [hep-lat] .
  * Jang  _et al._ [2020] Y.-C. Jang, R. Gupta, B. Yoon, and T. Bhattacharya, axial vector form factors from lattice QCD that satisfy the PCAC relation, Phys. Rev. Lett. 124, 072002 (2020), arXiv:1905.06470 [hep-lat] .
  * Alexandrou  _et al._ [2021] C. Alexandrou _et al._ (ETM), Nucleon axial and pseudoscalar form factors from lattice QCD at the physical point, Phys. Rev. D 103, 034509 (2021), arXiv:2011.13342 [hep-lat] .
  * Park  _et al._ [2022] S. Park, R. Gupta, B. Yoon, S. Mondal, T. Bhattacharya, Y.-C. Jang, B. Joó, and F. Winter (NME), Precision nucleon charges and form factors using (2+1)-flavor lattice QCD, Phys. Rev. D 105, 054505 (2022), arXiv:2103.05599 [hep-lat] .
  * Djukanovic  _et al._ [2022] D. Djukanovic, G. von Hippel, J. Koponen, H. B. Meyer, K. Ottnad, T. Schulz, and H. Wittig, Isovector axial form factor of the nucleon from lattice QCD, Phys. Rev. D 106, 074503 (2022), arXiv:2207.03440 [hep-lat] .
  * Meyer  _et al._ [2016] A. S. Meyer, M. Betancourt, R. Gran, and R. J. Hill, Deuterium target data for precision neutrino-nucleus cross sections, Phys. Rev. D 93, 113015 (2016), arXiv:1603.03048 [hep-ph] .
  * Bali  _et al._ [2015] G. S. Bali, S. Collins, B. Glässle, M. Göckeler, J. Najjar, R. H. Rödl, A. Schäfer, R. W. Schiel, W. Söldner, and A. Sternbeck (RQCD), Nucleon isovector couplings from $N_{f}=2$ lattice QCD, Phys. Rev. D 91, 054501 (2015), arXiv:1412.7336 [hep-lat] .
  * Capitani  _et al._ [2019] S. Capitani, M. Della Morte, D. Djukanovic, G. M. von Hippel, J. Hua, B. Jäger, P. M. Junnarkar, H. B. Meyer, T. D. Rae, and H. Wittig, Isovector axial form factors of the nucleon in two-flavor lattice QCD, Int. J. Mod. Phys. A 34, 1950009 (2019), arXiv:1705.06186 [hep-lat] .
  * Gupta  _et al._ [2017] R. Gupta, Y.-C. Jang, H.-W. Lin, B. Yoon, and T. Bhattacharya (PNDME), Axial vector form factors of the nucleon from lattice QCD, Phys. Rev. D 96, 114503 (2017), arXiv:1705.06834 [hep-lat] .
  * Tsukamoto  _et al._ [2018] N. Tsukamoto, K.-I. Ishikawa, Y. Kuramashi, S. Sasaki, and T. Yamazaki (PACS), Nucleon structure from $2+1$ flavor lattice QCD near the physical point, EPJ Web Conf. 175, 06007 (2018), arXiv:1710.10782 [hep-lat] .
  * Bali  _et al._ [2019] G. S. Bali, S. Collins, M. Gruber, A. Schäfer, P. Wein, and T. Wurm (RQCD), Solving the PCAC puzzle for nucleon axial and pseudoscalar form factors, Phys. Lett. B 789, 666 (2019), arXiv:1810.05569 [hep-lat] .
  * Bär [2019a] O. Bär, $N\pi$-state contamination in lattice calculations of the nucleon axial form factors, Phys. Rev. D 99, 054506 (2019a), arXiv:1812.09191 [hep-lat] .
  * Bär [2019b] O. Bär, $N\pi$-state contamination in lattice calculations of the nucleon pseudoscalar form factor, Phys. Rev. D 100, 054507 (2019b), arXiv:1906.03652 [hep-lat] .
  * Bär [2020] O. Bär, $N\pi$ states and the projection method for the nucleon axial and pseudoscalar form factors, Phys. Rev. D 101, 034515 (2020), arXiv:1912.05873 [hep-lat] .
  * Tiburzi [2015] B. C. Tiburzi, Chiral corrections to nucleon two- and three-point correlation functions, Phys. Rev. D 91, 094510 (2015), arXiv:1503.06329 [hep-lat] .
  * Hansen and Meyer [2017] M. T. Hansen and H. B. Meyer, On the effect of excited states in lattice calculations of the nucleon axial charge, Nucl. Phys. B 923, 558 (2017), arXiv:1610.03843 [hep-lat] .
  * Ruso  _et al._ [2022] L. A. Ruso _et al._ , Theoretical tools for neutrino scattering: interplay between lattice QCD, EFTs, nuclear physics, phenomenology, and neutrino event generators, (2022), arXiv:2203.09030 [hep-ph] .
  * Göckeler  _et al._ [2012] M. Göckeler, R. Horsley, M. Lage, U. Meißner, P. E. L. Rakow, A. Rusetsky, G. Schierholz, and J. M. Zanotti, Scattering phases for meson and baryon resonances on general moving-frame lattices, Phys. Rev. D 86, 094513 (2012), arXiv:1206.4141 [hep-lat] .
  * Lang  _et al._ [2017] C. B. Lang, L. Leskovec, M. Padmanath, and S. Prelovsek, Pion-nucleon scattering in the Roper channel from lattice QCD, Phys. Rev. D 95, 014510 (2017), arXiv:1610.01422 [hep-lat] .
  * Prelovsek  _et al._ [2017] S. Prelovsek, U. Skerbis, and C. B. Lang, Lattice operators for scattering of particles with spin, J. High Energy Phys. 01, 129, arXiv:1607.06738 [hep-lat] .
  * Maiani  _et al._ [1987] L. Maiani, G. Martinelli, M. L. Paciello, and B. Taglienti, Scalar densities and baryon mass differences in lattice QCD with Wilson fermions, Nucl. Phys. B 293, 420 (1987).
  * Sommer [1995] R. Sommer, Leptonic decays of $B$ and $D$ mesons,  _Proceedings, 12th International Symposium on Lattice Field Theory (Lattice 94), Bielefeld, Germany, September 27–October 1, 1994_ , Nucl. Phys. Proc. Suppl. 42, 186 (1995), arXiv:hep-lat/9411024 [hep-lat] .
  * Foster and Michael [1999a] M. Foster and C. Michael (UKQCD), Hadrons with a heavy color adjoint particle, Phys. Rev. D 59, 094509 (1999a), arXiv:hep-lat/9811010 [hep-lat] .
  * Foster and Michael [1999b] M. Foster and C. Michael (UKQCD Collaboration), Quark mass dependence of hadron masses from lattice QCD, Phys. Rev. D 59, 074503 (1999b).
  * Bruno  _et al._ [2015] M. Bruno _et al._ (CLS), Simulation of QCD with N${}_{f}=$ 2 $+$ 1 flavors of non-perturbatively improved Wilson fermions, J. High Energy Phys. 02, 043, arXiv:1411.3982 [hep-lat] .
  * Bali  _et al._ [2022] G. S. Bali, S. Collins, P. Georg, D. Jenkins, P. Korcyl, A. Schäfer, E. E. Scholz, J. Simeth, W. Söldner, and S. Weishäupl (RQCD), Scale setting and the light baryon spectrum in $N_{f}=2+1$ QCD with Wilson fermions, (2022), arXiv:2211.03744 [hep-lat] .
  * Berg [1982] B. Berg, Glueball calculations in lattice gauge theories, J. Phys. Colloq. 43, 272 (1982).
  * Michael [1985] C. Michael, Adjoint sources in lattice gauge theory, Nucl. Phys. B 259, 58 (1985).
  * Lüscher and Wolff [1990] M. Lüscher and U. Wolff, How to calculate the elastic scattering matrix in two-dimensional quantum field theories by numerical simulation, Nucl. Phys. B 339, 222 (1990).
  * Blossier  _et al._ [2009] B. Blossier, M. Della Morte, G. von Hippel, T. Mendes, and R. Sommer, On the generalized eigenvalue method for energies and matrix elements in lattice field theory, J. High Energy Phys. 04, 094, arXiv:0902.1265 [hep-lat] .
  * Bulava  _et al._ [2012] J. Bulava, M. Donnellan, and R. Sommer, On the computation of hadron-to-hadron transition matrix elements in lattice QCD, J. High Energy Phys. 01, 140, arXiv:1108.3774 [hep-lat] .
  * Owen  _et al._ [2013] B. J. Owen, J. Dragos, W. Kamleh, D. B. Leinweber, M. S. Mahbub, B. J. Menadue, and J. M. Zanotti, Variational approach to the calculation of $g_{A}$, Phys. Lett. B 723, 217 (2013), arXiv:1212.4668 [hep-lat] .
  * Dragos  _et al._ [2016] J. Dragos, R. Horsley, W. Kamleh, D. B. Leinweber, Y. Nakamura, P. E. L. Rakow, G. Schierholz, R. D. Young, and J. M. Zanotti, Nucleon matrix elements using the variational method in lattice QCD, Phys. Rev. D 94, 074505 (2016), arXiv:1606.03195 [hep-lat] .
  * Liang  _et al._ [2017] J. Liang, Y.-B. Yang, K.-F. Liu, A. Alexandru, T. Draper, and R. S. Sufian ($\chi$QCD), Lattice calculation of nucleon isovector axial charge with improved currents, Phys. Rev. D 96, 034519 (2017), arXiv:1612.04388 [hep-lat] .
  * Jülich Supercomputing Centre [2018] Jülich Supercomputing Centre, JURECA: modular supercomputer at Jülich Supercomputing Centre, J. of large-scale research facilities 4, A132 (2018).
  * Edwards and Joó [2005] R. G. Edwards and B. Joó (SciDAC, LHP and UKQCD), The Chroma software system for Lattice QCD, Nucl. Phys. B Proc. Suppl. 140, 832 (2005), arXiv:hep-lat/0409003 [hep-lat] .