Can the strong interactions between hadrons be determined using femtoscopy?

Abstract

Report issue for preceding element

In the last decades, femtoscopic measurements from heavy-ion collisions have become a popular tool to investigate the strong interactions between hadrons. The key observables measured in such experiments are the two-hadron momentum correlations, which depend on the production mechanism of hadron pairs and the final-state interactions. Given the complexity of ultra-relativistic collision experiments, the source term describing the production mechanism can only be modeled phenomenologically based on numerous assumptions. The commonly employed approach for analyzing femtoscopic data relies on the Koonin-Pratt formula, which relates the measured correlation functions with the relative wave function of an outgoing hadron pair and a source term that is assumed to be universal. Here, we critically examine this universality assumption and show that for strongly interacting particles such as nucleons, the interpretation of femtoscopic measurements suffers from a potentially large intrinsic uncertainty. We also comment on the ongoing efforts to explore three-body interactions using this experimental technique.

Report issue for preceding element

Introduction.—Femtoscopic measurements are becoming increasingly popular as a tool to study hadronic interactions and constitute an integral part of the ongoing and upcoming experimental programs at high-energy heavy-ion facilities such as RHIC, LHC, GSI/FAIR and J-PARC-HI. Originally developed as an imaging technique to extract spatial and temporal characteristics of particle production mechanisms in ultra-relativistic heavy-ion collisions [1, 2, 3, 4], femtoscopy has nowadays evolved into a method for measuring the strong interactions between hadrons [5, 6]. A traditional way of exploring hadron structure and dynamics by means of scattering experiments aims at the determination of the corresponding on-shell scattering amplitudes, which are unambiguously defined and experimentally measurable quantities. Unfortunately, the short lifetime of hadronic resonances makes it difficult to perform scattering experiments beyond the lightest hadrons. For example, there is a wealth of experimental data on nucleon-nucleon (NN) scattering, while only a handful of/no data are available for hyperon-nucleon/hyperon-hyperon scattering. Femtoscopic measurements from high-energy proton-proton (pp) or heavy-ion collisions do not require preparing a beam of unstable particles and thus provide an attractive alternative to scattering experiments. This, however, comes at the price of dealing with a complicated initial state after hadronization, whose modeling becomes an integral part of the analysis of the experimental data. After making several approximations and assumptions [4], the two-particle correlation function $C({\bf k})$ in the center-of-mass system (cms) is usually written in the form [1, 7, 4]

Report issue for preceding element

$$C({\bf k})=\int d{\bf r}\,S_{12}({\bf r})\,\big{|}\Psi({\bf r},\,{\bf k})\big{|}^{2}\,,$$ (1)

where ${\bf k}$ is the relative momentum, $S_{12}({\bf r})$ is the source function and $\Psi({\bf r},\,{\bf k})$ is the relative wave function of the outgoing two-body state, which coincides with the stationary solution of the scattering problem $\Psi({\bf r},\,{\bf k})=\langle{\bf r}|\Psi^{(+)}_{-\bf k}\rangle$ normalized to yield $\langle{\bf r}|\Psi^{(+)}_{-\bf k}\rangle\to e^{-i{\bf k}\cdot{\bf r}}$ in the absence of the interaction, see e.g. Ref. [8]. This is the celebrated Koonin-Pratt formula, which is widely utilized to analyze femtoscopic measurements. The pathway from experiment to interpretation can then be schematically summarized as follows (see Refs. [5, 9] for details):

Report issue for preceding element

• Measurement of the correlation functions $C({\bf k})$.
  

  Report issue for preceding element

• Modeling of the source function $S_{12}({\bf r})$ which is deemed to be universal. Here, one usually assumes a spherically symmetric Gaussian form, whose size $r_{0}$ is extracted from experimental data on $C({\bf k})$ using some model to describe the strong interaction. In particular, $r_{0}$ is often fitted to pp correlation functions using the Reid Soft-Core [10] or Argonne $v_{18}$ [11] potential models [12]. A more sophisticated modeling of the source is described in Ref. [13].
  

  Report issue for preceding element

• Once the source function $S_{12}({\bf r})$ is fixed, Eq. (1) allows one to probe hadronic interaction models using experimental data on $C({\bf k})$.

  Report issue for preceding element

It goes beyond the scope of this paper to address the validity of various assumptions and approximations used to arrive at the formula (1), which are, in fact, well documented in the literature [14, 15, 4, 8]. Rather, we would like to point out that the way Eq. (1) is used to probe hadronic interactions as outlined above suffers from a fundamental flaw: Combined with the universality assumption for the source function $S_{12}({\bf r})$, it implies the measurability of hadronic wave functions and thus also of the corresponding interaction potentials. Yet, hadronic potentials are well known not to represent observable quantities, see e.g. Ref. [16] for a recent discussion. While interaction potentials can, at least in principle, be determined ab initio using lattice QCD [17], their form depends on the choice of interpolating fields for parametrizing hadrons in terms of valence quarks. This inherent scheme dependence does, of course, not affect observable quantities like S-matrix, as can be shown by virtue of the LSZ reduction formalism [18].

Report issue for preceding element

To make the essence of the problem more explicit while keeping the presentation simple, we restrict ourselves to non-relativistic systems in the framework of quantum mechanics and start with rewriting Eq. (1) in a more general form

Report issue for preceding element

$$C({\bf k})=\langle\Psi_{{-\bf k}}^{(+)}|\hat{S}_{12}|\Psi_{{-\bf k}}^{(+)}\rangle\,.$$ (2)

The Koonin-Pratt formula (1) is the coordinate-space representation of Eq. (2), subject to the restriction that the source term is local, $\langle{\bf r}^{\prime}|\hat{S}_{12}|{\bf r}\rangle=\delta({\bf r}^{\prime}-{\bf r})S_{12}({\bf r})$. Eq. (2) makes it explicit that the correlation functions are invariant under a change of basis in the Hilbert space induced by unitary transformations (UTs) $\hat{U}$ as required for any observable quantity:

Report issue for preceding element

$$C({\bf k}) = \big{(}\langle\Psi_{{-\bf k}}^{(+)}|\hat{U}^{\dagger}\big{)}\big{(}\hat{U}\hat{S}_{12}\hat{U}^{\dagger}\big{)}\big{(}\hat{U}|\Psi_{{-\bf k}}^{(+)}\rangle\big{)}$$ (3)

$$\equiv \langle\Psi_{{-\bf k}}^{\prime(+)}|\hat{S}_{12}^{\prime}|\Psi_{{-\bf k}}^{\prime(+)}\rangle\,.$$

However, it also shows that the source term cannot be assumed to be universal across interaction models, and its form must be off-shell consistent with that of the interaction under consideration. Failure to maintain off-shell consistency between the source and the wave function (or interaction potential) by assuming $\hat{S}_{12}^{\prime}=\hat{S}_{12}$ leads to model dependence of the calculated correlation functions. Given the quest for precision in femtoscopy experiments [6], it is important to quantify model dependence in $C({\bf k})$ caused by the assumed universality of the source.

Report issue for preceding element

Gedankenexperiment.—For demonstration purposes, we focus here on two stable distinguishable spinless particles interacting via a short-range force (but all arguments made here remain valid for identical particles, nonvanishing spin and/or if the Coulomb interaction is included) and perform a Gedankenexperiment by letting Alice and Bob analyze two-body correlations.

Report issue for preceding element

The interaction between particles considered by Alice coincides with the spin-$0$ projection of the neutron-proton chiral effective field theory (EFT) potential at fifth order (N⁴LO⁺) [19] with the cutoff parameter $\Lambda=450$ MeV. That is,

Report issue for preceding element

$$\langle p^{\prime}l|\hat{V}_{\rm Alice}|pl\,\rangle\;=\;\langle{p^{\prime}}l|\hat{V}_{\rm np,\;N^{4}LO^{+}}^{s=0}|p\,l\rangle,$$ (4)

where $p\equiv|{\bf p}|$ and $p^{\prime}\equiv|{\bf p}^{\prime}|$ with ${\bf p}$ and ${\bf p}^{\prime}$ the initial and final momenta in the cms, respectively, while $l$ is the orbital angular momentum. The corresponding phase shifts in the $l\leq 3$ partial waves are shown by black open circles in Fig. 1. Note that we use natural units with $\hbar=c=1$ throughout this paper.

Report issue for preceding element

Figure 1: (Color online). S-, P-, D- and F-wave phase shifts for the interaction $V_{\rm Alice}$ are shown with open circles as functions of the cms momentum. The overlapping blue dotted and dashed lines show the phase shifts obtained from $V_{\rm Bob\mbox{-}I}$ and $V_{\rm Bob\mbox{-}II}$, respectively.

For the source term, Alice adopts the commonly chosen Gaussian form

Report issue for preceding element

$$S_{12}^{\rm Alice}({\bf r})\;=\;\frac{e^{-r^{2}/(4r_{0}^{2})}}{(4\pi r_{0}^{2})^{3/2}}\;\;\Rightarrow\;\langle{\bf p}^{\prime}|\hat{S}_{12}^{\rm Alice}|{\bf p}\rangle=e^{-q^{2}r_{0}^{2}}\,,$$ (5)

where $q\equiv|{\bf q}|$ with ${\bf q}={\bf p}^{\prime}-{\bf p}$ the momentum transfer. The radius $r_{0}$ is set to $r_{0}=1.5$ fm. The S-wave momentum-space matrix elements of the potential and the source term used by Alice and specified in Eqs. (4) and (5) are visualized in the upper row of Fig. 2.

Report issue for preceding element

Figure 2: (Color online). S-wave momentum-space matrix elements of the potential in units of GeV^-2 (left column) and of the dimensionless source term (right column) as chosen by Alice (upper row) and seen by Bob (middle and bottom rows).

The correlation functions are usually calculated in coordinate space. Under the commonly made assumption that particles interact only in the $l=0$ state, Eq. (1) can, for a spherically symmetric source function, be reduced to (see, e.g., Ref. [20]):

Report issue for preceding element

$$C(k)=1+\int_{0}^{\infty}4\pi r^{2}drS_{12}(r)\Big{[}|\Psi(r,\,k)|^{2}-|j_{0}(kr)|^{2}\Big{]},$$ (6)

with $j_{0}(kr)$ the spherical Bessel function and the S-wave scattering wave function $\Psi(r,\,k)$ is normalized according to

Report issue for preceding element

$$\Psi(r,\,k)\stackrel{{\scriptstyle r\to\infty}}{{\longrightarrow}}\frac{1}{2ikr}\Big{(}-e^{-2i\delta_{0}(k)}e^{-ikr}+e^{ikr}\Big{)}.$$ (7)

Here, $\delta_{0}(k)$ is the S-wave phase shift. The wave-function $\Psi(r,\,k)$ can be obtained from the half-shell T-matrix, see e.g. Ref. [21]. However, to avoid numerical integrations of strongly oscillating functions, it is more convenient to directly compute the matrix element in Eq. (2) in the partial-wave momentum-space basis. For a spherically symmetric source, we find (see also Ref. [22]):

Report issue for preceding element

$$C(k) = \sum_{l}\frac{2l+1}{4\pi}\,\cos^{2}\delta_{l}(k)\,\Big{[}S_{kk}^{l}+K_{kp}^{l}\circ S_{pk}^{l}$$ (8)

$$+ S_{kp}^{l}\circ K_{pk}^{l}+K_{kp}^{l}\circ S_{pp^{\prime}}^{l}\circ K_{p^{\prime}k}^{l}\Big{]},$$

where we have introduced a short-hand notation $S_{pp^{\prime}}^{l}\equiv\langle pl|\hat{S}_{12}|p^{\prime}l\rangle$, $K_{pk}^{l}=K_{kp}^{l}\equiv\langle pl|\hat{K}|kl\rangle$ and

Report issue for preceding element

$$A_{pp^{\prime}}^{l}\circ B_{p^{\prime}p^{\prime\prime}}^{l} \equiv \mathchoice{{\vbox{\hbox{$\textstyle-$}}\kern-4.86108pt}}{{\vbox{\hbox{$\scriptstyle-$}}\kern-3.25pt}}{{\vbox{\hbox{$\scriptscriptstyle-$}}\kern-2.29166pt}}{{\vbox{\hbox{$\scriptscriptstyle-$}}\kern-1.875pt}}\!\int_{0}^{\infty}p^{\prime 2}dp^{\prime}\,A_{pp^{\prime}}^{l}\,\frac{2\mu}{k^{2}-p^{\prime 2}}\,B_{p^{\prime}p^{\prime\prime}}^{l}\,.$$ (9)

Here, $\mu$ is the reduced mass and $\mathchoice{{\vbox{\hbox{$\textstyle-$}}\kern-4.86108pt}}{{\vbox{\hbox{$\scriptstyle-$}}\kern-3.25pt}}{{\vbox{\hbox{$\scriptscriptstyle-$}}\kern-2.29166pt}}{{\vbox{\hbox{$\scriptscriptstyle-$}}\kern-1.875pt}}\!\int$ denotes the Cauchy principal value integral. The half-shell K-matrix can be easily calculated for arbitrary short-range interactions by solving the Lippmann-Schwinger integral equation

Report issue for preceding element

$$K_{pk}^{l} = V_{pk}^{l}+V_{pp^{\prime}}^{l}\circ K_{p^{\prime}k}^{l}$$ (10)

using standard numerical methods. The on-shell K-matrix is related to the S-matrix via $S_{l}(k)=e^{2i\delta_{l}(k)}=1-2i\pi\mu kK^{l}_{kk}$. We note in passing that Eq. (8) can be brought to the form of Eq. (6) if the source term is local and the interaction vanishes in all $l\neq 0$ channels. The correlation function calculated by Alice using Eq. (8) is shown by open circles in Fig. 3. We have also checked these results by using Eq. (1) directly.

Report issue for preceding element

Figure 3: (Color online). Correlation functions calculated by Alice and Bob. The results of Alice are shown by open circles. Red dash-dotted and dash-double-dotted lines show the correlation functions calculated by Bob using the interactions $V_{\rm Bob\mbox{-}I}$ and $V_{\rm Bob\mbox{-}II}$, respectively, and the source term $S_{12}^{\rm Alice}$ (assumed to be universal). The overlapping blue dotted and dashed lines show the results of Bob using the properly transformed source terms $S_{12}^{\rm Bob\mbox{-}I}$ and $S_{12}^{\rm Bob\mbox{-}II}$, respectively.

On the other hand, Bob performs his analysis using a different basis in the Hilbert space with $|\Psi_{\rm Bob}\rangle=\hat{U}|\Psi_{\rm Alice}\rangle$, where the unitary operator is taken as a rank-1 separable form

Report issue for preceding element

$$\hat{U} = 1-2|g\rangle\langle g|,\quad\langle g|g\rangle=1\,,$$ (11)

and $|g\rangle$ is chosen following Ref. [23] as

Report issue for preceding element

$$g({\bf r}) \equiv \langle{\bf r}|g\rangle=Cr(1-\beta r)e^{-\alpha r}\,,$$ (12)

with $C$ the normalization constant fixed from $\int d{\bf r}\big{[}g({\bf r})\big{]}^{2}=1$. The parameter $\alpha$ specifies the inverse range of the transformation, while $\beta^{-1}$ gives the position of the node in the profile function $g(r)$. The interaction $V_{\rm Alice}$ in Bob’s convention has the form

Report issue for preceding element

$$\hat{V}_{\rm Bob} = \hat{U}\bigg{(}\frac{\hat{p}^{2}}{2\mu}+\hat{V}_{\rm Alice}\bigg{)}\hat{U}^{\dagger}-\frac{\hat{p}^{2}}{2\mu}\,,$$ (13)

where $\hat{\bf p}$ is the cms momentum operator. It can easily be calculated numerically in momentum space. Notice that the considered transformation acts only in the S-wave. To be specific, consider two different transformations corresponding to the parameters $\alpha=1.0\mbox{ fm}^{-1}$, $\beta=0.25\mbox{ fm}^{-1}$ (set-I) and $\alpha=0.7\mbox{ fm}^{-1}$, $\beta=2.0\mbox{ fm}^{-1}$ (set-II). These transformations significantly affect the S-wave momentum-space matrix elements of the interaction as shown in Fig. 2. Yet, physical observables are, of course, independent of the choice of basis in the Hilbert space. In particular, we have verified that the phase shifts calculated from $\hat{V}_{\rm Alice}$, $\hat{V}_{\rm Bob\mbox{-}I}$ and $\hat{V}_{\rm Bob\mbox{-}II}$ agree with each other to machine precision, see Fig. 1. In contrast, the correlation functions calculated by Bob under the assumption of the universal source term, $\hat{S}^{\rm Bob}_{12}=\hat{S}^{\rm Alice}_{12}$, are as expected different from the result found by Alice, as depicted by the red lines in Fig. 3. To restore the result for $C(k)$ obtained by Alice, the source term needs to be brought into Bob’s convention by applying the UT, $\hat{S}^{\rm Bob}_{12}=\hat{U}\hat{S}^{\rm Alice}_{12}\hat{U}^{\dagger}$, see the right column of Fig. 2 and blue lines in Fig. 3.

Report issue for preceding element

Scheme dependence in chiral EFT.—The above example shows that off-shell ambiguities of the interaction potentials can result in a large model dependence of the correlation function if the source term is assumed to be universal. On the other hand, realistic models of hadronic interactions are constrained by physical principles like, e.g., pion-exchange dominance at large distances and often share similarities when it comes to modeling of the short-distance behavior. The remaining off-shell ambiguities may therefore be expected to be less pronounced in practice. In this context, it is instructive to examine how scheme dependence manifests itself in chiral EFT for nuclear forces, [24, 25, 26], the most extensively studied and best understood hadronic interactions. Leaving aside the (significant) ambiguity from the choice of the regulator, the inherent scheme dependence in nuclear chiral EFT first shows up at fourth expansion order (N³LO), stemming from the long-range relativistic corrections that depend on two arbitrary phases [27, 28] and three short-range interactions contributing to the ¹S₀, ³S₁ and ³S₁-³D₁ potentials [19]. The corresponding five parameters can be chosen arbitrarily, subject to the naturalness constraint, and their values can be changed by applying suitable UTs, see Supplemental Material [29] for explicit expressions. Importantly, such UTs also induce many-body interactions and exchange currents. This shows, once again, that various scheme-dependent quantities such as two-body interactions, three-body forces (3BFs) and exchange currents must be used consistently with each other to avoid an uncontrollable model dependence.

Report issue for preceding element

The above power-counting-based arguments point towards a mild scheme-dependence of NN interactions in chiral EFT. However, the situation is different for 3BFs, which provide small but important corrections to the dominant pairwise forces and remain a challenging frontier in nuclear physics [30], see also Refs. [31, 32, 33] for recent attempts to explore 3BFs through femtoscopy. 3BFs first appear at third order (N²LO) in chiral EFT, and thus are much more sensitive to the above mentioned off-shell ambiguities. To illustrate this point we have generated a set of $27$ NN N⁴LO⁺ potentials corresponding to different choices of the off-shell short-range interactions, see [29] for details. These potentials are nearly phase-equivalent in the two-body sector and equally valid from the EFT point of view. However, they lead to vastly different predictions for three-nucleon observables as visualized in Fig. 4, which illustrates the well-known inherent scheme dependence of [34] 3BFs [35]. In particular, the required 3BF contributions to the ³H binding energy can be both attractive and repulsive depending on the off-shell behavior of the employed NN interaction.

Report issue for preceding element

Figure 4: (Color online). Neutron-deuteron total cross section (left) and the ³H binding energy (right) calculated using the N⁴LO⁺ NN potential of Ref. [36] (orange lines). Light-shaded blue bands show the results from the phase-equivalent but off-shell differing N⁴LO⁺ NN potentials as explained in the text.

Discussion and conclusions.—Hadronic correlations measured in ultra-relativistic collisions are sensitive to the strong interactions. However, probing final state interactions by means of the Koonin-Pratt formula violates basic principles of quantum mechanics if the source model is regarded to be universal. Using an example of two strongly interacting distinguishable particles, we have shown that off-shell ambiguities in the interaction can then translate into a significant model dependence for $C({\bf k})$. It is important to emphasize that the problematic universality assumption is essential as its relaxation sacrifices the predictive power of the femtoscopy approach.

Report issue for preceding element

The sensitivity of $C({\bf k})$ to the off-shell behavior of the strong force decreases for source radii $r_{0}$ being large compared to the interaction range since the large-distance behavior of the wave function $\Psi({\bf r},\,{\bf k})$ is unambiguously determined by the S-matrix [37, 8, 4]. The interpretation of $C({\bf k})$ in terms of the average $|\Psi({\bf r},\,{\bf k})|^{2}$ then becomes no more inherently problematic. Large $r_{0}$-values, however, also reduce the strength of femtoscopic signals. The source size in our example is, in fact, comparable to or even larger than those typically used in the literature [5, 6, 13, 9, 32]. Even smaller values for $r_{0}$ were obtained recently using a precise pion-kaon interaction [38].

Report issue for preceding element

Finally, we also discussed off-shell ambiguities of nuclear interactions in chiral EFT. While scheme-dependent contributions to the NN force are suppressed by power counting, the 3BF is highly sensitive to the off-shell components of the NN potentials. The lack of the off-shell consistency in femtoscopy analyses, therefore, raises concerns about the feasibility of extracting 3BFs from heavy-ion collisions as anticipated in Refs. [31, 32, 33].

Report issue for preceding element

I Supplemental Material

Report issue for preceding element

I.1 Scheme dependence in nuclear chiral EFT and phase equivalent NN potentials at N⁴LO⁺

Report issue for preceding element

In this Supplemental Material we focus on the off-shell ambiguities in the chiral nuclear interactions and describe the construction of the phase-equivalent but off-shell differing N⁴LO⁺ NN potentials.

Report issue for preceding element

As emphasized in the main text, scheme dependence is an inherent feature of nuclear interactions. Already the starting point for the derivation of the nuclear interactions, the effective chiral Lagrangian truncated at a given order in the derivative expansion, features a certain degree of ambiguity as reflected in the freedom to perform field redefinitions and apply equations of motion to eliminate various terms. The derivation of nuclear potentials requires off-the-energy-shell extensions of the (few-nucleon-irreducible) parts of the scattering amplitudes and introduces additional scheme dependence. For a description of various methods to derive interactions in chiral EFT see Ref. [26] and references therein. A novel approach using the path integral formalism, which is capable of dealing with regularized (non-local) effective Lagrangians [2], has been proposed in Ref. [3].

Report issue for preceding element

Here, we focus on the expressions for the NN potentials, 3NFs and current operators given in Refs. [4, 5, 6, 7, 8, 9, 10, 28, 11, 12, 13, 14, 15, 16, 17, 18], which were obtained using the method of unitary transformations (UTs). In this approach, one performes a UT of the effective pion-nucleon Hamiltonian to decouple the states on the Fock space involving pions from the purely nucleonic ones. This is achieved in a perturbative framework by utilizing the standard chiral expansion [4, 9]. To arrive at renormalized expressions for the 3BFs and current operators it was necessary, starting from N³LO, to systematically exploit the unitary ambiguity of the nuclear interactions by considering a broad class of transformations on the nucleonic subspace of the Fock space starting from N³LO [8, 10, 11, 14, 15, 16, 17, 18]. The corresponding unitary phases were then chosen in such a way that the resulting potentials remain finite after renormalization (using dimensional regularization). For the considered class of UTs, the resulting static expressions for the nuclear forces at N³LO and N⁴LO turn out to be determined unambiguously, while the leading relativistic corrections depend on two arbitrary real parameters $\bar{\beta}_{8}$, $\bar{\beta}_{9}$ corresponding to the UTs [14]

Report issue for preceding element

$$\hat{U}_{1/m}=e^{-\bar{\beta}_{8}\hat{S}_{8}-\bar{\beta}_{9}\hat{S}_{9}}\,,$$ (S1)

with the anti-Hermitian generators given by

Report issue for preceding element

$$\hat{S}_{8}=\eta\bigg{(}\hat{H}_{\rm kin}\,\eta\,\hat{H}_{g_{A}}\,\frac{\lambda_{1}}{\omega_{\pi}^{3}}\,\hat{H}_{g_{A}}-\mbox{h.c.}\bigg{)}\eta\,,\qquad\hat{S}_{9}=\eta\bigg{(}\hat{H}_{g_{A}}^{1/m}\,\frac{\lambda_{1}}{\omega_{\pi}^{2}}\,\hat{H}_{g_{A}}-\mbox{h.c.}\bigg{)}\eta\,,$$ (S2)

where $\eta$ and $\lambda_{1}$ are the projection operators on the purely nucleonic states and states involving a single pion, respectively, $\omega_{\pi}$ is the pion energy, $\hat{H}_{g_{A}}$ is the leading-order pion-nucleon vertex proportional to the nucleon axial-vector coupling constant $g_{A}$, $\hat{H}_{g_{A}}^{1/m}$ is the leading relativistic correction to that vertex with $m$ denoting the nucleon mass, while $\hat{H}_{\rm kin}$ refers to the kinetic energy operator of the nucleons. The UT in Eq. (S1) induces scheme-dependent contributions to the nuclear force of the form

Report issue for preceding element

$$\delta\hat{V}_{\rm nucl}^{1/m}=\hat{U}_{1/m}(\hat{H}_{\rm kin}+\hat{V}_{\rm nucl})\hat{U}_{1/m}^{\dagger}-\hat{H}_{\rm kin}=\Big{[}\big{(}\hat{H}_{\rm kin}+\hat{V}_{\rm nucl}^{\rm LO}\big{)},\;\big{(}\bar{\beta}_{8}\hat{S}_{8}+\bar{\beta}_{9}\hat{S}_{9}\big{)}\Big{]}+\ldots\,,$$ (S3)

where the ellipses refer to terms beyond N⁴LO and $\hat{V}_{\rm nucl}^{\rm LO}$ denotes the leading contribution to the NN force stemming from the one-pion exchange and derivative-less contact interactions, $\hat{V}_{\rm nucl}^{\rm LO}=\hat{V}_{1\pi}+\hat{V}_{\rm cont}^{\rm LO}$. Evidently, the commutator with the kinetic energy term $\hat{H}_{\rm kin}$ generates $1/m^{2}$ contributions to the one-pion exchange potential, while the one with $\hat{V}_{1\pi}$ induces the $1/m$ corrections to the NN and three-nucleon two-pion exchange potentials $\hat{V}_{2\pi}^{1/m}$ [28]. The resulting unitary ambiguities in the NN potential are well known and discussed in detail in Ref. [27]. Notice further that when considering reactions involving external electroweak probes, the unitary transformation in Eq. (S1) must be applied to the corresponding single-nucleon charge and current operators, yielding scheme-dependent two-nucleon contributions [19, 14, 15]. Approximate independence of electromagnetic NN observables on the unitary phases $\bar{\beta}_{8,9}$ requires using consistent expressions for nuclear potentials and current operators that correspond to the same values of $\bar{\beta}_{8,9}$ and has been explicitly verified in Refs. [20, 21].

Report issue for preceding element

In addition to the long-range UTs in Eq. (S1), nuclear potentials at N³LO feature the short-range ambiguity related to the UTs

Report issue for preceding element

$$\hat{U}_{\rm short-range}=e^{-\gamma_{1}\hat{T}_{1}-\gamma_{2}\hat{T}_{2}-\gamma_{3}\hat{T}_{3}}\,,$$ (S4)

where $\gamma_{i}$ are arbitrary real parameters and the short-range anti-Hermitian NN operators $\hat{T}_{i}$ are given by [19]

Report issue for preceding element

$$\langle{\bf p}_{1}^{\prime}{\bf p}_{2}^{\prime}|\hat{T}_{1}|{\bf p}_{1}{\bf p}_{2}\rangle = \frac{m}{8\Lambda_{\rm b}^{4}}\big{(}p_{1}^{\prime 2}+p_{2}^{\prime 2}-p_{1}^{2}-p_{2}^{2}\big{)}\,,$$

$$\langle{\bf p}_{1}^{\prime}{\bf p}_{2}^{\prime}|\hat{T}_{2}|{\bf p}_{1}{\bf p}_{2}\rangle = \frac{m}{8\Lambda_{\rm b}^{4}}\big{(}p_{1}^{\prime 2}+p_{2}^{\prime 2}-p_{1}^{2}-p_{2}^{2}\big{)}{\bm{\sigma}}_{1}\cdot{\bm{\sigma}}_{2}\,,$$

$$\langle{\bf p}_{1}^{\prime}{\bf p}_{2}^{\prime}|\hat{T}_{3}|{\bf p}_{1}{\bf p}_{2}\rangle = \frac{m}{16\Lambda_{\rm b}^{4}}\Big{(}{\bm{\sigma}}_{1}\cdot({\bf p}_{1}-{\bf p}_{2}+{\bf p}_{1}^{\prime}-{\bf p}_{2}^{\prime})\;{\bm{\sigma}}_{2}\cdot({\bf p}_{1}^{\prime}-{\bf p}_{2}^{\prime}-{\bf p}_{1}+{\bf p}_{2})$$ (S5)

$${}+{\bm{\sigma}}_{1}\cdot({\bf p}_{1}^{\prime}-{\bf p}_{2}^{\prime}-{\bf p}_{1}+{\bf p}_{2})\;{\bm{\sigma}}_{2}\cdot({\bf p}_{1}-{\bf p}_{2}+{\bf p}_{1}^{\prime}-{\bf p}_{2}^{\prime})\Big{)}\,.$$

Here, ${\bf p}_{i}$ and ${\bf p}_{i}^{\prime}$ denote the incoming and outgoing momenta of the nucleon $i$ while $\Lambda_{\rm b}$ is the breakdown scale of chiral EFT, respectively. Note that in line with the commonly used notation, we do not show in Eq. (I.1) the overall factor of $(2\pi)^{3}$ and the $\delta$-function corresponding to the total momentum conservation. We further emphasize that we do not consider here short-ranged UTs, whose generators vanish in the NN cms [22]. In a full analogy to Eq. (S3), the induced contributions to the nuclear forces have the form

Report issue for preceding element

$$\delta\hat{V}_{\rm short-range}=\Big{[}\big{(}\hat{H}_{\rm kin}+\hat{V}_{\rm nucl}^{\rm LO}\big{)},\;\big{(}\gamma_{1}\hat{T}_{1}+\gamma_{2}\hat{T}_{2}+\gamma_{3}\hat{T}_{3}\big{)}\Big{]}+\ldots\,.$$ (S6)

In the two-nucleon system, the commutator with $\hat{V}_{\rm nucl}^{\rm LO}$ leads to short-range NN operators that can be absorbed into a redefinition of the contact interactions. On the other hand, the commutator with $\hat{H}_{\rm kin}$ generates purely off-shell short-range interactions in the NN ¹S₀, ³S₁ and ³S₁- ³D₁ partial waves. Up to N⁴LO, the contact interactions in these channels have the form

Report issue for preceding element

$$\langle^{1}{\rm S}_{0},\,p^{\prime}|\hat{V}_{\rm cont}|^{1}{\rm S}_{0},\,p\rangle = \tilde{C}_{{}^{1}{\rm S}_{0}}+C_{{}^{1}{\rm S}_{0}}(p^{2}+p^{\prime 2})+D_{{}^{1}{\rm S}_{0}}p^{2}p^{\prime 2}+D_{{}^{1}{\rm S}_{0}}^{\rm off}(p^{2}-p^{\prime 2})^{2}\,,$$

$$\langle^{3}{\rm S}_{1},\,p^{\prime}|\hat{V}_{\rm cont}|^{3}{\rm S}_{1},\,p\rangle = \tilde{C}_{{}^{3}{\rm S}_{1}}+C_{{}^{3}{\rm S}_{1}}(p^{2}+p^{\prime 2})+D_{{}^{3}{\rm S}_{1}}p^{2}p^{\prime 2}+D_{{}^{3}{\rm S}_{1}}^{\rm off}(p^{2}-p^{\prime 2})^{2}\,,$$

$$\langle^{3}S_{1},\,p^{\prime}|\hat{V}_{\rm cont}|^{3}D_{1},\,p\rangle = C_{\epsilon_{1}}p^{2}+D_{\epsilon_{1}}p^{2}p^{\prime 2}+D_{\epsilon_{1}}^{\rm off}p^{2}(p^{\prime 2}-p^{2})\,,$$ (S7)

where $\tilde{C}_{i}$ are the LO low-energy constants (LECs), $C_{i}$ are the next-to-leading order (NLO) LECs while $D_{i}$ and $D_{i}^{\rm off}$ are the N³LO LECs. Notice that the contact interactions proportional to $D_{{}^{1}{\rm S}_{0}}^{\rm off}$, $D_{{}^{3}{\rm S}_{1}}^{\rm off}$ and $D_{\epsilon_{1}}^{\rm off}$ vanish in the on-shell kinematics with $p^{\prime}=p=k$ and thus cannot be determined reliably from fits to NN scattering data. This can also be understood from a different perspective by observing that the considered short-ranged UTs generate the NN contact interactions whose form is identical to those proportional to $D_{i}^{\rm off}$ when the commutator in Eq. (S6) is evaluated with $\hat{H}_{\rm kin}$. This implies that the contact terms proportional to the LECs $D_{i}^{\rm off}$ can be rotated away by means of a suitably chosen UT. Accordingly, the NN potentials of Refs. [19, 36] were constructed using the scheme with

Report issue for preceding element

$$D_{{}^{1}{\rm S}_{0}}^{\rm off}=D_{{}^{3}{\rm S}_{1}}^{\rm off}=D_{\epsilon_{1}}^{\rm off}=0\,.$$ (S8)

It is important to emphasize that the considered UTs also induce 3BFs when evaluating the commutator in Eq. (S6) with $\hat{V}_{\rm nucl}^{\rm LO}$ [19, 23]. This demonstrates, once again, that 3BFs can only be meaningfully defined in conjunction with the NN potential, i.e., within a given scheme fixed by a specific off-shell behavior of the NN interaction.

Report issue for preceding element

The choice specified in Eq. (S8) and adopted in Refs. [19, 36] is obviously just one possibility out of infinitely many. From the EFT point of view, any values of $D_{i}^{\rm off}$ are acceptable provided the LECs of the NN contact interactions are of natural size. To illustrate scheme dependence of the 3BF in chiral EFT, we have generated a set of semi-local momentum-space regularized (SMS) NN potentials [19, 36, 24], which are nearly phase equivalent but differ by the choices of the off-shell LECs $D_{{}^{1}{\rm S}_{0}}^{\rm off}$, $D_{{}^{3}{\rm S}_{1}}^{\rm off}$ and $D_{\epsilon_{1}}^{\rm off}$. Specifically, we consider the values

Report issue for preceding element

$$D_{{}^{1}{\rm S}_{0}}^{\rm off}=D_{{}^{3}{\rm S}_{1}}^{\rm off}=\pm 3,\,0;\qquad D_{\epsilon_{1}}^{\rm off}=\pm 1,\,0$$ (S9)

in units of $10^{4}$ GeV^-6 used in Ref. [19] and employ the cutoff value of $\Lambda=450$ MeV. These units coincide with the natural units of $4\pi/(F_{\pi}^{2}\Lambda_{\rm b}^{4})$ used in Refs. [25, 26], where $F_{\pi}$ and $\Lambda_{\rm b}$ denote the pion decay constant and the breakdown scale of chiral EFT, respectively, if one sets $\Lambda_{\rm b}\simeq 620$ MeV. To ensure approximate phase equivalence of the resulting potentials, we stay at the highest available EFT order N⁴LO⁺, which provides a sufficient flexibility to achieve a statistically perfect description of NN data up through the pion production threshold, see Refs. [19, 24] for details, and restrict ourselves to the cutoff value of $\Lambda=450$ MeV. The assumed range of values for $D_{{}^{1}{\rm S}_{0}}^{\rm off}$ and $D_{{}^{3}{\rm S}_{1}}^{\rm off}$ can be considered as fairly conservative with regard to the naturalness assumption. For example, with the off-shell LECs switched off, the resulting N³LO LECs $D_{i}$ fulfill $|D_{i}|\leq 2.7\times 10^{4}$ GeV^-6 for $\Lambda=450$ MeV [19]. For the LEC $D_{\epsilon_{1}}^{\rm off}$, the values $D_{\epsilon_{1}}^{\rm off}=\pm 3$ turn out to cause unrealistically strong modifications of the deuteron D-state probability, which would render the resulting potentials impractical for many-body applications [26]. We have therefore limited ourselves to the values $|D_{\epsilon_{1}}^{\rm off}|\leq 1$.

Report issue for preceding element

For each of the 26 combinations of the off-shell LECs in Eq. (S9) beyond the already available one corresponding to Eq. (S8), we have fitted the remaining LECs accompanying the NN contact interactions to the neutron-proton and proton-proton scattering data up to $E_{\rm lab}=280$ MeV, following the same protocol and using the same database as employed in Refs. [36, 27]. We emphasize that when changing the values of the off-shell LECs and/or the LECs beyond the considered expansion order, the LECs accompanying contact interactions at all orders must be refitted. This procedure corresponds to implicit renormalization within the considered EFT framework and ensures that all results are given in terms of physical parameters (i.e., experimental data for observables used in the fit) rather than bare LECs, see Ref. [26] for details. For all considered cases, we found an essentially perfect description of the neutron-proton (np) and proton-proton scattering data as reflected in the $\chi^{2}/{\rm datum}$-values in the range of $\chi^{2}/{\rm datum}=1.003\ldots 1.007$. These results provide an important consistency check of our calculations and show, in particular, that the contributions in Eq. (I.1) that have been neglected at the considered accuracy level based on the power-counting arguments are indeed numerically negligibly small.

Report issue for preceding element

Figure S1: (Color online). The results for the neutron-proton phase shifts and mixing angles calculated using the 27 N⁴LO⁺ potentials with different choices of the off-shell LECs specified in Eq. (S9) are shown by light-shaded blue bands. Phase shifts and mixing angles obtained using the potential from Ref. [36] with $D_{{}^{1}{\rm S}_{0}}^{\rm off}=D_{{}^{3}{\rm S}_{1}}^{\rm off}=D_{\epsilon_{1}}^{\rm off}=0$ are shown by orange lines. In all cases, the cutoff is set to $\Lambda=450$ MeV.

Report issue for preceding element

In Fig. S1, we show the resulting np phase shifts and mixing angles in low partial waves. These results demonstrate that the 27 potentials can indeed be regarded as essentially phase equivalent. The visible (but tiny) dependence of the mixing angle $\epsilon_{1}$ on the off-shell LECs reflects the fact that this particular observable is less well constrained by the experimental data as compared to other phase shifts. For an estimation of the EFT truncation uncertainty and other types of uncertainties see Refs. [26, 24, 27].

Report issue for preceding element

The predictions for the deuteron properties for the 27 N⁴LO⁺ potentials are collected in Table S1.

Report issue for preceding element

	N4LO+, Ref. [36]	N4LO+, off-shell LECs from Eq. (S9)	Empirical
$B_{d}$ [MeV] xxxxxxxxxxxxx	$2.2246$⋆xxxxxxxxxxxxxxx	$2.2246$⋆	$2.22456614$	$(41)$	[28]
$A_{S}$ [fm-1/2]	$0.8846$	$0.8845\ldots 0.8848$	$0.8845$	$(8)$	[29]
$\eta$	$0.0261$	$0.0260\ldots 0.0263$	$0.0256$	$(4)$	[30]
$r_{m}$ [fm]	$1.9662$	$1.9588\ldots 1.9709$	—
$Q_{0}$ [fm2]	$0.275$	$0.269\ldots 0.280$	—
$P_{D}$ [%]	$4.79$	$3.80\ldots 6.33$	—

⋆The deuteron binding energy has been taken as input in the fit.

Table S1: Deuteron binding energy $B_{d}$, asymptotic S-state normalization $A_{S}$, asymptotic D/S-state ratio $\eta$, matter radius $r_{m}$, leading contribution to the quadrupole moment $Q_{0}$ and D-state probability $P_{D}$ obtained using the SMS N⁴LO⁺ potential of Ref. [36] and the 27 N⁴LO⁺ potentials with the off-shell LECs $D_{i}^{\rm off}$ set according to Eq. (S9) are given in the second and third columns, respectively. The cutoff is set to $\Lambda=450$ MeV.

Report issue for preceding element

The resulting values for the S-state normalization observable $A_{S}$ and the asymptotic D/S-state ratio $\eta$ show, similarly to the phase shifts, a very small sensitivity to the off-shell LECs and agree with the experimental data within errors, see Ref. [24] for the error analysis. On the other hand, the matter radius $r_{m}$ (i.e., the expectation value of the relative distance between the nucleons), the quadrupole moment of the deuteron wave function $Q_{0}$ and the D-state probability $P_{D}$ are well known not to correspond to observable quantities, see e.g. [31]. Not surprisingly, our predictions for these quantities demonstrate a significant degree of scheme dependence. In contrast to $r_{m}$ and $Q_{0}$, the deuteron charge radius and quadrupole moment related to the electric and quadrupole form factors $G_{C}(q^{2})$ and $G_{Q}(q^{2})$, respectively, via $r_{c}=6\frac{1}{G_{C}(0)}\frac{dG_{C}(q^{2})}{dq^{2}}\big{|}_{q^{2}=0}$ and $Q=\frac{1}{m_{d}^{2}}G_{Q}(0)$ with $m_{d}$ denoting the deuteron mass are, of course, observables. The scheme-dependent quantities $r_{m}$ and $Q_{0}$ provide the bulk contributions to the charge radius and the quadrupole moment as they originate from the dominant single-nucleon charge density, but one also needs to take into account two-body contributions to the charge density operator, which are often referred to as “meson-exchange currents”. These (scheme-dependent) meson-exchange currents must be chosen consistently with the nuclear interactions to ensure that the resulting form factors are unambiguous. Indeed, it is easy to see that the unitary transformation in Eq. (S4) generates short-range meson-exchange currents when acting on the single-nucleon charge density [21], rendering the calculated form factors independent of the arbitrary phases $\gamma_{i}$. Note that the spread in the obtained results for $r_{m}$ and $Q_{0}$ agrees qualitatively with the expected size of meson-exchange contributions from the power-counting, see Refs. [21] for details.

Report issue for preceding element

Finally, in Fig. S2, we show the results for the np and neutron-deuteron (nd) total cross section in comparison with the experimental data. As already evident from the results for phase shifts, the np cross section comes out nearly identical for all 27 potentials with the maximal relative differences below two-permille level. On the other hand, the nd cross section calculated using the NN forces only shows a significant scheme dependence, as already visible from Fig. 4.

Report issue for preceding element

Figure S2: (Color online). Neutron-proton and neutron-deuteron total cross section as a function of the cms momentum $k_{\rm cms}$. The upper and lower inset plots show the relative cross section differences for nd and np scattering, respectively, defined as $\delta\sigma\equiv(\sigma-\sigma_{\rm V_{NN,\,[36]}})/\sigma_{\rm V_{NN,\,[36]}}$. Light-shaded blue bands (not visible for neutron-proton scattering) are obtained from the 27 nearly phase-equivalent N⁴LO⁺ NN potentials and illustrate inherent scheme dependence of nuclear interactions in chiral EFT. Experimental data for the nd and np total cross section are taken from Refs. [34] and [34], respectively.

Report issue for preceding element

The difference between the nd experimental data and calculations utilizing a particular NN potential gives the contribution of the 3NF needed to reproduce the cross section data. The large spread in the predictions for three-nucleon observables from different but nearly phase equivalent N⁴LO⁺ NN potentials, reflected in the light-shaded blue bands in Figs. 4 and S2, demonstrates the strong inherent scheme dependence of the 3BF in the framework of chiral EFT.

Report issue for preceding element

The results for the nd total cross section and the ³H binding energy shown in Figs. 4 and S2 agree well with estimations based on the power counting, see Ref. [32] for a related discussion. For example, given the expansion parameter in chiral EFT of the order of $Q\sim 1/3$ in chiral EFT and the expectation value of the NN interaction for ³H of $\langle V_{\rm NN}\rangle\sim 40\ldots 50$ MeV, one may expect the 3BF contribution to the ³H binding energy roughly of the order of $\langle V_{\rm 3NF}\rangle\sim Q^{3}\,\langle V_{\rm NN}\rangle\sim 1.5$ MeV, where we took into account the first appearance of the 3BF at N²LO ($Q^{3}$). Furthermore, the approximately constant deviations between the experimental data for the nd total cross sections and the blue band comprising the results based on the 27 N⁴LO⁺ NN potentials in Fig. 4 translate into the relative deviation that grows with energy and reaches about $\sim 10\ldots 20$% at $E_{\rm lab}=200$ MeV, which corresponds to $k_{\rm cms}\sim 400$ MeV. This pattern and the amount of underprediction of $\sigma^{\rm nd}$ are also well in line with the expected size of the 3BF contributions based on the chiral power counting [32]. In particular, the growing contributions of the 3NF at increasing energies, also seen in the context of the nuclear equation of state, are consistent with the assumed EFT expansion parameter $Q\sim k_{\rm cms}/\Lambda_{\rm b}$. Last but not least, we emphasize that the leading 3BF at N²LO was found to increase the nd total cross section in Ref. [33], thereby improving the agreement between theory and experiment.

Report issue for preceding element