# Phase discovery with active learning:
Application to structural phase transitions in equiatomic NiTi

Jonathan Vandermause  Department of Physics, Harvard University, Cambridge,
Massachusetts 02138, USA  John A. Paulson School of Engineering and Applied
Sciences, Harvard University, Cambridge, MA 02138, USA  Anders Johansson  John
A. Paulson School of Engineering and Applied Sciences, Harvard University,
Cambridge, MA 02138, USA  Yucong Miao  John A. Paulson School of Engineering
and Applied Sciences, Harvard University, Cambridge, MA 02138, USA  Joost J.
Vlassak  John A. Paulson School of Engineering and Applied Sciences, Harvard
University, Cambridge, MA 02138, USA  Boris Kozinsky<EMAIL_ADDRESS>John A. Paulson School of Engineering and Applied Sciences, Harvard
University, Cambridge, MA 02138, USA  Robert Bosch LLC, Research and
Technology Center, Cambridge, MA 02139, USA

###### Abstract

Nickel titanium (NiTi) is a protypical shape-memory alloy used in a range of
biomedical and engineering devices, but direct molecular dynamics simulations
of the martensitic $B19^{\prime}\rightarrow B2$ phase transition driving its
shape-memory behavior are rare and have relied on classical force fields with
limited accuracy. Here, we train four machine-learned force fields for
equiatomic NiTi based on the LDA, PBE, PBEsol, and SCAN DFT functionals. The
models are trained on the fly during NPT molecular dynamics, with DFT
calculations and model updates performed automatically whenever the
uncertainty of a local energy prediction exceeds a chosen threshold. The
models achieve accuracies of 1-2 meV/atom during training and are shown to
closely track DFT predictions of $B2$ and $B19^{\prime}$ elastic constants and
phonon frequencies. Surprisingly, in large-scale molecular dynamics
simulations, only the SCAN model predicts a reversible
$B19^{\prime}\leftrightarrow B2$ phase transition, with the LDA, PBE, and
PBEsol models predicting a reversible transition to a previously
uncharacterized low-volume phase, which we hypothesize to be a new stable
high-pressure phase. We examine the structure of the new phase and estimate
its stability on the temperature-pressure phase diagram. This work establishes
an automated active learning protocol for studying displacive transformations,
reveals important differences between DFT functionals that can only be
detected in large-scale simulations, provides an accurate force field for
NiTi, and identifies a new phase.

##  I Introduction

The bimetallic alloy nickel titanium (NiTi) exhibits a range of fascinating
behaviors, including superelasticity and the shape-memory effect. Because of
its unique properties, NiTi is used in a range of biomedical, aerospace, and
robotics applications Jani  _et al._ (2014). Its shape-memory behavior is
driven by a structural phase transition between a high-temperature $B2$
austenite phase and a low-temperature $B19^{\prime}$ martensite phase at $\sim
340\text{ K}$ Frenzel  _et al._ (2010). Direct molecular dynamics simulations
of the martensitic phase transition in NiTi could shed light on the
microscopic mechanisms of the shape-memory effect and even guide the rational
design of shape-memory alloys with tunable properties, but an accurate
atomistic description of the stable phases of equiatomic NiTi has posed a
number of theoretical challenges over the last two decades. The principal
computational tool available for studying the microsopic properties of NiTi,
density functional theory (DFT), has revealed a complex picture of its
potential energy surface, with several distinct phases lying close in energy
to the convex hull of the phase diagram Huang  _et al._ (2001, 2003); Holec
_et al._ (2011); Wang and Sehitoglu (2012); Jain  _et al._ (2013); Haskins
_et al._ (2016). For example, DFT predicts the minimum-energy phase at zero
Kelvin to be the orthorombic $B33$ phase (space group $C_{mcm}$) rather than
the $B19^{\prime}$ phase (space group $P2_{1}/m$) reported in experiments,
with the energy gap between the two as small as $8\text{ meV/atom}$ Huang  _et
al._ (2003). A recent ab initio molecular dynamics study demonstrated that
$B19^{\prime}$ has lower free energy than $B33$ above $50\text{ K}$, providing
evidence that $B19^{\prime}$ is stabilized by entropic effects Haskins  _et
al._ (2016). However, there has to date been no direct ab initio simulation of
the martensitic phase transition, due to the extreme computational cost that a
sufficiently large simulation would require.

A more computationally feasible approach to simulating the martensitic phase
transition is to use an empirical interatomic potential or force field to
drive the simulation, as the cost of energy and force prediction with a
classical force field scales only linearly with the size of the system.
However, parameterizing an accurate classical force field that captures the
rich and subtle potential energy surface of NiTi is a formidable challenge. A
number of classical force fields for NiTi have been developed over the years
Farkas  _et al._ (1996); Lai and Liu (2000); Ishida and Hiwatari (2007);
Mutter and Nielaba (2010); Saitoh  _et al._ (2010); Zhong  _et al._ (2011); Ko
_et al._ (2015); Kim  _et al._ (2017); Kavousi  _et al._ (2019), with the 2NN
MEAM model of Ko  _et al._ (2015) specifically targeting the martensitic phase
transition. These classical models rely on strict functional forms, with a
fairly small number of parameters fit to a manually constructed DFT database
or a small set of experimentally measured properties.

In the last decade, machine learning (ML) based approaches have dramatically
increased the flexibility of classical force fields, opening up the
possibility of linear-scaling models that achieve very high accuracy relative
to their ab initio reference data Behler and Parrinello (2007); Handley  _et
al._ (2009); Bartók  _et al._ (2010); Thompson  _et al._ (2015); Wood and
Thompson (2018). A highly accurate ML force field for NiTi would help address
the fundamental challenge of simulating the martensitic phase transition at ab
initio accuracy, offering both quantum accuracy and efficiency of a classical
force field. While many early ML force fields were restricted to relatively
simple single-element systems Szlachta  _et al._ (2014); Bartók  _et al._
(2018); Deringer  _et al._ (2020), recently a neural network based force field
was developed for NiTi based on the PBE DFT functional Tang  _et al._ (2022).
However, a reversible martensitic phase transition was reported with this
force field only in a very small 6x4x4 supercell, with manually selected
phases included in the training set and stress labels excluded.

Here, we train four highly accurate ML force fields for equiatomic NiTi based
on the LDA, PBE, PBEsol, and SCAN DFT functionals. We use our recently
developed active learning methodology Vandermause  _et al._ (2020); Xie  _et
al._ (2021); Vandermause  _et al._ (2021) to autonomously train the models “on
the fly” during molecular dynamics. Structural phase transitions are directly
observed during training, with all models predicting the
$B19^{\prime}\rightarrow B2$ transition on heating at zero pressure. On
cooling, however, the LDA, PBE, and PBEsol models predict transitions from
$B2$ to a previously uncharacterized phase, which we call $M2$. Only the SCAN
model predicts the martensitic $B2\rightarrow B19^{\prime}$ phase transition.
The results are confirmed in much larger cooling/heating simulations of
$12,168$ atoms, which show reversible $B2\rightarrow M2$ predicted by the LDA,
PBE, and PBEsol models and a reversible $B2\rightarrow B19^{\prime}$ predicted
by the SCAN model. We hypothesize that the new $M2$ phase is stable at high
pressure.

##  II Results

###  II.1 On-the-fly training of NiTi force fields

Fig. 1 illustrates our autonomous active learning approach to training NiTi
force fields during molecular dynamics simulations. We train many-body sparse
Gaussian process (SGP) models on four DFT functionals: LDA, PBE, PBEsol, and
SCAN (see Sec. A, Methods for details). To sample different phases of NiTi in
an unbiased fashion, we perform two independent $50$-ps on-the-fly training
simulations of $144$-atom NiTi structures in the NPT ensemble at zero
pressure, i.e. with cell lengths and angles allowed to freely change. The
first is a cooling simulation that begins in the $B2$ austenite phase at
$1000\text{ K}$ and terminates at $50\text{ K}$ (blue lines, Fig. 1(b-f)); the
second is a heating simulation that begins in the $B19^{\prime}$ martensite
phase at low temperature and terminates at $1500\text{ K}$ (red lines, Fig.
1(b-f)). In both the cooling and heating simulations, temperature shocks are
introduced every $5$ ps during the first half of the simulation (vertical
lines, Fig. 1(b-f)) to induce structural phase transitions during training.
Example snapshots from the beginning and end of both simulations are shown in
Fig. 1(a, d) for the PBE and SCAN models, illustrating the changes in cell
shape and size observed during the simulations.

Figure 1: On-the-fly training of NiTi force fields. (a, d) Snapshots from PBE
and SCAN training simulations. Atoms are colored by uncertainty in the local
energy prediction. (b, e) Potential energy predictions during the cooling and
heating simulations. DFT potential energies are shown in black, with vertical
lines indicating changes in temperature. (c, f) Atomic volumes during
training. For clarity, a moving average is plotted in dark red and blue.

Our training protocol is based on our previously introduced active learning
method Vandermause  _et al._ (2020, 2021), with individual MD steps driven by
a Bayesian SGP model that is automatically updated when high-uncertainty local
environments are encountered in the simulation (see Sec. B, Methods for
details). In each training simulation, the model starts from scratch, with no
data initially stored in the training set. DFT calculations (shown in black in
Fig. 1(b,e)) are performed whenever the uncertainty of a local energy
prediction exceeds a fixed prediction threshold. When a DFT calculation is
performed, the training set of the model is updated with the potential energy,
atomic forces and virial stress tensor obtained from the calculation, and
subsequent MD steps are taken with the refined model. In general, the models
display excellent agreement with DFT during training, with typical errors
around $1\text{ meV/atom}$ for energies, $50\text{ meV/\AA}$ for forces, and
$0.1\text{ GPa}$ for stresses in all training simulations (see Figs. S5-S8).

Structural phase transitions are observed in all zero pressure training
simulations, with heating simulations displaying the expected
$B19^{\prime}\rightarrow B2$ transition. This is most clearly seen in plots of
the cell lengths and angles as a function of time (Figs. S1-S4), with sharp
changes in the cell observed between $t=15$ and $t=20\text{ ps}$ as the system
is heated above $1000\text{ K}$. Surprisingly, in the cooling simulations, the
expected $B2\rightarrow B19^{\prime}$ was observed with the SCAN model only,
as indicated by the sharp increase in atomic volume (blue line, Fig. 1(f)) and
the formation of a monoclinic distortion in the unit cell (Fig. S2) observed
at $t\approx 18\text{ ps}$. Structural phase transitions are also observed
with the PBE, LDA, and PBEsol models during cooling, but to lower-volume
structures that are inconsistent with the expected $B19^{\prime}$ martensite
phase.

Visualizations of the structures obtained in each of the four cooling
simulations are presented in Fig. S9. The structures are drawn from the final
picosecond of the simulations, and subsequently relaxed with the final trained
models to eliminate residual forces and stresses. The relaxed SCAN structure
has space group $P2_{1}/m$ and exhibits a monoclinic distortion and elongated
$\vec{c}$ lattice vector, consistent with $B19^{\prime}$ martensite (Fig.
S9(b)). The PBE, LDA, and PBEsol relaxed structures, on the other hand, have a
rhombic prism unit cell and exhibit the orthorombic space group $C_{mcm}$.
Curiously, this is the same space group as the ground state $B33$ phase of
NiTi first characterized in Huang  _et al._ (2003), but the 144-atom
orthorombic structures we observe are not consistent with $B33$ as they cannot
be reduced to a simple 4-atom primitive cell and have significantly smaller
atomic volume. To our best knowledge, this low-volume phase has not been
characterized before, but is the stable phase predicted by our LDA, PBE, and
PBEsol models upon cooling of $B2$ NiTi below the transition temperature.

Figure 2: Cohesive energy versus volume of $B2$ and $B19^{\prime}$ NiTi
predicted by the four models trained in this work (LDA, black; PBEsol, green;
PBE, blue; and SCAN, magenta) and by the MEAM model of Ref. Ko  _et al._
(2015) (red). DFT cohesive energies are shown as dots for $B2$ and triangles
for $B19^{\prime}$. Experimental estimates of the cohesive energy and atomic
volume of $B2$ NiTi derived from Refs. Kubaschewski  _et al._ (1956); Kittel
_et al._ (1996); Prokoshkin  _et al._ (2004) are shown as a gold star.

###  II.2 Model validation

We validate our models by computing a range of challenging properties and
comparing against DFT and a recent 2NN-MEAM NiTi force field Ko  _et al._
(2015). We first consider predictions of cohesive energies as a function of
atomic volume for the $B2$ and $B19^{\prime}$ phases (Fig. 2). The plotted
volumes range between $86\%$ and $114\%$ of the equilibrium volumes and are
obtained by applying isotropic strains to the fully relaxed $B2$ and
$B19^{\prime}$ structures predicted by each model. This allows us to examine
model performance over a wide range of pressures between approximately $-20$
and $40\text{ GPa}$ that were not encountered in the zero-pressure training
simulations. DFT cohesive energies are plotted for a subset of these strained
structures and show excellent agreement with our models over the entire range
of volumes. For reference, we include experimental estimates of the cohesive
energy and atomic volume of $B2$ NiTi derived from Refs. Kubaschewski  _et
al._ (1956); Kittel  _et al._ (1996); Prokoshkin  _et al._ (2004) (gold star,
Fig. 2). The four DFT functionals overestimate the cohesive energy of $B2$ and
underestimate the atomic volume, with LDA showing the most severe disagreement
with experiment. The MEAM force field is in the best agreement with experiment
of the five force field models, in part because the experimentally measured
cohesive energies of pure Ni and Ti were directly incorporated into the fit.

We next consider the elastic constants of $B2$ and $B19^{\prime}$ NiTi.
Predictions of the bulk modulus and three independent elastic constants of the
$B2$ phase are plotted in Fig. S10, with generally good agreement between DFT
and the machine-learned force fields, despite the absence of any relaxed or
low-temperature $B2$ structures in the training set of the models. The MEAM
model is in better agreement with experiment than DFT for the $C_{11}$ and
$C_{12}$ constants, but overestimates $C_{44}$ by more than $130\%$.

Figure 3: Elastic constants of $B19^{\prime}$ NiTi predicted with the PBE-
trained model (dark blue), with PBE DFT (light blue), and with the MEAM model
of Ref. Ko  _et al._ (2015). Figure 4: $B2$ and $B19^{\prime}$ phonon
frequencies predicted with the PBE-trained model (blue), with PBE DFT
(dotted), and with the MEAM model of Ref. Ko  _et al._ (2015).

Predicting the thirteen independent elastic constants of the monoclinic
$B19^{\prime}$ phase represents an even more stringent test of the machine-
learned force fields, and to our knowledge these constants have not been
reported with any NiTi force field in the literature. The predictions of our
PBE-trained model are plotted in Fig. 3 and compared against PBE DFT,
displaying excellent agreement to within $\sim 10\text{ GPa}$ for all thirteen
elastic constants. Similar agreement is observed with the SCAN, LDA, and
PBEsol models (Figs. S11-S13). Our models achieve higher accuracy on the
$B19^{\prime}$ elastic constants than on the $B2$ constants, likely because
many low-temperature $B19^{\prime}$ structures were sampled during the on-the-
fly heating simulations. We also include in Fig. 3 predictions made with the
MEAM model of Ref. Ko  _et al._ (2015), which was also trained on PBE
reference data. While the MEAM model gives accurate predictions of four of the
thirteen $B19^{\prime}$ elastic constants ($C_{12}$, $C_{44}$, $C_{55}$, and
$C_{45}$), it disagrees severely with DFT on the remaining nine, with errors
of more than $100\text{ GPa}$ for $C_{11},C_{33},C_{66},$ and $C_{16}$.

As an additional validation of the PBE-trained force field, we plot the
vibrational phonon frequencies of the $B2$ and $B19^{\prime}$ phases of NiTi
in Fig. 4. Correct description of the $B2$ phonon frequencies is an
interesting test of model transferability, as the model was trained on finite-
temperature $B2$ structures above $300\text{ K}$ while the phonon frequencies
are derived from forces on structures only slightly perturbed from fully
relaxed $B2$. We see good agreement between the ML-predicted phonons and PBE
DFT (blue and black lines, Fig. 4(a)), with negative phonon frequencies near
the M point in the Brillouin zone correctly reproduced. The MEAM model of Ref.
Ko  _et al._ (2015) predicts qualitatively accurate phonon frequencies
relative to PBE DFT, but is less accurate than the PBE-trained ML model. We
plot the phonon frequencies of $B19^{\prime}$ martensite in Fig. 4b, and
observe excellent agreement between the PBE-trained ML model and DFT
throughout the Brillouin zone. The MEAM $B19^{\prime}$ phonon spectrum is in
severe disagreement with PBE DFT, with frequencies overestimated by more than
$100\%$ for the highest-frequency optical branch (see Fig. S14).

###  II.3 Large-scale MD

Finally, we perform large-scale molecular dynamics simulations with the ML
models to probe the structural phase transition in equiatomic NiTi supercells
of $12,168$ atoms (see Fig. 5(a)). We first check the stability of
$B19^{\prime}$ martensite by performing heating simulations of pristine
$B19^{\prime}$ NiTi from $100$ to $1500\text{ K}$ with a heating rate of
$1\text{ K/ps}$ (light colors, Figs. 5 and S15). All models predict that
$B19^{\prime}$ remains stable until transforming to $B2$, indicated by sharp
drops in the atomic volume. The $B19^{\prime}\rightarrow B2$ phase transition
is observed at $750\text{ K}$ for the LDA-trained model, $930\text{ K}$ for
SCAN, $780\text{ K}$ for PBE, and $750\text{ K}$ for PBEsol.

Figure 5: Large-scale MD simulations of phase transitions in equiatomic NiTi
at zero pressure. (a) Snapshots of $B2$ NiTi from the SCAN cooling/heating
simulation (left) and of $B19^{\prime}$ NiTi from the SCAN heating simulation
(right). (b) Cooling/heating simulation of $B2$ NiTi (dark) and heating
simulation of $B19^{\prime}$ NiTi (light colors). Zero-Kelvin volumes of
$B19^{\prime}$ (triangle) and $M2$ (square) are shown for reference. To avoid
overlap in the plot, the PBEsol simulations are reported in Fig. S15.

Next, we perform cooling/heating simulations with each model to examine
whether the $B19^{\prime}\rightarrow B2$ phase transition is reversible (dark
colors, Fig. 5(b)). For each of the four models, the simulation begins in a
pristine $B2$ supercell and is equilibrated for $10\text{ ps}$ at $1000\text{
K}$. The system is then gradually cooled to $100\text{ K}$ and heated to
$1500\text{ K}$ at a rate of $1\text{ K/ps}$. All models predict reversible
structural phase transitions between $\sim 200$ and $700\text{ K}$. Consistent
with the training simulations, the LDA, SCAN, and PBEsol models predict a
transition to the $M2$ phase, indicated by sharp drops in the atomic volume.
Only the SCAN model predicts a reversible transition to $B19^{\prime}$
martensite, with the $B2\rightarrow B19^{\prime}$ transition observed at
$400\text{ K}$ and the $B19^{\prime}\rightarrow B2$ transition observed at
$680\text{ K}$.

We investigate the possibility that $M2$ is stabilized at elevated pressures
by performing cooling/heating and compression/decompression simulations at a
range of pressures and temperatures with the SCAN-trained force field. The
results of these simulations are summarized in the phase diagram shown in Fig.
6, with a clear transition to a low-volume $M2$ phase observed between 12 and
16 GPa (see Figs. S16 and S17 for representative plots of the system volume
during these pressure and temperature sweeps). Computed XRD spectra establish
that the observed low-volume phase is distinct from $B19^{\prime}$ martensite
(Fig. S18), but the spectra of this phase differ across the three force
fields, indicating that the functionals agree in predicting a low volume phase
but produce different atomic packing. Furthermore, symmetry analysis of frames
extracted from molecular dynamics simulations of the $M2$ phase failed to
reveal a space group larger than $P_{1}$ or $P_{-1}$, suggesting $M2$ is a
low-symmetry, possibly incommensurate, phase with properties that are highly
sensitive to the DFT functional.

Figure 6: Phase diagram predicted by the SCAN-trained force field. Vertical
gray lines indicate hysteresis windows of fixed-pressure temperature sweeps
performed between 0 and 20 GPa, while horizontal lines indicate approximate
phase-transition pressures derived from fixed-temperature pressure sweeps
performed between 100 and 600 K (see Figs. S16-17 for example simulations).
Estimated phase boundaries, shown in black, are derived from linear fits to
the centers of the transition windows.

##  III Discussion

We have developed four force fields for equiatomic NiTi based on the LDA, PBE,
PBEsol, and SCAN DFT functionals. The models were trained autonomously during
molecular dynamics, without any constraints applied to the unit cell as the
144-atom model systems were heated and cooled. During these training
simulations, the models are in remarkable agreement with the DFT reference to
within about $1\text{ meV/atom}$ on potential energy predictions, allowing us
to generate MD trajectories at ab initio accuracy at a fraction of the cost of
ab initio molecular dynamics. Nearly all previous computational studies of
equiatomic NiTi have relied on a single DFT functional, most often PBE, and
explored either the zero-Kelvin potential energy surface or the finite-
temperature free energies of specific phases. Our training simulations
eliminate these restrictions, allowing us to directly observe which phases
form as the system is cooled or heated past the phase-transition temperature.

Surprisingly, on cooling of the $B2$ austenite phase, our LDA, PBE, and PBEsol
models predict a transition to a lower-volume phase that does not match any
known phases of equiatomic NiTi. We have labelled this phase $M2$. Only our
SCAN model predicts the $B_{2}\rightarrow B_{19}^{\prime}$ expected from
experiments. Similar behavior is observed in much larger cooling and heating
$12,168$-atom simulations of the phase transition, with the LDA, PBE, and
PBEsol models predicting a sharp decrease in atomic volume on cooling,
consistent with the zero-Kelvin volume of the $M2$ phase.

We have validated the models on a range of challenging properties, and shown
excellent agreement with DFT in predictions of the elastic constants and
phonon frequencies of $B2$ and $B19^{\prime}$ NiTi, showing dramatic
improvement over a previous MEAM force field. The accuracy of our models
suggests that the departure from experiment observed at low temperature is not
due to inaccuracies of the force field models in reproducing the DFT
reference, but rather due to inaccuracies in DFT itself, at least with the
LDA, PBE, and PBEsol functionals. This has the surprising implication that the
most common DFT functionals used in the past to study NiTi, most often PBE,
are in qualitative disagreement with experiment at zero pressure, predicting
the more densely packed $M2$ phase rather than $B19^{\prime}$ martensite.
This, however, eluded detection because large-scale DFT-accurate simulations
were not possible before the possibility to develop machine learning force
fields.

It is possible that the $M2$ phase predicted by our LDA, PBE, and PBEsol
models is stabilized at higher pressure. The phase diagram of equiatomic NiTi
has not been extensively explored in experiments. A recent DFT study has
proposed candidate phases of NiTi at high pressures Liu  _et al._ (2018), but
the analysis was based on zero-temperature enthalpies and therefore excluded
entropic effects. Because $M2$ is more densely packed than $B19^{\prime}$, it
may become favored as pressure increases. We intend to analyze both the
enthalpies and finite-temperature stability of the $M2$ phase at elevated
pressures, to investigate whether $M2$ may form in experiments in the high-
pressure limit.

##  IV Methods

###  IV.1 Sparse Gaussian process models

The models discussed in this work are sparse Gaussian processes trained with
the FLARE++ code, as discussed in greater detail in Ref. Vandermause  _et al._
(2021). Here, we briefly summarize the key features of the SGP model and
report our choice of hyperparameters. The SGP predicts the local energy
$\varepsilon$ of a local environment $\rho_{i}$ centered at atom $i$ by
summing over a set of representative local environments $S$,

$\varepsilon(\rho_{i})=\sum_{s\in
S}k(\mathbf{d}_{i},\mathbf{d}_{s})\alpha_{s},$ (1)

where $k$ is a kernel function quantifying the similarity of two inputs,
$\mathbf{d}_{i}$ and $\mathbf{d}_{s}$ are descriptors of the local
environments $\rho_{i}$ and $\rho_{s}$, and $\bm{\alpha}$ is a vector of
training coefficients. Atomic forces and virial stresses are computed from
gradients of local energy predictions with respect to the atomic positions and
cell vectors. In this work, we use a normalized dot product kernel with
integer power $\xi=2$,

$k(\mathbf{d}_{1},\mathbf{d}_{2})=\sigma^{2}\left(\frac{\mathbf{d}_{1}\cdot\mathbf{d}_{2}}{d_{1}d_{2}}\right)^{\xi}.$
(2)

The signal variance $\sigma^{2}$ is a continuous hyperparameter that is
optimized together with the energy, force and stress noise parameters
$\sigma_{E},\sigma_{F}$ and $\sigma_{S}$ by maximizing the log marginal
likelihood of the SGP. For the descriptor vectors $\mathbf{d}_{i}$, we use the
$B_{2}$ term of the many-body Atomic Cluster Expansion Drautz (2019) with
$N_{\text{rad}}=8$ Chebyshev radial basis functions, $\ell_{\text{max}}=3$ in
the spherical harmonics expansion, and a quadratic cutoff with radius
$r_{\text{cut}}=5\text{ \AA}$.

Uncertainties in the local energy $\varepsilon(\rho_{i})$ are computed as the
square root of a normalized predictive variance $V_{\epsilon}$, defined as

$V_{\varepsilon}=1-\frac{1}{\sigma^{2}}\mathbf{k}_{\varepsilon
S}\mathbf{K}_{SS}^{-1}\mathbf{k}_{\varepsilon S}^{\intercal},$ (3)

where $\mathbf{k}_{\varepsilon S}$ is the vector of kernels between $\rho_{i}$
and the environments $\rho_{s}$ in the sparse set, and $\mathbf{K}_{SS}$ is
the matrix of kernels between pairs of sparse environments.

Final models for each DFT functional are constructed by pooling together the
structures collected during the heating and cooling training simulations and
training a single SGP on all structures. The optimized hyperparameters
$\sigma$, $\sigma_{E}$, $\sigma_{F}$, and $\sigma_{S}$ from the heating
simulations are used in the final trained models. The final models are trained
on cohesive energies rather than potential energies by correcting raw
potential energy labels with isolated-atom energies of Ni and Ti. To
accelerate predictions, the final models are mapped onto equivalent parametric
models that are quadratic in the descriptor $\mathbf{d}_{i}$ and that have a
prediction cost that is independent of the size of the sparse set Vandermause
_et al._ (2021).

###  IV.2 On-the-fly training

Training simulations are performed in the NPT ensemble using the Atomic
Simulation Environment (ASE) Larsen  _et al._ (2017). We use a timestep of
$5\text{ fs}$ and characteristic thermostat and barostat timescales of $25$
and $75\text{ fs}$, respectively. Each training simulation is initialized with
a $144$-atom supercell of $B2$ or $B19^{\prime}$ NiTi with the experimentally
measured lattice parameters from Refs. Goo and Sinclair (1985); Kudoh  _et
al._ (1985). Atomic positions of these initial structures are randomly
jittered by up to $0.1\text{ \AA}$ to give nonzero forces on the first frame.
The temperature schedule for the training simulations is reported in Table I.

Time $t$ (ps) |  Cooling $T$ (K) |  Heating $T$ (K)
---|---|---
$0$ |  1000 |  250
$5$ |  750 |  500
$10$ |  500 |  750
$15$ |  250 |  1000
$20$ |  100 |  1250
$25$ |  50 |  1500
Table 1: Temperature schedule for the $50$-ps on-the-fly training simulations.

During the training simulations, hyperparameter optimization of
$\sigma,\sigma_{E},\sigma_{F},$ and $\sigma_{S}$ is performed with the
L-BFGS-B algorithm implemented in SciPy Virtanen  _et al._ (2020) with the
energy noise $\sigma_{E}$ constrained to be greater than $1\text{ meV/atom}$.
Hyperparameter optimization is performed after each of the first
$N_{\text{hyp}}=10$ updates to the training set, after which the
hyperparameters are held fixed. The sparse set of the SGP is initialized with
$7$ arbitrarily chosen environments from the first frame. DFT is called
whenever the predictive uncertainty on a local energy prediction exceeds
$1\%$. After the DFT calculation is performed, the potential energy, atomic
forces, and six independent components of the virial stress tensor are added
to the training set of the SGP, and local environments with predictive
uncertainty greater than $0.5\%$ are added to the sparse set. After each
update to the model, a minimum of $10$ MD steps are taken with the SGP before
the next DFT calculation is performed.

###  IV.3 DFT details

DFT calculations are performed with plane-wave basis sets and the projector
augmented-wave (PAW) method implemented in the Vienna Ab Initio Simulation
Package (VASP) Kresse (1995). We use a plane-wave kinetic energy cutoff of
$336.9\text{ eV}$ and a k-point mesh of $\sim 4000$ k points per reciprocal
atom. The Methfessel-Paxton smearing method is used with a broadening value of
$0.2\text{ eV}$. The $3\text{d}^{9}4\text{s}^{1}$ and
$3\text{d}^{3}4\text{s}^{1}$ electrons are included in the valence for Ni and
Ti, respectively. For elastic constant and phonon frequency calculations, an
electronic energy convergence criteria of $10^{-7}\text{ eV}$ is used. Phonon
frequencies are calculated using the frozen phonon method as implemented in
Phonopy Togo and Tanaka (2015). For on-the-fly DFT calculations, the default
convergence criteria of $10^{-4}\text{ eV}$ is used. Except for isolated-atom
calculations for Ni and Ti, the calculations are not spin polarized. For the
isolated-atom calculations, a vacuum size of $8\text{ \AA}$ was used.

###  IV.4 MD details

Large-scale MD simulations are performed with LAMMPS using a custom pair style
available in the FLARE++ repository. Simulations are performed in the NPT
ensemble at zero pressure. We use a timestep of $2\text{ fs}$, a temperature
damping parameter of $0.2\text{ ps}$, and a stress damping parameter of
$2\text{ ps}$ in all simulations. Atomic volumes reported in Fig. 5 are
averaged over $500$ consecutive MD steps.

##  V Acknowledgments

This work was supported by the National Science Foundation under Award No.
DMR-1808162 and the Harvard University Materials Research Science and
Engineering Center Grant No. DMR-2011754.

## References

  * Jani  _et al._ (2014) Jaronie Mohd Jani, Martin Leary, Aleksandar Subic, and Mark A Gibson, “A review of shape memory alloy research, applications and opportunities,” Materials & Design (1980-2015) 56, 1078–1113 (2014).
  * Frenzel  _et al._ (2010) J Frenzel, Easo P George, A Dlouhy, Ch Somsen, MF-X Wagner, and G Eggeler, “Influence of Ni on martensitic phase transformations in NiTi shape memory alloys,” Acta Materialia 58, 3444–3458 (2010).
  * Huang  _et al._ (2001) Xiangyang Huang, Claudia Bungaro, Vitaliy Godlevsky, and Karin M Rabe, “Lattice instabilities of cubic NiTi from first principles,” Physical Review B 65, 014108 (2001).
  * Huang  _et al._ (2003) Xiangyang Huang, Graeme J Ackland, and Karin M Rabe, “Crystal structures and shape-memory behaviour of NiTi,” Nature materials 2, 307–311 (2003).
  * Holec  _et al._ (2011) David Holec, Martin Friák, Antonín Dlouhỳ, and Jörg Neugebauer, “Ab initio study of pressure stabilized NiTi allotropes: Pressure-induced transformations and hysteresis loops,” Physical Review B 84, 224119 (2011).
  * Wang and Sehitoglu (2012) Jifeng Wang and Huseyin Sehitoglu, “Resolving quandaries surrounding NiTi,” Applied Physics Letters 101, 081907 (2012).
  * Jain  _et al._ (2013) Anubhav Jain, Shyue Ping Ong, Geoffroy Hautier, Wei Chen, William Davidson Richards, Stephen Dacek, Shreyas Cholia, Dan Gunter, David Skinner, Gerbrand Ceder, _et al._ , “Commentary: The materials project: A materials genome approach to accelerating materials innovation,” APL materials 1, 011002 (2013).
  * Haskins  _et al._ (2016) Justin B Haskins, Alexander E Thompson, and John W Lawson, “Ab initio simulations of phase stability and martensitic transitions in NiTi,” Physical Review B 94, 214110 (2016).
  * Farkas  _et al._ (1996) D Farkas, D Roqueta, A Vilette, and K Ternes, “Atomistic simulations in ternary Ni-Ti-al alloys,” Modelling and Simulation in Materials Science and Engineering 4, 359 (1996).
  * Lai and Liu (2000) WS Lai and BX Liu, “Lattice stability of some Ni-Ti alloy phases versus their chemical composition and disordering,” Journal of Physics: Condensed Matter 12, L53 (2000).
  * Ishida and Hiwatari (2007) H Ishida and Y Hiwatari, “Md simulation of martensitic transformations in TiNi alloys with meam,” Molecular Simulation 33, 459–461 (2007).
  * Mutter and Nielaba (2010) Daniel Mutter and Peter Nielaba, “Simulation of structural phase transitions in NiTi,” Physical Review B 82, 224201 (2010).
  * Saitoh  _et al._ (2010) Ken-ichi Saitoh, Keisuke Kubota, and Tomohiro Sato, “Atomic-level structural change in Ni-Ti alloys under martensite and amorphous transformations,” Technische Mechanik-European Journal of Engineering Mechanics 30, 269–279 (2010).
  * Zhong  _et al._ (2011) Yuan Zhong, Ken Gall, and Ting Zhu, “Atomistic study of nanotwins in NiTi shape memory alloys,” Journal of Applied Physics 110, 033532 (2011).
  * Ko  _et al._ (2015) Won-Seok Ko, Blazej Grabowski, and Jörg Neugebauer, “Development and application of a Ni-Ti interatomic potential with high predictive accuracy of the martensitic phase transition,” Physical Review B 92, 134107 (2015).
  * Kim  _et al._ (2017) Young-Kwang Kim, Hong-Kyu Kim, Woo-Sang Jung, and Byeong-Joo Lee, “Development and application of Ni-Ti and Ni-al-Ti 2nn-meam interatomic potentials for Ni-base superalloys,” Computational Materials Science 139, 225–233 (2017).
  * Kavousi  _et al._ (2019) Sepideh Kavousi, Brian R Novak, Michael I Baskes, Mohsen Asle Zaeem, and Dorel Moldovan, “Modified embedded-atom method potential for high-temperature crystal-melt properties of Ti–Ni alloys and its application to phase field simulation of solidification,” Modelling and Simulation in Materials Science and Engineering 28, 015006 (2019).
  * Behler and Parrinello (2007) Jörg Behler and Michele Parrinello, “Generalized neural-network representation of high-dimensional potential-energy surfaces,” Physical review letters 98, 146401 (2007).
  * Handley  _et al._ (2009) Chris M Handley, Glenn I Hawe, Douglas B Kell, and Paul LA Popelier, “Optimal construction of a fast and accurate polarisable water potential based on multipole moments trained by machine learning,” Physical Chemistry Chemical Physics 11, 6365–6376 (2009).
  * Bartók  _et al._ (2010) Albert P Bartók, Mike C Payne, Risi Kondor, and Gábor Csányi, “Gaussian approximation potentials: The accuracy of quantum mechanics, without the electrons,” Physical review letters 104, 136403 (2010).
  * Thompson  _et al._ (2015) Aidan P Thompson, Laura P Swiler, Christian R Trott, Stephen M Foiles, and Garritt J Tucker, “Spectral neighbor analysis method for automated generation of quantum-accurate interatomic potentials,” Journal of Computational Physics 285, 316–330 (2015).
  * Wood and Thompson (2018) Mitchell A Wood and Aidan P Thompson, “Extending the accuracy of the snap interatomic potential form,” The Journal of chemical physics 148, 241721 (2018).
  * Szlachta  _et al._ (2014) Wojciech J Szlachta, Albert P Bartók, and Gábor Csányi, “Accuracy and transferability of gaussian approximation potential models for tungsten,” Physical Review B 90, 104108 (2014).
  * Bartók  _et al._ (2018) Albert P Bartók, James Kermode, Noam Bernstein, and Gábor Csányi, “Machine learning a general-purpose interatomic potential for silicon,” Physical Review X 8, 041048 (2018).
  * Deringer  _et al._ (2020) Volker L Deringer, Miguel A Caro, and Gábor Csányi, “A general-purpose machine-learning force field for bulk and nanostructured phosphorus,” Nature communications 11, 1–11 (2020).
  * Tang  _et al._ (2022) Hao Tang, Yin Zhang, Qing-Jie Li, Haowei Xu, Yuchi Wang, Yunzhi Wang, and Ju Li, “High accuracy neural network interatomic potential for NiTi shape memory alloy,” Acta Materialia 238, 118217 (2022).
  * Vandermause  _et al._ (2020) Jonathan Vandermause, Steven B Torrisi, Simon Batzner, Yu Xie, Lixin Sun, Alexie M Kolpak, and Boris Kozinsky, “On-the-fly active learning of interpretable bayesian force fields for atomistic rare events,” npj Computational Materials 6, 1–11 (2020).
  * Xie  _et al._ (2021) Yu Xie, Jonathan Vandermause, Lixin Sun, Andrea Cepellotti, and Boris Kozinsky, “Bayesian force fields from active learning for simulation of inter-dimensional transformation of stanene,” npj Computational Materials 7, 1–10 (2021).
  * Vandermause  _et al._ (2021) Jonathan Vandermause, Yu Xie, Jin Soo Lim, Cameron J Owen, and Boris Kozinsky, “Active learning of reactive bayesian force fields: Application to heterogeneous hydrogen-platinum catalysis dynamics,” arXiv preprint arXiv:2106.01949 (2021).
  * Kubaschewski  _et al._ (1956) O Kubaschewski, H Villa, and WA Dench, “The reaction of titanium tetrachloride with hydrogen in contact with various refractories,” Transactions of the Faraday Society 52, 214–222 (1956).
  * Kittel  _et al._ (1996) Charles Kittel, Paul McEuen, and Paul McEuen,  _Introduction to solid state physics_ , Vol. 8 (Wiley New York, 1996).
  * Prokoshkin  _et al._ (2004) SD Prokoshkin, AV Korotitskiy, Vladimir Brailovski, S Turenne, I Yu Khmelevskaya, and IB Trubitsyna, “On the lattice parameters of phases in binary Ti–Ni shape memory alloys,” Acta Materialia 52, 4479–4492 (2004).
  * Liu  _et al._ (2018) Guangtao Liu, Hanyu Liu, Xiaolei Feng, and Simon AT Redfern, “High-pressure phase transitions of nitinol NiTi to a semiconductor with an unusual topological structure,” Physical Review B 97, 140104 (2018).
  * Drautz (2019) Ralf Drautz, “Atomic cluster expansion for accurate and transferable interatomic potentials,” Physical Review B 99, 014104 (2019).
  * Larsen  _et al._ (2017) Ask Hjorth Larsen, Jens Jørgen Mortensen, Jakob Blomqvist, Ivano E Castelli, Rune Christensen, Marcin Dułak, Jesper Friis, Michael N Groves, Bjørk Hammer, Cory Hargus, _et al._ , “The atomic simulation environment — a Python library for working with atoms,” Journal of Physics: Condensed Matter 29, 273002 (2017).
  * Goo and Sinclair (1985) E Goo and R Sinclair, “The b2 to r transformation in Ti50Ni47fe3 and Ti49. 5Ni50. 5 alloys,” Acta Metallurgica 33, 1717–1723 (1985).
  * Kudoh  _et al._ (1985) Yasuhiro Kudoh, Masayasu Tokonami, S Miyazaki, and Kazuhiro Otsuka, “Crystal structure of the martensite in Ti-49.2 at.% Ni alloy analyzed by the single crystal x-ray diffraction method,” Acta Metallurgica 33, 2049–2056 (1985).
  * Virtanen  _et al._ (2020) Pauli Virtanen _et al._ , “SciPy 1.0: Fundamental Algorithms for Scientific Computing in Python,” Nature Methods 17, 261–272 (2020).
  * Kresse (1995) G. Kresse, “Ab-initio molecular-dynamics for liquid-metals,” J. Non-Cryst. Solids 193, 222–229 (1995).
  * Togo and Tanaka (2015) Atsushi Togo and Isao Tanaka, “First principles phonon calculations in materials science,” Scripta Materialia 108, 1–5 (2015).


for the Baryon Scattering (BaSc) Collaboration

# Two-pole nature of the $\Lambda(1405)$ from lattice QCD

John Bulava  Deutsches Elektronen-Synchrotron (DESY), Platanenallee 6, 15738
Zeuthen, Germany  Bárbara Cid-Mora  GSI Helmholtz Centre for Heavy Ion
Research, 64291 Darmstadt, Germany  Andrew D. Hanlon  Physics Department,
Brookhaven National Laboratory, Upton, New York 11973, USA  Ben Hörz  Intel
Deutschland GmbH, Dornacher Str. 1, 85622 Feldkirchen, Germany  Daniel Mohler
Institut für Kernphysik, Technische Universität Darmstadt,
Schlossgartenstrasse 2, 64289 Darmstadt, Germany  GSI Helmholtz Centre for
Heavy Ion Research, 64291 Darmstadt, Germany  Colin Morningstar  Department of
Physics, Carnegie Mellon University, Pittsburgh, Pennsylvania 15213, USA
Joseph Moscoso  Department of Physics and Astronomy, University of North
Carolina, Chapel Hill, NC 27516-3255, USA  Amy Nicholson  Department of
Physics and Astronomy, University of North Carolina, Chapel Hill, NC
27516-3255, USA  Fernando Romero-López  Center for Theoretical Physics,
Massachusetts Institute of Technology, Cambridge, MA 02139, USA  Sarah Skinner
Department of Physics, Carnegie Mellon University, Pittsburgh, Pennsylvania
15213, USA  André Walker-Loud  Nuclear Science Division, Lawrence Berkeley
National Laboratory, Berkeley, CA 94720, USA

###### Abstract

This letter presents the first lattice QCD computation of the coupled channel
$\pi\Sigma$–$\bar{K}N$ scattering amplitudes at energies near $1405\,{\rm
MeV}$. These amplitudes contain the resonance $\Lambda(1405)$ with strangeness
$S=-1$ and isospin, spin, and parity quantum numbers $I(J^{P})=0(1/2^{-})$.
However, whether there is a single resonance or two nearby resonance poles in
this region is controversial theoretically and experimentally. Using single-
baryon and meson-baryon operators to extract the finite-volume stationary-
state energies to obtain the scattering amplitudes at slightly unphysical
quark masses corresponding to $m_{\pi}\approx 200$ MeV and $m_{K}\approx 487$
MeV, this study finds the amplitudes exhibit a virtual bound state below the
$\pi\Sigma$ threshold in addition to the established resonance pole just below
the $\bar{K}N$ threshold. Several parametrizations of the two-channel
$K$-matrix are employed to fit the lattice QCD results, all of which support
the two-pole picture suggested by $SU(3)$ chiral symmetry and unitarity.

lattice QCD, scattering amplitudes

††preprint:  MIT-CTP/5579

Introduction.—The strong nuclear force is described by quantum chromodynamics
(QCD) which governs the dynamics and interactions of quarks and gluons. Due to
an important property of QCD known as asymptotic freedom, the use of
perturbation theory is useful for QCD scattering calculations at very high
energies. The binding of quarks and gluons to form hadrons, such as protons
and neutrons, is a low-energy phenomenon of QCD, requiring a nonperturbative
calculational technique. Such techniques are difficult to apply, so that
understanding the hadron spectrum of QCD remains an important outstanding
issue for the Standard Model of particle physics. In particular, resonances
such as the $\Lambda(1405)$ defy the naive quark-model picture of baryons and
mesons. In this letter, a Markov-chain Monte Carlo method using QCD formulated
on a space-time lattice is applied to shed light on the puzzling hadron
resonance structure in the region of the $\Lambda(1405)$.

The history of the $\Lambda(1405)$ began in Refs. [1, 2] which suggested that
the low-energy $K^{-}p$ amplitude measured in bubble chamber experiments
implies a resonance in the $\pi^{-}\Sigma^{+}$ spectrum just below the
$K^{-}p$ threshold. The intervening decades have witnessed considerable
experimental progress in this system [3, 4], but a consensus about whether
there is a single resonance or two nearby resonance poles in this energy
region has not yet been reached. This is evidenced by the most recent Particle
Data Group review [5] which lists an additional $\Lambda(1380)$ resonance pole
with lower confidence. An improved determination of the $K^{-}p$ scattering
length [6] was enabled by measurements of the energy shift and width of kaonic
hydrogen by the SIDDHARTHA collaboration at DA$\Phi$NE [7]. The angular
analysis of the process $\gamma+p\rightarrow K^{+}+\Sigma+\pi$ by the CLAS
collaboration at JLab determined the line shapes [8] and the spin parity
quantum numbers [9] $J^{P}=1/2^{-}$. The CLAS data has been analyzed in Refs.
[10, 11], which suggest the existence of two isoscalar poles. Recent data from
the BGOOD collaboration [12] and a preliminary study by the GlueX
collaboration [13] also support the two-pole scenario. Similarly, inter-hadron
potentials determined by the ALICE collaboration using the ‘femtoscopy’
approach favor the two-pole picture [14]. However, recent data from J-PARC is
successfully described by a single pole [15] and a combined analysis in Ref.
[16] concludes that a single resonance sufficiently describes the experimental
data, without ruling out the two-pole description.

On the theoretical side, the relatively low mass and quantum numbers of the
$\Lambda(1405)$ are difficult to accommodate in constituent quark models [17].
However, some insight is gained from $SU(3)$ chiral effective theory [18, 19,
20]. The leading-order interaction between the octet of Goldstone bosons and
(ground-state) octet baryons [21, 22] predicts an attractive interaction in
both the flavor-$SU(3)$ singlet and octet combinations. After employing a
unitarization procedure, this attraction leads to two poles in the scattering
matrix analytically continued to complex center-of-mass energies [23]. Despite
the agreement of nearly all chiral approaches (which are reviewed in Refs.
[24, 25, 26]) on the two-pole scenario, the position of the lower pole remains
somewhat poorly constrained [5]. Recent theoretical works about the
$\Lambda(1405)$ can be found in Refs. [27, 28, 29, 30, 31, 32, 33, 34].

Lattice QCD is a first-principles method that can be used to unambiguously
determine the nature of the $\Lambda(1405)$ and provide two unique insights.
First, the elastic $\pi\Sigma$ scattering amplitude can be computed directly
below the $\bar{K}N$ threshold. This process is difficult to access
experimentally and lattice results may help identify and constrain a second
lower pole. Second, the motion of the poles in the complex plane upon varying
the $u$, $d$, and $s$ quark masses away from their physical values provides
additional input to future chiral effective theory analyses [35, 36].

The computation of real-time two-to-two scattering amplitudes below three-
hadron thresholds from imaginary-time lattice QCD calculations is well-
developed and relies on the finite-volume spectrum of interacting two-hadron
states [37, 38, 39, 40, 41, 42, 43, 44]. Previous lattice QCD computations of
the $\Lambda(1405)$ have not computed scattering amplitudes and instead aimed
only to isolate the lowest finite-volume energy eigenstate using single-baryon
three-quark interpolating fields [45, 46, 47, 48, 49, 50, 51, 52, 53]. Using
only such operators is known to be insufficient to extract scattering
information, such as scattering amplitudes and pole locations. The $\bar{K}N$
scattering length for $I=0$ has been computed long ago using the quenched
approximation [54], but mixing with the kinematically-open $\pi\Sigma$ channel
was neglected. The $\pi\Sigma$ and $\bar{K}N$ scattering lengths in other
(non-singlet) flavor and isospin combinations not directly relevant for the
$\Lambda(1405)$ have been computed in Refs. [55, 56, 57].

This work computes the isospin $I=0$ and strangeness $S=-1$ coupled-channel
$\pi\Sigma-\bar{K}N$ scattering amplitudes below the $\pi\pi\Lambda$ threshold
from lattice QCD for the first time. A single ensemble of gauge configurations
with dynamical $u$, $d$, and $s$ quarks is employed with pion and kaon masses
of $m_{\pi}\approx 200\,{\rm MeV}$ and $m_{\rm K}\approx 487$ MeV,
respectively, which deviate slightly from their physical values $m_{\pi}^{\rm
phys}\approx 140$ MeV and $m_{\rm K}^{\rm phys}\approx 495$ MeV. The $u$ and
$d$ quark masses are set to be equal and electroweak interactions are
neglected, so isospin is a good quantum number. The main result of this work
is a set of parametrizations of the amplitudes which are constrained by fits
to the finite-volume energy spectrum. These parametrizations can accommodate
zero, one, or two poles, but when fit to the lattice results and analytically
continued to the complex-energy plane, they all confirm the presence of two
poles, the positions of which vary little and are consistent with predictions
from chiral effective theory. Our use of $m_{\pi}>m_{\pi}^{\rm phys}$ moves
the lower pole just below the $\pi\Sigma$-threshold leading to its unambiguous
identification as a virtual bound state. The higher pole near the
$\bar{K}N$-threshold is also clearly present.

This letter provides a summary of the computation while technical details are
left to a companion paper [58]. The main result is Fig. 1, which shows fits to
the finite-volume spectrum using all parametrizations of the coupled-channel
amplitude and the associated pole positions. Statistical errors are shown for
the parametrization with the lowest Akaike information criterion (AIC) value.

Figure 1: The $I=0$ and $S=-1$ coupled-channel $\pi\Sigma-\bar{K}N$ amplitude
computed on a single lattice QCD gauge-field ensemble with $m_{\pi}\approx
200\,{\rm MeV}$ as a function of the energy difference to the $\pi\Sigma$
threshold in the center-of-mass frame. The upper panel shows the transition
matrix elements, defined in Eq. (4), using the $K$-matrix parametrization with
the lowest AIC constrained by the finite-volume spectrum in the bottom panel.
The second panel shows the model variation for the same quantities using
several parametrizations. The third and fourth panels show the position of
poles in the complex center-of-mass energy ($E_{\rm cm}$) plane on the sheets
closest to the physical one: using the parametrization with lowest AIC (third
panel), and for several parametrizations (fourth panel). In the second and
fourth panel, the transparency of each line and corresponding pair of pole
positions is proportional to
$\exp{[-\left(\text{AIC}-\text{AIC}_{\text{min}}\right)/2]}$, where
$\text{AIC}_{\text{min}}$ is the lowest AIC corresponding to the fit in Eq.
(3), which is also shown in the top panel. The subscripts $i,j$ index the two
open scattering channels. In the lowest panel, the lattice QCD energy levels
that serve as input to the amplitude analyses are displayed. For clarity,
these energy levels are displaced vertically by the total spatial momentum
$\bm{d}^{2}$ defined below Eq. (1).

Determination of finite-volume energies.— The ensemble of gauge configurations
and algorithm for evaluating correlation functions are briefly reviewed here
and discussed more deeply in the companion paper. The $N_{\rm f}=2+1$ QCD
gauge configurations comprise the ‘D200’ ensemble generated by the Coordinated
Lattice Simulations (CLS) consortium [59] which is detailed in Table 1. The
lattice spacing is determined in Ref. [60] and updated in Ref. [61]. All
correlation matrices are computed using the stochastic-LapH [62]
implementation of Ref. [63]. The flexibility afforded by the source-sink
factorization and subsequent computation of correlators via optimized tensor
contractions [64] is particularly advantageous for large Hermitian correlation
matrices containing single-baryon, $\pi\Sigma$, and $\bar{K}N$ interpolating
operators.

Table 1:  Parameters of the D200 ensemble [59]. The lattice dimensions in space and time ($L$ and $T$), as well as the mass of the pion ($m_{\pi}$) and kaon ($m_{\rm K}$) are given in units of the lattice spacing $a$. In pion mass units, the size of the box is $m_{\pi}L=4.181(16)$.  $a[\textup{fm}]$ | $(L/a)^{3}\times T/a$ | $am_{\pi}$ | $am_{\rm K}$
---|---|---|---
0.0633(4)(6) | $64^{3}\times 128$ | 0.06533(25) | 0.15602(16)

The determination of finite-volume stationary-state energies is also discussed
in detail in the companion paper and summarized here. The interaction shift
$\Delta E_{\rm lab}$ of a lab-frame energy from a nearby non-interacting
energy is extracted from a single-state fit to the ratio of a diagonalized
correlation function over the product of correlators for the individual
constituents of the nearby non-interacting energy. The diagonalization of the
correlation matrices is done by solving a generalized eigenvalue problem
(GEVP) as described in Ref. [63]. We have verified the insensitivity of our
extracted energies to reasonable variations of the GEVP parameters, the use of
different nearby non-interacting levels, and different fit forms. An example
energy determination is shown in Fig. 2 for the ground state of the $G_{\rm
1u}(0)$ irreducible representation (irrep), which predominantly contains the
parity-odd, $s$-wave scattering system. All levels used to constrain the
amplitude are shown in Fig. 3.

Scattering amplitude determination.— In lattice QCD, scattering amplitudes
below three-hadron thresholds are inferred from finite-volume spectra [37, 38,
39, 40, 41, 42, 43, 44] using the relationship [65]

$\displaystyle{\rm det}[\widetilde{K}^{-1}(E_{\rm cm})-B^{\bm{P}}(E_{\rm
cm})]=0\,.$ (1)

The matrix $\widetilde{K}$ is related to the usual scattering $K$-matrix
(normalized such that the single-channel equivalent of $\widetilde{K}^{-1}$ is
the $s$-wave $k\cot\delta_{0}$), and the ‘box matrix’ $B^{\bm{P}}$ for a
particular total momentum $\bm{P}=(2\pi/L)\bm{d}$ (with
$\bm{d}\in\mathbb{Z}^{3}$) encodes the reduction in symmetry due to the finite
toroidal spatial volume. $E_{\rm cm}$ denotes center-of-mass energy. Eq. (1)
ignores terms which are suppressed exponentially with the spatial extent $L$.
For scattering between baryons and pseudoscalar mesons, the $K$-matrix does
not mix different $J^{P}$, but does couple the $\pi\Sigma$ and $\bar{K}N$
channels. By contrast, the box matrix is diagonal in the two scattering
channels, but mixes partial waves. $B^{\bm{P}}$ is, however, block diagonal in
the basis given by irreps of the finite-volume little group of momentum
$\bm{P}$. Demanding a vanishing determinant in one of these infinite-
dimensional blocks provides a relationship between the $K$-matrix and the
finite-volume spectrum in a particular irrep. In practice, partial waves with
orbital angular momentum $\ell>\ell_{\rm max}$ are neglected; here $\ell_{\rm
max}=0$ is chosen in both the $\pi\Sigma$ and $\bar{K}N$ channels. The
systematic error due to this is estimated by considering $\ell_{\rm max}=1$,
and found to be insignificant for the near-threshold energy region relevant
here. Specifically, the effect of including additional waves with $\ell=1$
leads to shifts that are significantly smaller that the statistical
uncertainties in fit results for the $\ell=0$ $K$-matrix. Levels from all
irreps in Table 1 of Ref. [63] to which the $J^{P}=1/2^{-}$ partial wave
contributes are employed, as well as one level each from the $G_{\rm 1g}(0)$,
$F_{1}(3)$, and $F_{2}(3)$ irreps for the $\ell_{\rm max}=1$ check. All
elements of $B^{\bm{P}}$ required for this work are given in Ref. [65].

Figure 2: Example determination of a finite-volume stationary-state energy,
illustrating the sensitivity of the fitted energy to the lower end of the fit
range ($t_{\rm min}$) for the lowest level of the $G_{\mathrm{1u}}(0)$ irrep.
Each set of points corresponds to a different fit form. The two-exponential
and geometric [63] fits are performed to the diagonalized correlation function
only. The single-exponential ratio fits are performed to the same correlator
divided by either the product $\bar{K}(0)N(0)$ or $\pi(0)\Sigma(0)$ of
correlators, and the lab frame energy $aE_{\rm lab}$ is reconstructed from the
interaction shifts. The dark horizontal band and filled symbol denote the
chosen fit.  Figure 3: Finite-volume spectrum in the center-of-mass frame used
as input data to constrain parametrizations of the coupled-channel
$\pi\Sigma-\bar{K}N$ scattering amplitude. Each column corresponds to a
particular irrep $\Lambda(\bm{d}^{2})$ of the little group of total momentum
$\bm{P}^{2}=(2\pi/L)^{2}\bm{d}^{2}$. Only irreps where the $\ell=0$ partial
wave contributes are included. Dashed lines indicate various thresholds, as
labeled. Model energies from the resultant scattering-amplitude fit are given
by blue squares.

For $\ell_{\rm max}=0$, the finite-volume spectrum shown in Fig. 3 constrains
the coupled-channel scattering amplitude via Eq. (1) at center-of-mass
energies near the $\pi\Sigma$ and $\bar{K}N$ thresholds. The effective range
expansion (ERE) is used to parametrize the $K$-matrix

$\frac{E_{\rm
cm}}{m_{\pi}}\widetilde{K}_{ij}=A_{ij}+B_{ij}\Delta_{\pi\Sigma},$ (2)

where $A_{ij}$ and $B_{ij}$ are symmetric and real coefficients with $i$ and
$j$ denoting either of the two scattering channels (channel 0 is $\pi\Sigma$
and channel 1 is $\bar{K}N$). Moreover, $\Delta_{\pi\Sigma}=(E_{\rm
cm}^{2}-(m_{\pi}+m_{\Sigma})^{2})/(m_{\pi}+m_{\Sigma})^{2}$ labels the
distance to the $\pi\Sigma$ threshold. The parameters, which are the elements
of the $A$ and $B$ matrices, are determined from fits to the lattice QCD
results using the spectrum method [66]. Similar fits are performed with
variations of the above parametrization: an ERE for $\widetilde{K}^{-1}$,
removing the factor of $E_{\rm cm}$ in Eq. (2), parametrizations inspired by
the Weinberg-Tomozawa potential [19], or using the Blatt-Biedenharn [67] form.
The effect of fixing some (or all) of the elements of $B$ to zero is also
explored.

The correlated-$\chi^{2}$ of the above fits is defined by comparing the
center-of-mass interaction shifts $\Delta E_{\rm cm}$ obtained from the model
with those determined from the ratio fits with a particular choice of the non-
interacting levels. The fit with the lowest Akaike Information Criterion (AIC)
value is a four-parameter fit to Eq. (2). The result is

$\displaystyle\begin{split}&A_{00}=4.1(1.8),\quad\ \,A_{11}=-10.5(1.1),\\\
&A_{01}=10.3(1.5),\quad B_{01}=-29(18),\end{split}$ (3)

with fixed $B_{00}=B_{11}=0$ and $\chi^{2}=10.52$ for 11 degrees of freedom.
This fit is shown in Fig. 1. All statistical uncertainties and correlations
are taken into account using the bootstrap method with $800$ samples.

Analytic structure of the amplitude.—The various parametrizations discussed
above constrain the energy dependence of the amplitudes near the finite-volume
energies, even if they do not accommodate left-hand (cross-channel) cuts.
Knowledge over this limited range enables the analytic continuation of the
scattering amplitude (denoted $\mathcal{T}$) to complex $E_{\rm cm}$ and the
identification of poles close to the real axis on sheets adjacent to the
physical one.

The $K$-matrix, the $J^{P}=1/2^{-}$ scattering amplitude $\mathcal{T}$, and
the normalized amplitude $t$ shown in Fig. 1 are related by

$t^{-1}=\frac{8\pi E_{\rm
cm}}{m_{\pi}}\,\mathcal{T}^{-1}=\widetilde{K}^{-1}-i\widehat{k},\quad$ (4)

where $m_{\pi}\widehat{k}=\text{diag }(k_{\pi\Sigma},k_{\bar{K}N})$,

$\displaystyle k_{\pi\Sigma}^{2}$ $\displaystyle=\frac{1}{4E^{2}_{\rm
cm}}\lambda(E^{2}_{\rm cm},m_{\pi}^{2},m_{\Sigma}^{2})\,.$

Here, $\lambda(x,y,z)$ is the Källén function [68] and $k_{\bar{K}N}$ is
defined similarly. Analytic continuation of the coupled channel
$\pi\Sigma-\bar{K}N$ amplitude involves four different Riemann sheets, each
labelled by the sign of the imaginary parts of $(k_{\pi\Sigma},k_{\bar{K}N})$,
with $(+,+)$ denoting the physical sheet. Complex poles in the scattering
amplitude correspond to vanishing eigenvalues in the right-hand side of Eq.
(4), and are determined numerically. In the vicinity of a pole, the divergent
part of the amplitude is

$t=\frac{m_{\pi}}{E_{\rm cm}-E_{\rm
pole}}\begin{pmatrix}c_{\pi\Sigma}^{2}&c_{\pi\Sigma}\,c_{\bar{K}N}\\\
c_{\pi\Sigma}\,c_{\bar{K}N}&c_{\bar{K}N}^{2}\end{pmatrix}\,+\,\dots,$ (5)

where the (complex) residues $c_{\pi\Sigma}$ and $c_{\bar{K}N}$ denote the
coupling of the resonance pole to each channel.

Two poles are found on the $(-,+)$ sheet, which is the one closest to physical
scattering in the region between the two thresholds. Their locations are

$\displaystyle\begin{split}E_{1}=&1392(9)(2)(16)\text{ MeV},\\\
E_{2}=&[1455(13)(2)(17)-i11.5(4.4)(4)(0.1)]\text{ MeV},\end{split}$ (6)

and their couplings

$\left|\frac{c^{(1)}_{\pi\Sigma}}{c^{(1)}_{\bar{K}N}}\right|=1.9(4)(6),\qquad\left|\frac{c^{(2)}_{\pi\Sigma}}{c^{(2)}_{\bar{K}N}}\right|=0.53(9)(10).$
(7)

The first uncertainty is statistical, the second accounts for parametrization
dependence, and for the pole positions, the third comes from the uncertainty
in the lattice spacing in Table 1. Two poles are present for all
parametrizations of the $K$-matrix. The pole at $E_{1}$ is likely a virtual
bound state, except in 0.5% of bootstrap samples where it is located on the
physical sheet and thus a bound state, while the one at $E_{2}$ is a
resonance. The first pole has a stronger coupling to the $\pi\Sigma$ channel,
while for the second, the hierarchy is reversed, a pattern also predicted by
chiral unitary models. Further confirmation of the existence of the lower pole
as a virtual bound state comes from a single-channel analysis of the energy
levels near the $\pi\Sigma$ threshold, as shown in Fig. 4.

Figure 4: The elastic $\pi\Sigma$ amplitude near threshold. The points are
obtained from Eq. (1) using a single channel and $\ell_{\rm max}=0$. The
shaded band denotes a fit of the four levels shown to a two-parameter
effective range expansion. A pole on the real axis in the second Riemann sheet
(a virtual bound state) occurs when
$k_{\pi\Sigma}\cot\delta_{\pi\Sigma}-ik_{\pi\Sigma}=0$ below threshold. This
is where the black dashed line intersects the fit.

Conclusion.— This study of $\pi\Sigma-\bar{K}N$ scattering in the
$\Lambda(1405)$ region is the first coupled-channel meson-baryon scattering
amplitude determined from lattice QCD. Hermitian correlation matrices using
both single-baryon and meson-baryon interpolating operators for a variety of
different total momenta and irreducible representations were used. The
analytic continuation of the amplitudes into the complex center-of-mass energy
plane is stabilized by finite-volume energies just below the $\pi\Sigma$ and
$\bar{K}N$ thresholds and clearly exhibits two poles. At a slightly heavier-
than-physical pion mass of $m_{\pi}\approx 200\,{\rm MeV}$, the lower pole is
a virtual bound state below the $\pi\Sigma$ threshold and the higher a
resonance just below the $\bar{K}N$ threshold. Due to our use of
$m_{\pi}>m_{\pi}^{\rm phys}$, the real parts of the pole positions in Eq. (6)
are somewhat larger than those determined at the physical point from
experiment using chiral approaches [5], which lie within the ranges ${\rm
Re}\,E_{1}=1325-1380\,{\rm MeV}$ and ${\rm Re}\,E_{2}=1421-1434\,{\rm MeV}$.
Importantly, this qualitative consistency supports the two-pole picture
predicted by chiral symmetry and unitarity.

Future work with this system includes moving to physical quark masses which
requires the consideration of three particle effects, but this should not
present a major problem in the region relevant for the $\Lambda(1405)$.
Estimating residual finite-volume and lattice spacing effects are also
planned. Studying this system along the quark-mass trajectory toward the
SU(3)-symmetric point will also test the motion of the pole positions
predicted by chiral effective theories. Finally, this work opens the door to
investigations of other baryon resonances, such as the $N(1535)$,
$\Lambda(1670)$, $\Sigma(1620)$, and $\Xi(1620)$.

Acknowledgements: We acknowledge helpful discussions with R.J. Hudspith, A.
Jackura, M. Mai, and M. Lutz. We are also grateful to S. Kuberski for
providing re-weighting factors for the D200 ensemble, as computed according to
Ref. [69]. We thank our colleagues within the CLS consortium for sharing
ensembles. Computations were carried out on Frontera [70] at the Texas
Advanced Computing Center (TACC), and at the National Energy Research
Scientific Computing Center (NERSC), a U.S. Department of Energy (DOE) Office
of Science User Facility located at Lawrence Berkeley National Laboratory,
operated under Contract No. DE-AC02-05CH11231 using NERSC awards NP-
ERCAP0005287, NP-ERCAP0010836 and NP-ERCAP0015497. This work was supported in
part by the U.S. National Science Foundation under awards PHY-1913158 and
PHY-2209167 (C.M., S.S.), the Faculty Early Career Development Program
(CAREER) under award PHY-2047185 (A.N.) and by the Graduate Research
Fellowship Program under Grant No. DGE-2040435 (J.M.), the U.S. Department of
Energy, Office of Science, Office of Nuclear Physics, under grant contract
numbers DE-SC0011090 and DE-SC0021006 (F.R.L.), DE-SC0012704 (A.D.H.), DE-
AC02-05CH11231 (A.W.L.) and within the framework of Scientific Discovery
through Advanced Computing (SciDAC) award “Fundamental Nuclear Physics at the
Exascale and Beyond” (A.D.H.), the Mauricio and Carlota Botton Fellowship
(F.R.L.), and the Heisenberg Programme of the Deutsche Forschungsgemeinschaft
(DFG, German Research Foundation) project number 454605793 (D.M.). NumPy [71],
matplotlib [72], and the CHROMA software suite [73] were used for analysis,
plotting, and correlator evaluation.

## References

  * Dalitz and Tuan [1959] R. H. Dalitz and S. F. Tuan, Possible resonant state in pion-hyperon scattering, Phys. Rev. Lett. 2, 425 (1959).
  * Dalitz and Tuan [1960] R. H. Dalitz and S. F. Tuan, The phenomenological representation of $\overline{K}$-nucleon scattering and reaction amplitudes, Annals Phys. 10, 307 (1960).
  * Hyodo and Jido [2012] T. Hyodo and D. Jido, The nature of the $\Lambda(1405)$ resonance in chiral dynamics, Prog. Part. Nucl. Phys. 67, 55 (2012), arXiv:1104.4474 [nucl-th] .
  * Hyodo and Niiyama [2021] T. Hyodo and M. Niiyama, QCD and the strange baryon spectrum, Prog. Part. Nucl. Phys. 120, 103868 (2021), arXiv:2010.07592 [hep-ph] .
  * Workman  _et al._ [2022] R. L. Workman _et al._ (Particle Data Group), Review of Particle Physics, PTEP 2022, 083C01 (2022).
  * Meissner  _et al._ [2004] U. G. Meissner, U. Raha, and A. Rusetsky, Spectrum and decays of kaonic hydrogen, Eur. Phys. J. C 35, 349 (2004), arXiv:hep-ph/0402261 .
  * Bazzi  _et al._ [2011] M. Bazzi _et al._ (SIDDHARTA), A New Measurement of Kaonic Hydrogen X-rays, Phys. Lett. B 704, 113 (2011), arXiv:1105.3090 [nucl-ex] .
  * Moriya  _et al._ [2013] K. Moriya _et al._ (CLAS), Measurement of the $\Sigma$$\pi$ photoproduction line shapes near the $\Lambda$(1405), Phys. Rev. C 87, 035206 (2013), arXiv:1301.5000 [nucl-ex] .
  * Moriya  _et al._ [2014] K. Moriya _et al._ (CLAS Collaboration), Spin and parity measurement of the $\mathrm{\Lambda}(1405)$ baryon, Phys. Rev. Lett. 112, 082004 (2014).
  * Mai and Meißner [2015] M. Mai and U.-G. Meißner, Constraints on the chiral unitary $\bar{K}N$ amplitude from $\pi\Sigma K^{+}$ photoproduction data, Eur. Phys. J. A 51, 30 (2015), arXiv:1411.7884 [hep-ph] .
  * Roca and Oset [2013] L. Roca and E. Oset, Isospin 0 and 1 resonances from $\pi\Sigma$ photoproduction data, Phys. Rev. C 88, 055206 (2013), arXiv:1307.5752 [nucl-th] .
  * Scheluchin  _et al._ [2022] G. Scheluchin _et al._ (BGOOD), Photoproduction of $K^{+}\Lambda(1405)\rightarrow K^{+}\pi^{0}\Sigma^{0}$ extending to forward angles and low momentum transfer, Phys. Lett. B 833, 137375 (2022), arXiv:2108.12235 [nucl-ex] .
  * Wickramaarachchi  _et al._ [2022] N. Wickramaarachchi, R. A. Schumacher, and G. Kalicy (GlueX), Decay of the $\Lambda$(1405) hyperon to $\Sigma^{0}$$\pi^{0}$ measured at GlueX, EPJ Web Conf. 271, 07005 (2022), arXiv:2209.06230 [nucl-ex] .
  * Acharya  _et al._ [2023] S. Acharya _et al._ (ALICE), Constraining the $\overline{K}\,N$ coupled channel dynamics using femtoscopic correlations at the LHC, Eur. Phys. J. C 83, 340 (2023), arXiv:2205.15176 [nucl-ex] .
  * Aikawa  _et al._ [2023] S. Aikawa _et al._ (J-PARC E31), Pole position of $\Lambda$(1405) measured in d(K-,n)$\pi$$\Sigma$ reactions, Phys. Lett. B 837, 137637 (2023), arXiv:2209.08254 [nucl-ex] .
  * Anisovich  _et al._ [2020] A. V. Anisovich, A. V. Sarantsev, V. A. Nikonov, V. Burkert, R. A. Schumacher, U. Thoma, and E. Klempt, Hyperon III: $K^{-}p-\pi\Sigma$ coupled-channel dynamics in the $\Lambda(1405)$ mass region, Eur. Phys. J. A 56, 139 (2020).
  * Isgur and Karl [1978] N. Isgur and G. Karl, P Wave Baryons in the Quark Model, Phys. Rev. D 18, 4187 (1978).
  * Kaiser  _et al._ [1995] N. Kaiser, P. B. Siegel, and W. Weise, Chiral dynamics and the low-energy kaon - nucleon interaction, Nucl. Phys. A 594, 325 (1995), arXiv:nucl-th/9505043 .
  * Oset and Ramos [1998] E. Oset and A. Ramos, Nonperturbative chiral approach to $s$ wave $\overline{K}N$ interactions, Nucl. Phys. A 635, 99 (1998), arXiv:nucl-th/9711022 .
  * Martinez Torres  _et al._ [2012] A. Martinez Torres, M. Bayar, D. Jido, and E. Oset, Strategy to find the two $\Lambda(1405)$ states from lattice QCD simulations, Phys. Rev. C 86, 055201 (2012), arXiv:1202.4297 [hep-lat] .
  * Weinberg [1966] S. Weinberg, Pion scattering lengths, Phys. Rev. Lett. 17, 616 (1966).
  * Tomozawa [1966] Y. Tomozawa, Axial vector coupling renormalization and the meson baryon scattering lengths, Nuovo Cim. A 46, 707 (1966).
  * Oller and Meissner [2001] J. A. Oller and U. G. Meissner, Chiral dynamics in the presence of bound states: Kaon nucleon interactions revisited, Phys. Lett. B 500, 263 (2001), arXiv:hep-ph/0011146 .
  * Mai [2018] M. Mai, Status of the $\Lambda(1405)$, Few Body Syst. 59, 61 (2018).
  * Mai [2021] M. Mai, Review of the ${\Lambda}$(1405) A curious case of a strangeness resonance, Eur. Phys. J. ST 230, 1593 (2021), arXiv:2010.00056 [nucl-th] .
  * Meißner [2020] U.-G. Meißner, Two-pole structures in QCD: Facts, not fantasy!, Symmetry 12, 981 (2020), arXiv:2005.06909 [hep-ph] .
  * Ezoe and Hosaka [2020] T. Ezoe and A. Hosaka, $\Lambda$(1405) as a $\overline{K}N$ Feshbach resonance in the Skyrme model, Phys. Rev. D 102, 014046 (2020), arXiv:2006.03788 [hep-ph] .
  * Myint  _et al._ [2018] K. S. Myint, Y. Akaishi, M. Hassanvand, and T. Yamazaki, Single-pole Nature of the Detectable Lambda(1405), PTEP 2018, 073D01 (2018), arXiv:1804.08240 [nucl-th] .
  * Miyahara  _et al._ [2018] K. Miyahara, T. Hyodo, and W. Weise, Construction of a local $\bar{K}N-\pi\Sigma-\pi\Lambda$ potential and composition of the $\Lambda(1405)$, Phys. Rev. C 98, 025201 (2018), arXiv:1804.08269 [nucl-th] .
  * Miyahara and Hyodo [2018] K. Miyahara and T. Hyodo, Theoretical study of $\Lambda(1405)$ resonance in $\Xi_{b}^{0}\to D^{0}(\pi\Sigma)$ decay, Phys. Rev. C 98, 025202 (2018), arXiv:1803.05572 [nucl-th] .
  * Azizi  _et al._ [2023] K. Azizi, Y. Sarac, and H. Sundu, Investigation of $\Lambda(1405)$ as a molecular pentaquark state, (2023), arXiv:2306.07393 [hep-ph] .
  * Xie  _et al._ [2023] J.-M. Xie, J.-X. Lu, L.-S. Geng, and B.-S. Zou, Two-pole structures demystified: chiral dynamics at work, (unpublished) (2023), arXiv:2307.11631 [hep-ph] .
  * Lu  _et al._ [2023] J.-X. Lu, L.-S. Geng, M. Doering, and M. Mai, Cross-Channel Constraints on Resonant Antikaon-Nucleon Scattering, Phys. Rev. Lett. 130, 071902 (2023), arXiv:2209.02471 [hep-ph] .
  * Hyodo and Weise [2022] T. Hyodo and W. Weise, Theory of Kaon-Nuclear Systems, in  _Handbook of Nuclear Physics_, edited by I. Tanihata, H. Toki, and T. Kajino (2022) pp. 1–34, arXiv:2202.06181 [nucl-th] .
  * Jido  _et al._ [2003] D. Jido, J. A. Oller, E. Oset, A. Ramos, and U. G. Meissner, Chiral dynamics of the two $\Lambda(1405)$ states, Nucl. Phys. A 725, 181 (2003), arXiv:nucl-th/0303062 .
  * Molina and Döring [2016] R. Molina and M. Döring, Pole structure of the $\Lambda$(1405) in a recent QCD simulation, Phys. Rev. D 94, 056010 (2016), [Addendum: Phys.Rev.D 94, 079901 (2016)], arXiv:1512.05831 [hep-lat] .
  * Lüscher [1991] M. Lüscher, Two particle states on a torus and their relation to the scattering matrix, Nucl. Phys. B354, 531 (1991).
  * Rummukainen and Gottlieb [1995] K. Rummukainen and S. A. Gottlieb, Resonance scattering phase shifts on a nonrest frame lattice, Nucl. Phys. B 450, 397 (1995), arXiv:hep-lat/9503028 .
  * Kim  _et al._ [2005] C. H. Kim, C. T. Sachrajda, and S. R. Sharpe, Finite-volume effects for two-hadron states in moving frames, Nucl. Phys. B727, 218 (2005), arXiv:hep-lat/0507006 [hep-lat] .
  * He  _et al._ [2005] S. He, X. Feng, and C. Liu, Two particle states and the S-matrix elements in multi-channel scattering, JHEP 07, 011, arXiv:hep-lat/0504019 [hep-lat] .
  * Bernard  _et al._ [2011] V. Bernard, M. Lage, U. G. Meissner, and A. Rusetsky, Scalar mesons in a finite volume, JHEP 01, 019, arXiv:1010.6018 [hep-lat] .
  * Göckeler  _et al._ [2012] M. Göckeler, R. Horsley, M. Lage, U. G. 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. D86, 094513 (2012), arXiv:1206.4141 [hep-lat] .
  * Briceno and Davoudi [2013] R. A. Briceno and Z. Davoudi, Moving multichannel systems in a finite volume with application to proton-proton fusion, Phys. Rev. D88, 094507 (2013), arXiv:1204.1110 [hep-lat] .
  * Briceno [2014] R. A. Briceno, Two-particle multichannel systems in a finite volume with arbitrary spin, Phys. Rev. D89, 074507 (2014), arXiv:1401.3312 [hep-lat] .
  * Gubler  _et al._ [2016] P. Gubler, T. T. Takahashi, and M. Oka, Flavor structure of $\Lambda$ baryons from lattice QCD: From strange to charm quarks, Phys. Rev. D 94, 114518 (2016), arXiv:1609.01889 [hep-lat] .
  * Menadue  _et al._ [2012] B. J. Menadue, W. Kamleh, D. B. Leinweber, and M. S. Mahbub, Isolating the $\Lambda(1405)$ in Lattice QCD, Phys. Rev. Lett. 108, 112001 (2012), arXiv:1109.6716 [hep-lat] .
  * Engel  _et al._ [2013a] G. P. Engel, C. B. Lang, and A. Schäfer (BGR (Bern-Graz-Regensburg)), Low-lying $\Lambda$ baryons from the lattice, Phys. Rev. D 87, 034502 (2013a), arXiv:1212.2032 [hep-lat] .
  * Engel  _et al._ [2013b] G. P. Engel, C. B. Lang, D. Mohler, and A. Schäfer (BGR), QCD with Two Light Dynamical Chirally Improved Quarks: Baryons, Phys. Rev. D 87, 074504 (2013b), arXiv:1301.4318 [hep-lat] .
  * Nemoto  _et al._ [2003] Y. Nemoto, N. Nakajima, H. Matsufuru, and H. Suganuma, Negative parity baryons in quenched anisotropic lattice QCD, Phys. Rev. D 68, 094505 (2003), arXiv:hep-lat/0302013 .
  * Burch  _et al._ [2006] T. Burch, C. Gattringer, L. Y. Glozman, C. Hagen, D. Hierl, C. B. Lang, and A. Schafer, Excited hadrons on the lattice: Baryons, Phys. Rev. D 74, 014504 (2006), arXiv:hep-lat/0604019 .
  * Takahashi and Oka [2010] T. T. Takahashi and M. Oka, Low-lying Lambda Baryons with spin 1/2 in Two-flavor Lattice QCD, Phys. Rev. D 81, 034505 (2010), arXiv:0910.0686 [hep-lat] .
  * Meinel and Rendon [2022] S. Meinel and G. Rendon, Charm-baryon semileptonic decays and the strange $\Lambda$* resonances: New insights from lattice QCD, Phys. Rev. D 105, L051505 (2022), arXiv:2107.13084 [hep-ph] .
  * Hall  _et al._ [2015] J. M. M. Hall, W. Kamleh, D. B. Leinweber, B. J. Menadue, B. J. Owen, A. W. Thomas, and R. D. Young, Lattice QCD Evidence that the $\Lambda$(1405) Resonance is an Antikaon-Nucleon Molecule, Phys. Rev. Lett. 114, 132002 (2015), arXiv:1411.3402 [hep-lat] .
  * Fukugita  _et al._ [1995] M. Fukugita, Y. Kuramashi, M. Okawa, H. Mino, and A. Ukawa, Hadron scattering lengths in lattice QCD, Phys. Rev. D52, 3003 (1995), arXiv:hep-lat/9501024 [hep-lat] .
  * Detmold and Nicholson [2016] W. Detmold and A. Nicholson, Low energy scattering phase shifts for meson-baryon systems, Phys. Rev. D93, 114511 (2016), arXiv:1511.02275 [hep-lat] .
  * Torok  _et al._ [2010] A. Torok, S. R. Beane, W. Detmold, T. C. Luu, K. Orginos, A. Parreno, M. J. Savage, and A. Walker-Loud, Meson-Baryon Scattering Lengths from Mixed-Action Lattice QCD, Phys. Rev. D81, 074506 (2010), arXiv:0907.1913 [hep-lat] .
  * Meng  _et al._ [2004] G.-W. Meng, C. Miao, X.-N. Du, and C. Liu, Lattice study on kaon nucleon scattering length in the $I=1$ channel, Int. J. Mod. Phys. A 19, 4401 (2004), arXiv:hep-lat/0309048 .
  * Bulava  _et al._ [2024] J. Bulava, B. Cid-Mora, A. D. Hanlon, B. Hörz, D. Mohler, C. Morningstar, J. Moscoso, A. Nicholson, F. Romero-López, S. Skinner, and A. Walker-Loud (Baryon Scattering (BaSc)), Lattice QCD study of $\pi\Sigma-\bar{K}N$ scattering and the $\Lambda(1405)$ resonance, Phys. Rev. D 109, 014511 (2024), arXiv:2307.13471 [hep-lat] .
  * Bruno  _et al._ [2015] M. Bruno _et al._ , Simulation of QCD with N${}_{f}=$ 2 $+$ 1 flavors of non-perturbatively improved Wilson fermions, JHEP 02, 043, arXiv:1411.3982 [hep-lat] .
  * Bruno  _et al._ [2017] M. Bruno, T. Korzec, and S. Schaefer, Setting the scale for the CLS $2+1$ flavor ensembles, Phys. Rev. D95, 074504 (2017), arXiv:1608.08900 [hep-lat] .
  * Strassberger  _et al._ [2022] B. Strassberger _et al._ , Scale setting for CLS 2+1 simulations, PoS LATTICE2021, 135 (2022), arXiv:2112.06696 [hep-lat] .
  * Morningstar  _et al._ [2011] C. Morningstar, J. Bulava, J. Foley, K. J. Juge, D. Lenkner, M. Peardon, and C. H. Wong, Improved stochastic estimation of quark propagation with Laplacian Heaviside smearing in lattice QCD, Phys. Rev. D83, 114505 (2011), arXiv:1104.3870 [hep-lat] .
  * Bulava  _et al._ [2023] J. Bulava, A. D. Hanlon, B. Hörz, C. Morningstar, A. Nicholson, F. Romero-López, S. Skinner, P. Vranas, and A. Walker-Loud, Elastic nucleon-pion scattering at $m_{\pi}=200$ MeV from lattice QCD, Nucl. Phys. B 987, 116105 (2023), arXiv:2208.03867 [hep-lat] .
  * Hörz and Hanlon [2019] B. Hörz and A. Hanlon, Two- and three-pion finite-volume spectra at maximal isospin from lattice QCD, Phys. Rev. Lett. 123, 142002 (2019), arXiv:1905.04277 [hep-lat] .
  * Morningstar  _et al._ [2017] C. Morningstar, J. Bulava, B. Singha, R. Brett, J. Fallica, A. Hanlon, and B. Hörz, Estimating the two-particle $K$-matrix for multiple partial waves and decay channels from finite-volume energies, Nucl. Phys. B924, 477 (2017), arXiv:1707.05817 [hep-lat] .
  * Guo  _et al._ [2013] P. Guo, J. Dudek, R. Edwards, and A. P. Szczepaniak, Coupled-channel scattering on a torus, Phys. Rev. D 88, 014501 (2013), arXiv:1211.0929 [hep-lat] .
  * Blatt and Biedenharn [1952] J. M. Blatt and L. C. Biedenharn, Neutron-Proton Scattering with Spin-Orbit Coupling. 1. General Expressions, Phys. Rev. 86, 399 (1952).
  * G. Källén [1964] G. Källén,  _Elementary Particle Physics_ (Addison-Wesley, 1964).
  * Kuberski [2023] S. Kuberski, Low-mode deflation for twisted-mass and RHMC reweighting in lattice QCD, (2023), arXiv:2306.02385 [hep-lat] .
  * Stanzione  _et al._ [2020] D. Stanzione, J. West, R. Evans, T. Minyard, O. Ghattas, and D. Panda, Frontera: The evolution of leadership computing at the national science foundation, in  _Proceedings of Practice and Experience in Advanced Research Computing (PEARC ’20)_ (2020).
  * Harris  _et al._ [2020] C. R. Harris, K. J. Millman, S. J. Van Der Walt, R. Gommers, P. Virtanen, D. Cournapeau, E. Wieser, J. Taylor, S. Berg, N. J. Smith, _et al._ , Array programming with numpy, Nature 585, 357 (2020).
  * Hunter [2007] J. D. Hunter, Matplotlib: A 2d graphics environment, Computing in Science & Engineering 9, 90 (2007).
  * Edwards and Joo [2005] R. G. Edwards and B. Joo (SciDAC), The Chroma software system for lattice QCD, Nucl. Phys. Proc. Suppl. 140, 832 (2005).