#  Properties of chiral nucleon-nucleon interaction at N3LO with high cutoffs
studied by local projection

H.Y. Shang  State Key Laboratory of Nuclear Physics and Technology, School of
Physics, Peking University, Beijing 100871, China  R.Z. Hu  State Key
Laboratory of Nuclear Physics and Technology, School of Physics, Peking
University, Beijing 100871, China  J.C. Pei<EMAIL_ADDRESS>State Key
Laboratory of Nuclear Physics and Technology, School of Physics, Peking
University, Beijing 100871, China  Southern Center for Nuclear-Science Theory
(SCNT), Institute of Modern Physics, Chinese Academy of Sciences, Huizhou
516000, China  F.R. Xu<EMAIL_ADDRESS>State Key Laboratory of Nuclear
Physics and Technology, School of Physics, Peking University, Beijing 100871,
China  Southern Center for Nuclear-Science Theory (SCNT), Institute of Modern
Physics, Chinese Academy of Sciences, Huizhou 516000, China

###### Abstract

The chiral nucleon-nucleon ($\displaystyle NN$) interaction at high cutoffs
has been plagued by the presence of spurious bound states. In this work, the
chiral $\displaystyle NN$ interaction at N3LO is studied by the local
projection method as the cutoff increases. The evolution of short-range
behaviors of pion-exchange interactions and contact interactions is
intuitively demonstrated. The $\displaystyle P$-channel potentials toward high
cutoffs appear to be erratic at short ranges to compromise with phase shifts,
while such erratic behaviors can be avoided in $\displaystyle S$ and
$\displaystyle D$ channels. Furthermore, a chiral $\displaystyle NN$
interaction at N3LO is studied at a cutoff of 700 MeV. The properties of
deuteron and triton are testified with this interaction. Such a hard
interaction is expected to provide an alternative choice for studies of short-
range correlations and high density nuclear matter.

##  I Introduction

The chiral effective field theory ($\displaystyle\chi$EFT) presents a major
breakthrough to derive modern nuclear forces [1, 2, 3]. $\displaystyle\chi$EFT
provides a perturbation solution of low energy QCD and maintains its crucial
feature of chiral symmetry. Within $\displaystyle\chi$EFT, the nuclear force
invokes the pion degrees of freedoms and many-body forces can be generated in
a systematic expansion. The three-nucleon forces appear at N2LO and are
important for descriptions of nuclear properteis [4]. The applications of
$\displaystyle\chi$EFT nuclear forces have been very successful in ab initio
calculations of nuclear structures [5, 6, 7, 8, 9, 10, 11, 12, 13], few-
nucleon reactions [14, 15, 16, 17, 13] and nuclear matter [18, 19, 20, 21].

The nuclear force based on $\displaystyle\chi$EFT has been conventionally
written in terms of ($\displaystyle Q/\Lambda_{\chi}$)v, where $\displaystyle
Q$ denotes the momentum and $\displaystyle\Lambda_{\chi}$ is the chiral
symmetry broken scale, and $\displaystyle v$ is the order of the Weinberg
power counting [22, 23]. In addition, an explicit cutoff
$\displaystyle\Lambda$ in chiral interaction is realized by multiplying a
regulator function. Currently the high quality chiral nuclear interaction up
to N4LO has been developed [24, 25]. The chiral nuclear interactions usually
adopt a finite cutoff around 500 MeV and the $\displaystyle NN$ scattering
data are fitted below 300 MeV [26, 27, 28, 24]. In the framework of Weinberg
power counting, the chiral interactions are known to have regulator dependence
[29, 25].

In Weinberg’s power counting, the potential is obtained by summing irreducible
diagrams in orders of $\displaystyle Q/\Lambda_{\chi}$, followed by a non-
perturbative resummation in the Lippmann-Schwinger equation to infinite order
to yield the amplitude [30]. This can be a problem for renormalization of
chiral interactions [31, 32]. With a high cutoff, the limit-cycle behaviors of
parameters [33] and unphysical deep bound states can emerge in some partial
waves, which are related to the presence of attractive singular potentials.
Such behaviors have been studied at LO and N2LO before [30, 34]. However, the
emergence of spurious states at high cutoffs at N3LO have not been studied
comprehensively yet.

The chiral interaction is naturally developed in momentum space. The local
projection method provides a direct visualization of the range-dependent
features of the non-local momentum interactions. This has been successfully
applied in studying the softening behavior of nuclear interactions in
similarity renormalization group evolution [35]. This method is expected to be
useful to study the cutoff dependence of chiral interactions, since it is too
complicated to study the limit-cycle behaviors and deep bound states at N3LO
analytically. For nuclear interactions, not only the phase shift but also the
behaviors of nuclear interactions in coordinate space should be considered.

It is always curious to explore realistic nuclear forces at high cutoffs
although soft interactions are favorable in many-body calculations. Actually
$\displaystyle NN$ scattering data in the elastic channel is available up to 3
GeV [36]. The $\displaystyle\Delta$ degree of freedom would play a role at
high cutoffs but the role of full chiral interactions with
$\displaystyle\Delta$ is not conclusive in actual calculations [37, 38, 1].
The role of three-nucleon forces is different in hard $\displaystyle NN$
interactions [24]. The high energy electron scattering measurement indicates
that hard nuclear interactions such as Argonne-V18 [39] can better describe
short-range correlations [40, 41]. In addition, it would be more realistic to
describe nuclear matter at higher densities using nuclear interactions with
high cutoffs.

In this work, we study the behaviors of chiral nuclear interactions at N3LO
using the local projection method while the cutoff varies. Thus the emergence
of spurious deep bound states can be intuitively demonstrated, which is
difficult to be studied analytically at N3LO. The optimization of chiral
$\displaystyle NN$ interactions at N3LO is also performed with a high cutoff
of 700 MeV. The short range features of this interaction are also studied with
local projection. The properties of deuteron and triton are calculated with
this hard interaction.

##  II Methods

###  II.1 The chiral NN potential

Generally the nuclear potential is the sum of contact terms $\displaystyle
V_{\rm ct}$ and pion-exchange terms $\displaystyle V_{\rm pe}$ [1]:

$V\left({\vec{p}}~{}^{\prime},\vec{p}\right)=V_{\rm
ct}\left({\vec{p}}~{}^{\prime},\vec{p}\right)+V_{\rm
pe}\left({\vec{p}}~{}^{\prime},\vec{p}\right),$ (1)

where $\displaystyle{\vec{p}}~{}^{\prime}$ and $\displaystyle\vec{p}$ denote
the final and initial nucleon momenta in the center-of-mass (CMS) system.

In Weinberg’s power counting, the potential is expanded in terms of powers of
small external momenta over a large scale as
$\displaystyle(Q/\Lambda_{\chi})^{\nu}$. The pion-exchange potential is the
sum of various pion-exchange contributions. The leading order is abbreviated
as LO, and next-to leading order as NLO, etc. The contact potential is the sum
of even polynomials of the momentum as:

$\displaystyle V_{\rm ct}=$ $\displaystyle
C_{S}+C_{1}q^{2}+C_{2}k^{2}+D_{1}q^{4}+D_{2}k^{4}$ (2)
$\displaystyle+D_{3}q^{2}k^{2}+D_{4}\left(\vec{q}\times\vec{k}\right)^{2}+\cdots,$

where the momentum transfer is denoted by
$\displaystyle\vec{q}\equiv{\vec{p}}~{}^{\prime}-\vec{p}$ and the average
momentum by
$\displaystyle\vec{k}\equiv\left({\vec{p}}~{}^{\prime}+\vec{p}\right)/2$. In
coordinate spaces, the contact potential is associated with delta function and
derivatives of delta function.

The specific form of the pion-exchange component is given in Ref. [42, 24].
For parameters of the axial-vector coupling constant $\displaystyle g_{A}$,
pion-decay constant $\displaystyle f_{\pi}$, pion masses, nucleon masses and
$\displaystyle\pi N$ low-energy constants (LECs), we adopt the values in Table
I and II of Ref. [24]. To study the contributions of different types of
diagrams, the two-pion exchange (2PE) terms within
$\displaystyle\mathrm{N^{3}LO}$ ($\displaystyle V_{2\pi}^{(4)}$) are further
decomposed [42] to facilitate discussions. To be specific, the $\displaystyle
2\pi$ contribution is decomposed into one-loop and two-loop (tl)
contributions. The one-loop contribution includes terms of $\displaystyle
c_{i}^{2}$ (football diagram or fd), corrections related to $\displaystyle
c_{i}/M_{N}$ (cm), and leading relativistic corrections (rc), as given in
Appendix A. The contact terms under the partial wave representation are
summerized in Appendix B. The specific charge dependency will be discussed in
subsequent sections. The chiral $\displaystyle NN$ interaction can be
expressed as

$\displaystyle V_{\mathrm{LO}}$
$\displaystyle=V_{1\pi}+V_{\mathrm{ct}}^{(0)},$ (3) $\displaystyle
V_{\mathrm{NLO}}$
$\displaystyle=V_{\mathrm{LO}}+V_{2\pi}^{(2)}+V_{\mathrm{ct}}^{(2)},$
$\displaystyle V_{\mathrm{N^{2}LO}}$
$\displaystyle=V_{\mathrm{NLO}}+V_{2\pi}^{(3)},$ $\displaystyle
V_{\mathrm{N}^{3}\mathrm{LO}}$
$\displaystyle=V_{\mathrm{N^{2}LO}}+V_{\mathrm{ct}}^{(4)}+V_{2\pi,\mathrm{tl}}$
$\displaystyle+V_{2\pi,\mathrm{rc}}+V_{2\pi,\mathrm{fd}}+V_{2\pi,\rm cm}.$

###  II.2 The NN scattering

In principle, we can substitute the potential $\displaystyle V$ into the
Blankenbeclar-Sugar (BbS) equation to obtain the scattering matrix
$\displaystyle T$, as given in Ref. [43].

The Bbs equation can be transformed by defining

$\widehat{V}\left({\vec{p}}~{}^{\prime},\vec{p}\right)\equiv\frac{1}{(2\pi)^{3}}\sqrt{\frac{M_{N}}{E_{p^{\prime}}}}V\left({\vec{p}}~{}^{\prime},\vec{p}\right)\sqrt{\frac{M_{N}}{E_{p}}}$
(4)

and

$\widehat{T}\left({\vec{p}}~{}^{\prime},\vec{p}\right)\equiv\frac{1}{(2\pi)^{3}}\sqrt{\frac{M_{N}}{E_{p^{\prime}}}}T\left({\vec{p}}~{}^{\prime},\vec{p}\right)\sqrt{\frac{M_{N}}{E_{p}}},$
(5)

where $\displaystyle M_{N}$ is nucleon mass and $\displaystyle E_{p^{\prime}}$
is nucleon energy $\displaystyle\sqrt{M_{N}^{2}+p^{\prime 2}}$. Actually we
solve the non-relativistic Lippmann-Schwinger (LS) equation:

$\displaystyle\widehat{T}\left({\vec{p}}~{}^{\prime},\vec{p}\right)=$
$\displaystyle\widehat{V}\left({\vec{p}}~{}^{\prime},\vec{p}\right)+\int
d^{3}p^{\prime\prime}\widehat{V}\left({\vec{p}}~{}^{\prime},{\vec{p}}~{}^{\prime\prime}\right)$
(6) $\displaystyle\times\frac{M_{N}}{p^{2}-p^{\prime
2}+i\epsilon}\widehat{T}\left({\vec{p}}~{}^{\prime\prime},\vec{p}\right).$

$\displaystyle\widehat{V}$ and $\displaystyle\widehat{T}$ can be regarded as
the non-relativistic potential and scattering matrix, respectively. By
parameterizing the on-shell ($\displaystyle p$=$\displaystyle p^{\prime}$)
$\displaystyle T$-matrix elements, we can obtain the phase shifts. In this
work, $\displaystyle\widehat{V}$ and $\displaystyle\widehat{T}$ are denoted as
$\displaystyle V$ and $\displaystyle T$.

When solving the LS equation, a cutoff on $\displaystyle V$ at a high momenta
has to be applied to avoid infinities in the integral. It is convenient to
choose the regulator function as [1]

$f\left(p^{\prime},p\right)=\exp\left[-\left(p^{\prime}/\Lambda\right)^{2n}-(p/\Lambda)^{2n}\right],$
(7)

and it remains unaffected in the partial wave decomposition. For the pion-
exchange component, we adopt $\displaystyle n=2$ for 2PE and $\displaystyle
n=4$ for one-pion exchange (1PE). The choice of powers $\displaystyle n$ in
contact terms depends on the partial wave and the order.

The numerical methods in solving the LS equation can be found in Ref. [44] and
the specific formulas for the calculation of phase shifts are given in Ref.
[45]. Note that the Stapp convention is adopted in calculations of phase
shifts [45].

###  II.3 Local projection

In partial wave basis, we denote the potential matrix elements
$\displaystyle\mathrm{i}^{L-L^{\prime}}\left\langle
p^{\prime}L^{\prime}SJM|V|pLSJM\right\rangle$ by $\displaystyle
V^{JS}_{L^{\prime}L}\left(p^{\prime},p\right)$ and $\displaystyle\left\langle
r^{\prime}L^{\prime}SJM|V|rLSJM\right\rangle$ by $\displaystyle
V^{JS}_{L^{\prime}L}\left(r^{\prime},r\right)$. The coordinate-space and
momentum-space matrix elements are associated by the Fourier transform

$V^{JS}_{L^{\prime}L}\left(r^{\prime},r\right)=\frac{2}{\pi}\int_{0}^{\infty}p^{\prime
2}dp^{\prime}p^{2}dpj_{L^{\prime}}\left(p^{\prime}r^{\prime}\right)j_{L}(pr)V^{JS}_{L^{\prime}L}(p^{\prime},p),$
(8)

where $\displaystyle j_{L^{\prime}}$ and $\displaystyle j_{L}$ are spherical
Bessel functions. If the matrix elements in coordinate space can be expressed
by

$V^{JS}_{L^{\prime}L}\left(r^{\prime},r\right)=V^{JS}_{L^{\prime}L}(r)\frac{\delta\left(r-r^{\prime}\right)}{r^{\prime
2}},$ (9)

such a potential is a local potential. Otherwise the potential is nonlocal.

Up to N2LO, the pion-exchange terms obtained from calculations of Feynman
diagrams are local. However, the relativistic corrections and the
$\displaystyle c_{i}/M_{N}$ terms appearing at N3LO, the chosen regulator, and
the $\displaystyle\vec{k}$-dependent contact terms, would introduce non-
localities. These non-localities prohibit direct visualizations and it is
difficult to observe range-dependent features of the potential. In this work,
a local projection of the potential is applied, which is described in Ref.
[35].

The local projection can be straightforwardly obtained by integrating out one
of the coordinates $\displaystyle r^{\prime}$ in the coordinate space
potential,

$\displaystyle\mathrm{L}\left[V\left({\vec{r}}~{}^{\prime},\vec{r}\right)\right]$
$\displaystyle=\delta\left({\vec{r}}~{}^{\prime}-\vec{r}\right)\int
d^{3}r^{\prime\prime}V\left({\vec{r}}~{}^{\prime\prime},\vec{r}\right)$ (10)
$\displaystyle\equiv\delta\left({\vec{r}}~{}^{\prime}-\vec{r}\right)\bar{V}(\vec{r}),$

where $\displaystyle\mathrm{L}$ denotes the local projection, and
$\displaystyle\bar{V}$ is the resulted local potential. In the partial wave
basis, the equation above is written as

$\bar{V}^{JS}_{L^{\prime}L}(r)=\int_{0}^{\infty}r^{\prime
2}dr^{\prime}V^{JS}_{L^{\prime}L}\left(r^{\prime},r\right).$ (11)

By substituting this expression into Eq. (8) and integrating over the free
coordinates analytically, we obtain

$\bar{V}^{JS}_{L0}(r)=\int_{0}^{\infty}k^{2}dkj_{L}(kr)V^{JS}_{L0}(k,0)$ (12)

and (for $\displaystyle L,L^{\prime}>0$)

$\bar{V}^{JS}_{L^{\prime}L}(r)=N_{L^{\prime}}\int_{0}^{\infty}dkdk^{\prime}\frac{k^{2}}{k^{\prime}}j_{L}(kr)V^{JS}_{L^{\prime}L}(k^{\prime},k),$
(13)

where 111Compared to Ref. [35], a factor of $\displaystyle i^{L-L^{\prime}}$
is absorbed in the definiation of matrix elements in the momentumm space.

$N_{L^{\prime}}=\frac{4}{\sqrt{\pi}}\frac{\Gamma\left(\frac{L^{\prime}+3}{2}\right)}{\Gamma\left(\frac{L^{\prime}}{2}\right)}$
(14)

and $\displaystyle\Gamma$ denotes the gamma function.

This approach essentially averages the non-locality, retaining the main
features of the potential. By plotting $\displaystyle\bar{V}$ as a function of
$\displaystyle r$, we can provide a visual representation of the potential.
This allows us to study the variations and the counterbalance of contact terms
and pion-exchange contributions as the cutoff increases.

##  III results and discussions

###  III.1 Potential with different cutoffs

A high cutoff in momentum space can be roughly interpreted as a small-radius
cutoff in coordinate space. As elucidated in Refs. [33, 31], an attractive
potential of the form

$V_{\mathrm{L}}(r)=-\frac{\alpha}{2\mu r^{n}}$ (15)

together with the repulsive centrifugal barrier $\displaystyle
l(l+1)/\left(2\mu r^{2}\right)$, compose the radial potential. For
$\displaystyle n=2$ when $\displaystyle\alpha$ is not very small and
$\displaystyle n\geq 3$, $\displaystyle V_{\mathrm{L}}(r)$ dominates at small
distances in low partial waves, leading to a ‘singular’ potential. The
potential in the coordinate space can be separated into a short-range part and
a long-range part with a sharp cutoff $\displaystyle R$,

$V(r)=V_{\mathrm{S}}(R)\theta(R-r)+V_{\mathrm{L}}(r)\theta(r-R).$ (16)

One can ensure the low-energy observables by adjusting $\displaystyle V_{S}$
when $\displaystyle R$ is fixed. As $\displaystyle R\to 0$
($\displaystyle\Lambda\to\infty$), the parameter displays a limit-cycle
behavior [33] and unphysical bound states appear.

For $\displaystyle NN$ potentials, similar situation occurs in low partial
waves. At LO, Ref. [30] demonstrates that the tensor interaction of 1PE is
‘singular’ and spurious bound states occur in ${}^{3}P_{0}$,
${}^{3}S_{1}-{}^{3}D_{1}$ and ${}^{3}D_{2}$ channels. At N2LO, Ref.[34]
illustrates the running behavior of parameters in ${}^{3}P_{1}$ channel and
indicates the emergence of a spurious bound state at
$\displaystyle\Lambda=730$ MeV.

Figure 1:  Local projections of the contact and pion-exchange components of
the N3LO potential in ${}^{3}P_{0}$ channel. (a) Local projections of
different contributions with a cutoff of $\displaystyle\Lambda$=500 MeV, as
described in Eq. (3). (b) Local projection of total nuclear potential with
$\displaystyle\Lambda$=500 MeV. (c) Local projection of total nuclear
potential with $\displaystyle\Lambda$=850 MeV, in which the shadow region is
responsible for the appearance of deep bound states.

At N3LO, based on the chiral $\displaystyle NN$ potential, we analyze the
parameter running and the emergence of spurious bound states in partial waves
with $\displaystyle J\leq 2$, as $\displaystyle\Lambda$ varies from 500 MeV to
1000 MeV. In the N3LO potential, there are four parameters related to
${}^{1}S_{0}$ channel, and eight parameters related to the coupled
${}^{3}S_{1}-{}^{3}D_{1}$ channels. With such parameters in these channels, it
is challenging to fully determine them via phase shifts. Actually, spurious
bound states can always be eliminated by adjusting the parameters in
${}^{1}S_{0}$, ${}^{3}S_{1}-{}^{3}D_{1}$ channels at N3LO below 1 GeV. For
${}^{3}P_{0}$, ${}^{1}P_{1}$, ${}^{3}P_{1}$ channels, each have two relating
parameters. Then we perform a fitting to GWU SP07 partial wave phase shifts
[36] of the $\displaystyle np$ channel at laboratory kinetic energy
($\displaystyle T_{\mathrm{lab}}$) values of 1, 5, 10, 25, 50, 100, 150, 200,
250, and 300 MeV. We determine the two parameters for each $\displaystyle P$
wave by minimizing
$\displaystyle\chi^{2}=(\delta_{\mathrm{N^{3}LO}}-\delta_{\mathrm{GWU}})^{2}$
and calculate the bound states as the cutoff varies. We find that the first
spurious bound states appear in these three $\displaystyle P$ waves around
$\displaystyle\Lambda$ at 760, 800, and 850 MeV, respectively. Due to the
sensitive dependence of the parameters on $\displaystyle\Lambda$ near the
emergence of spurious bound states, it is difficult to precisely determine
parameters in such cases. For $\displaystyle J$=2 waves, no spurious bound
states appear as $\displaystyle\Lambda$ varies from 500 to 1000 MeV. Note
that, for the ${}^{3}D_{2}$ channel which has one related parameter, a
spurious bound state appears at $\displaystyle\Lambda$ around 1100 MeV. For
$\displaystyle\Lambda>1200$ MeV, phase shifts of the ${}^{3}D_{2}$ channel can
not be described at N3LO by adjusting the parameters. For higher partial
waves, spurious bound states do not appear within the considered cutoff range
due to a strong centrifugal barrier.

Figure 2:  Local projections of contributions of the spectral function
representation $\displaystyle V^{\mathrm{NLO}}_{T,\mathrm{SF}}$ (solid lines)
and the remaining polynomial $\displaystyle V^{\mathrm{NLO}}_{T,\mathrm{pol}}$
(dashed lines), which belong to the 2PE tensor term in ${}^{3}P_{0}$ channel
at NLO, see the text for definitions. The red, orange and blue lines
correspond to local projections of potentials with cutoffs of 500 MeV, 1500
MeV, and 7500 MeV, respectively. In the insert, the results with
$\displaystyle\Lambda=7500$ MeV are shown in the range from 0 to 1 fm.  Figure
3:  Local projections of the contact and pion-exchange components of the N3LO
potential in ${}^{1}P_{0}$ channel, similar to Fig. 1. Figure 4:  Local
projections of the contact and pion-exchange components of the N3LO potential
in ${}^{3}P_{1}$ channel, similar to Fig. 1.

###  III.2 Visualization

Based on the local projection method, we expect to observe the variations in
the pion-exchange and contact-term components as the cutoff changes. In Fig.
1, local projections are applied to N3LO potential [24] as the cutoff
$\displaystyle\Lambda$ varies in ${}^{3}P_{0}$ channel. The regulator powers
related to contact terms $\displaystyle C_{{}^{3}P_{0}}$ and $\displaystyle
D_{{}^{3}P_{0}}$ are $\displaystyle n=2$ and $\displaystyle n=3$. Potentials
with $\displaystyle\Lambda=$500 MeV and $\displaystyle\Lambda=$850 MeV are
selected to illustrate cases without and with spurious bound states,
respectively. Note that the potential with $\displaystyle\Lambda=$500 MeV
shown in Fig. 1 refers to EMN500 from Ref. [24].

To verify the local projection, the bound states are calculated in both the
momentum space and the projected local potential. At $\displaystyle\Lambda$ of
850 MeV for ${}^{3}P_{0}$, the first eigenvalues before and after local
projection are -44964 MeV and -41180 MeV, respectively. At
$\displaystyle\Lambda$ of 820 MeV, the eigenvalues before and after local
projection are -40431 MeV and -37531 MeV. This indicates that the local
projection includes the main features of the potential. Note that the fitting
is difficult as $\displaystyle\Lambda$ near the instable region.

At $\displaystyle\Lambda=500$ MeV, which is a safe choice for chiral nuclear
forces, the pion exchange term is observed to predominantly contribute the
attraction part, with a width of about 2 fm and a depth of approximately 400
MeV, as shown in Fig. 1(b). The contribution from the contact term is
primarily repulsive, with a width and scale comparable to the pion exchange
term. The contact and pion exchange terms have opposite signs and
counterbalance each other. The resulting total potential has slight
oscillations with an attraction part in the short range and a slight
attraction around 2.7 fm. The attraction potential in the core would be
overshadowed by the strong centrifugal barrier, and has no bound-state
solutions.

In Fig. 1(c), as the cutoff $\displaystyle\Lambda$ increases, the oscillation
width of the pion exchange potential gradually decreases, while its
oscillation amplitude increases. At $\displaystyle\Lambda=850$ MeV, the width
is about 1.3 fm. The description of phase shift requires that the contact term
is also attractive at short range. Note that the phase shift at this cutoff is
increased by $\displaystyle\pi$ compared to that at 500 MeV. The shadow region
is related to the appearance of bound states, which is also shown but is small
at $\displaystyle\Lambda=$500 MeV. In this case, the total potential is
significantly attractive at short range and results in spurious deep bound
states.

The decomposed potentials from different terms at N3LO are displayed in
Fig.1(a). The pion-exchange contributions, contact contributions are shown
separately. The separated terms are given in Eq. (3). The oscillations in
different terms up to N2LO are similar at short ranges although their signs
and amplitudes are different. All $\displaystyle P$-wave contributions after
the local projection are zero at $\displaystyle r=0.0$ fm due to the Bessel
function in Eq. (13) . At N3LO, the $\displaystyle c_{i}^{2}$ terms(fd) have
strong attraction contribution, while the $\displaystyle c_{i}/M_{N}$
terms(cm) are significantly repulsive. These two contributions are much larger
than all other terms due to the large value of $\displaystyle c_{i}$. The
contributions from the tl and rc terms are relatively small. The different
terms are given in details in the Appendixes.

If the potential in the moment space is not truncated by the regulator, the
Fourier transform of the irreducible 2PE terms is problematic [47]. We can
obtain the potential in coordinate space through its spectral-function
representation [48], as given in Appendix C. As an example, the tensor term of
the pion exchange contribution at NLO can be divided into $\displaystyle
V^{\mathrm{NLO}}_{T,\mathrm{SF}}$ and $\displaystyle
V^{\mathrm{NLO}}_{T,\mathrm{pol}}$, as given in Eq. (57). $\displaystyle
V^{\mathrm{NLO}}_{T,\mathrm{SF}}$ is related to the spectral-function and
converges in the Fourier transformation as the cutoff increases. The remaining
contribution $\displaystyle V^{\mathrm{NLO}}_{T,\mathrm{pol}}$ is a constant
at NLO in the momentum space. Fig. 2 displays the local projections of these
two terms as the cutoff $\displaystyle\Lambda$ varies. The $\displaystyle
V^{\mathrm{NLO}}_{T,\mathrm{SF}}$ term converges numerically in coordinate
space for $\displaystyle r$$\displaystyle>$1 fm at a cutoff of 7500 MeV. The
regulator with a lower cutoff causes this term to oscillate around the
converged value. The dashed-lines denote the remaining contribution
$\displaystyle V^{\mathrm{NLO}}_{T,\mathrm{pol}}$ with the regulator. With
increasing cutoffs, the oscillations in both $\displaystyle
V^{\mathrm{NLO}}_{T,\mathrm{SF}}$ and $\displaystyle
V^{\mathrm{NLO}}_{T,\mathrm{pol}}$ become more significant. In addition, the
oscillations quickly decays at long ranges as $\displaystyle\Lambda$
increases. This is consistent with Fig.1(c). At $\displaystyle\Lambda$=7500
MeV, both terms have significant oscillations only at very short ranges. Or
say the long range part of the potential becomes convergent as
$\displaystyle\Lambda$ increases. Note that the cutoff $\displaystyle\Lambda$
should be smaller than $\displaystyle\Lambda_{\chi}$ $\displaystyle\sim$1 GeV
for physical studies. Generally the oscillation amplitudes in $\displaystyle
V^{\mathrm{NLO}}_{T,\mathrm{pol}}$ is much more significant than that of
$\displaystyle V^{\mathrm{NLO}}_{T,\mathrm{SF}}$. The remaining contribution
is also a polynomial at N3LO, which causes problems in the Fourier transform
as $\displaystyle q$ increases. The remaining contributions can also be
absorbed into contact terms, so that the oscillations can be much suppressed.
In the local chiral potential, only the contribution of spectral-function
representation is included [49].

Similar studies of ${}^{1}P_{1}$ and ${}^{3}P_{1}$ by the local projection
method are shown in Fig. 3 and Fig. 4. We see again that the contributions of
tl and rc are relatively small. For ${}^{1}P_{1}$ channel, the contact term
and the pion-exchange term have the same oscillation width. This is because
the power in the regulator being adopted is the same. From
$\displaystyle\Lambda$=500 to 900 MeV, the oscillation width decreases and the
oscillation amplitude increases. The situation of ${}^{3}P_{1}$ is similar to
that of ${}^{3}P_{0}$. However, the emergence of spurious states in
${}^{1}P_{1}$ is different from ${}^{3}P_{0}$ and ${}^{3}P_{1}$. At a high
cutoff, the contact potential is still in an oppositive oscillation compared
to the pion-exchange potential.

###  III.3 Potential at cutoff $\displaystyle\Lambda=700$ MeV

Figure 5:  The fitting of $\displaystyle np$ phase shifts in ${}^{3}P_{0}$
channel with a cutoff of 700 MeV. The phase shifts up to 300 or 500 MeV are
taken in the fitting respectively. The red lines denote results with
$\displaystyle n$=2 in regulators for both contact and pion-exchange terms.
The green lines refer to results with different $\displaystyle n$ in
regulators for contact terms. The circles denote the experimental data.
Figure 6:  Local projections of the potentials shown in Fig. 5.

Up to date, chiral nuclear forces are usually developed with a cutoff
$\displaystyle\Lambda$ around 500 MeV, and the phase shifts up to 300 MeV are
taken into the fitting procedure [29, 25, 26, 27, 28, 24]. Recently there are
increasing strong interests to probe the short range properties of realistic
nuclear forces [50]. For example, high energy electron scattering is a very
useful tool to study short range correlations and for the verification of the
point-like nucleon models [41]. The cutoff of conventional chiral nuclear
forces is not sufficiently high to describe such measurements. In addition, ab
initio descriptions of dense nuclear matter suffers large uncertainties above
the saturation density [19]. With more observations of gravitational waves
[51, 52], it is of strong interests to explore nuclear forces at a higher
cutoff. In this context, we are developing chiral $\displaystyle NN$ forces at
N3LO by adopting a cutoff at $\displaystyle\Lambda=700$ MeV and taking into
account the phase shift data up to 500 MeV.

In this work, we only adjust the $\displaystyle NN$ contact terms to fit the
low partial wave phase shifts. The scattering experimental phase shifts up to
500 MeV are included. For $\displaystyle T_{\mathrm{lab}}$ energy between 0
and 10 MeV, we adopt data from the Nijmegen multienergy $\displaystyle np$
phase-shift analysis PWA93 [53]. For $\displaystyle T_{\mathrm{lab}}$ above 10
MeV, we adopt the phase shift analysis of GWU SP07. In the range where both
sets of analysis are available, the data differences are very small. We denote
the cutoff in the phase-shift database as $\displaystyle
T^{\mathrm{cut}}_{\mathrm{lab}}$. The optimization is achieved by using the
genetic algorithm [54], which is changed into an OpenMP parallel scheme in
this work.

Figure 7:  Same as Fig. 5, but for ${}^{1}P_{1}$ channel. Figure 8:  Same as
Fig. 6, but for ${}^{1}P_{1}$ channel. Figure 9:  Same as Fig. 5,but for
coupled $\displaystyle{}^{3}S_{1}-{}^{3}D_{1}$ channels. Figure 10:  The
descriptions of $\displaystyle np$ phase shifts for $\displaystyle S$,
$\displaystyle P$ and $\displaystyle D$ channels and mixing parameters
$\displaystyle\varepsilon_{1}$ and $\displaystyle\varepsilon_{2}$. The red
solid, blue dashed lines and filled circles correspond to $\displaystyle NN$
potential with $\displaystyle\Lambda$=700 MeV, EMN500 and experimental data,
respectively. The green solid lines are fitting with different $\displaystyle
n$ in regulators for contact terms, as in Fig. 5.  Figure 11:  Local
projections of $\displaystyle NN$ potentials for $\displaystyle S$,
$\displaystyle P$ and $\displaystyle D$ channels. The red lines correspond to
$\displaystyle NN$ potentials fitted with $\displaystyle\Lambda$=700 MeV. The
blue lines, black lines and green lines correspond to EMN500, N2LOopt, and
local chiral potential, respectively.

Charge dependence in $\displaystyle NN$ interactions must be taken into
account [1]. The charge dependence due to pion-mass splitting in the one-pion
exchange (1PE) and due to nucleon-mass splitting are considered. Besides, the
irreducible $\displaystyle\pi-\gamma$ exchange contribution in the
$\displaystyle np$ channel is included [55]. The charge-dependent
nonderivative $\displaystyle\tilde{C}^{\mathrm{pp}}_{{}^{1}S_{0}}$,
$\displaystyle\tilde{C}^{\mathrm{np}}_{{}^{1}S_{0}}$ and
$\displaystyle\tilde{C}^{\mathrm{nn}}_{{}^{1}S_{0}}$ contact terms are
employed to absorb the remaining corrections. And they are determined by
$\displaystyle NN$ scattering lengths.

Generally, in minimizing the chi-square
($\displaystyle\chi^{2}=(w\times(\delta_{\mathrm{N^{3}LO}}-\delta_{\mathrm{expt}}))^{2}$)
for fitting parameters, we adopt phase shifts at $\displaystyle
T_{\mathrm{lab}}$=1, 5, 10, 25, 50, 100, 150, 200, 250, 300, 350, 400, 450,
and 500 MeV. Due to large uncertainties in phase shifts between 300 and 500
MeV, and the influences of inelastic scattering, the contribution of this part
to the chi-square is multiplied by a small weight coefficient as
$\displaystyle w$=0.3. For $\displaystyle S$-waves with a strong attraction
and significant phase shift values at low energies, additional phase shifts at
$\displaystyle T_{\mathrm{lab}}$=0.01, 0.02,…, 0.1, 0.2,…, 1, 2,…, 10 MeV are
selected in the low-energy range to better describe bound states and
scattering lengths. We adjust charge-dependent parameters
$\displaystyle\tilde{C}^{\mathrm{pp}}_{{}^{1}S_{0}}$ and
$\displaystyle\tilde{C}^{\mathrm{nn}}_{{}^{1}S_{0}}$ to reproduce
$\displaystyle nn$ and $\displaystyle pp$ scattering lengths in
$\displaystyle{}^{1}S_{0}$. Note that in calculations of the $\displaystyle
pp$ scattering length, the Coulomb interaction is taken into account
appropriately via the Vincent-Phatak method [56].

The obtained low-energy scattering parameters are shown in Table. 1. These
scattering parameters are obtained by using the effective range expansion, as
given in Ref. [45]. $\displaystyle a^{C}$ and $\displaystyle r^{C}$ refer to
parameters by concluding Coulomb force, and $\displaystyle a^{N}$ and
$\displaystyle r^{N}$ denote parameters determined with nuclear interactions
only.

Table 1:  Scattering lengths ($\displaystyle a$) and effective ranges ($\displaystyle r$) in units of fm given by the $\displaystyle NN$ potential with $\displaystyle\Lambda$=700 MeV. $\displaystyle a_{pp}^{C}$ and $\displaystyle r_{pp}^{C}$ refer to $\displaystyle pp$ parameters with Coulomb force. $\displaystyle a^{N}$ and $\displaystyle r^{N}$ denote parameters determined without Coulomb potential. The ${}^{3}S_{1}$ effective range parameters $\displaystyle a_{t}$ and $\displaystyle r_{t}$ are also listed.  | $\displaystyle\mathrm{N}^{3}\mathrm{LO}_{700}$ | Empirical | References
---|---|---|---
$\displaystyle a_{pp}^{C}$ | -7.82 | -7.8196(26) | [57]
$\displaystyle r_{pp}^{C}$ | 2.735 | 2.790(14) | [57]
$\displaystyle a_{nn}^{N}$ | -17.1 |  |
$\displaystyle r_{nn}^{N}$ | 2.87 |  |
$\displaystyle a_{nn}^{N}$ | -18.95 | -18.95(40) | [58]
$\displaystyle r_{nn}^{N}$ | 2.80 | 2.75(11) | [59]
$\displaystyle a_{np}$ | -23.83 | -23.740(20) | [45]
$\displaystyle r_{np}$ | 2.70 | 2.77(5) | [45]
$\displaystyle a_{t}$ | 5.424 | 5.419(7) | [45]
$\displaystyle r_{t}$ | 1.769 | 1.753(8) | [45]

The $\displaystyle np$ scattering phase shifts in the ${}^{3}P_{0}$ channel
are shown in Fig. 5. We fit the parameters using phase shifts database cutoff
$\displaystyle T^{\mathrm{cut}}_{\mathrm{lab}}$=300 and 500 MeV with a cutoff
in the regulator of $\displaystyle\Lambda=$ 700 MeV and $\displaystyle n$=2,
obtaining chi-squared per datum values ($\displaystyle\chi^{2}$/datum) of 2.4
and 1.6, respectively. However, by setting $\displaystyle n$=1 for
$\displaystyle C_{{}^{3}P_{0}}$ and $\displaystyle n$=2 for $\displaystyle
D_{{}^{3}P_{0}}$ in contact terms, the fitting with phase shifts between 0 and
300 MeV can achieve a $\displaystyle\chi^{2}$/datum of 0.54. By setting
$\displaystyle n$=1 for $\displaystyle C_{{}^{3}P_{0}}$ and $\displaystyle
n$=3 for $\displaystyle D_{{}^{3}P_{0}}$, the fitting with phase shifts
between 0 and 500 MeV can achieve a $\displaystyle\chi^{2}$/datum of 0.44. The
fitting of phase shifts can be much improved by adopting different regulators
in contact terms and pion-exchange terms.

The local projection of the fitted potential in ${}^{3}P_{0}$ channel is shown
in Fig. 6. We find that by adjusting $\displaystyle n$ the phase shift data
can be well reproduced up to 500 MeV, however, the corresponding potential
exhibits strong oscillations. There is a notable deep attractive part around 2
fm. This is because of the different powers $\displaystyle n$ in contact and
pion-exchange terms, the resulted oscillation widths in the coordinate space
are different, and the oscillations in pion-exchange terms can not be well
absorbed by the contact terms. In this case, although the phase shift can be
well described, the resulting potential has unrealistic oscillations in the
medium range, which is not explicit in the momentum space. In Ref. [1], some
channels are also fitted with different powers $\displaystyle n$ in regulators
for contact and pion-exchange terms, but the oscillation is not a serious
problem with a lower cutoff.

A similar situation occurs in the ${}^{1}P_{1}$ channel. We present the phase
shifts in Fig. 7 and local projections in Fig. 8. For fittings with phase
shifts up to $\displaystyle T^{\mathrm{cut}}_{\mathrm{lab}}$=300 and 500 MeV,
by adopting $\displaystyle n$=1 for $\displaystyle C_{{}^{1}P_{1}}$ and
$\displaystyle n$=3 for $\displaystyle D_{{}^{1}P_{1}}$ in regulators for
contact terms, the phase shifts can be well described. However, the cases with
better phase shifts in Fig. 7 correspond to strong oscillations in Fig. 8.
Note that there is no spurious bound states with such oscillations. For
fitting with data up to 300 MeV, there is a strong attraction at short ranges
and a steep-rising phase shift above 300 MeV. In both ${}^{1}P_{1}$ and
${}^{3}P_{1}$ at high energies, the phase shifts are too repulsive and this is
also shown in the local projection around 1.2 fm.

The phase shifts for the coupled $\displaystyle{}^{3}S_{1}-{}^{3}D_{1}$
channels are shown in Fig. 9. It is found that the fitting up to
$\displaystyle T_{\mathrm{lab}}^{\rm cut}$=300 MeV in ${}^{3}S_{1}$ channel
extrapolates well to 500 MeV. In ${}^{3}D_{1}$ channel, the fitting up to 500
MeV results in a substantial improvement compared to EMN500. Both sets of
parameters show improvements over EMN500 for the mixing parameter
$\displaystyle\varepsilon_{1}$. The ground state of the deuteron is a
$\displaystyle 1^{+}$ state, which depends solely on the coupled
$\displaystyle{}^{3}S_{1}-{}^{3}D_{1}$ channels. In Table 2, we calculate the
deuteron properties using the nuclear forces in Figure 9, with comparison with
the EMN500 results and empirical values. The binding energy is further
optimized by adjusting the nonderivative ${}^{3}S_{1}$ term. Compared to
EMN500, a higher $\displaystyle D$-state probability $\displaystyle P_{D}$
reflects a stronger tensor force in this potential. This is consistent with a
larger mixing parameter $\displaystyle\varepsilon_{1}$ obtained between 50 MeV
and 200 MeV in Fig. 9. Consequently, this results in a larger asymptotic
$\displaystyle D/S$ state $\displaystyle\eta$ and a larger quadrupole moment
$\displaystyle Q_{d}$.

Figure 12:  Calculations of binding energies of ${}^{3}\mathrm{H}$ by NCSM with different HO basis. The $\displaystyle NN$ potentials of EMN500 and our potential with $\displaystyle\Lambda$=700 MeV are used. The convergence of results with different basis truncations $\displaystyle N_{\mathrm{max}}$ as a function of HO frequency is shown.  Table 2: Deuteron properties as predicted by $\displaystyle NN$ potentials at N3LO with $\displaystyle\Lambda$=700 MeV. The phase shifts up to 300 or 500 MeV are taken in the fitting respectively. The properties include binding energy $\displaystyle B_{d}$, asymptotic $\displaystyle S$ state $\displaystyle A_{S}$, asymptotic $\displaystyle D/S$ state $\displaystyle\eta$, structure radius $\displaystyle r_{\mathrm{str}}$, quadrupole moment $\displaystyle Q_{d}$, $\displaystyle D$-state probability $\displaystyle P_{D}$. | 500fit | 300fit | EMN500 | Empirical 222See Table XVIII of Ref. [45] for references; the empirical value for $\displaystyle r_{\mathrm{str}}$ is from Ref. [60].
---|---|---|---|---
$\displaystyle B_{d}$ (MeV) | 2.224602 | 2.224625 | 2.224615 | 2.224575(9)
$\displaystyle A_{S}$ (fm-1/2) | 0.8879 | 0.8873 | 0.8853 | 0.8846(9)
$\displaystyle\eta$ | 0.0260 | 0.0263 | 0.0257 | 0.0256(4)
$\displaystyle r_{\mathrm{str}}$ (fm) | 1.975 | 1.974 | 1.970 | 1.97507(78)
$\displaystyle Q_{d}$ (fm2) | 0.2733 | 0.2768 | 0.271 | 0.2859(3)
$\displaystyle P_{D}$ (%) | 4.35 | 4.60 | 4.15 | -

The phase shifts of $\displaystyle S$, $\displaystyle P$, and $\displaystyle
D$ channels and mixing parameters $\displaystyle\varepsilon_{1}$ and
$\displaystyle\varepsilon_{2}$ are displayed in red lines in Fig. 10, which
are obtained with $\displaystyle\Lambda=700$MeV and $\displaystyle n=2$ in the
regulators. The phase shifts in ${}^{1}D_{2}$ channel can be improved by
adjusting $\displaystyle n$=2.5 in the contact term $\displaystyle
D_{{}^{1}D_{2}}$, but this would result in strong oscillations. In
${}^{3}D_{2}$ channel, the phase shifts at high energies are improved by
adopting $\displaystyle\Lambda$=700 MeV compared to EMN500. For the
${}^{3}P_{2}$ channel, it is difficult to reproduce the phase shifts at low
and high energies simultaneously. The ${}^{3}D_{3}$ partial wave can be very
well described with $\displaystyle\Lambda$=700 MeV. The phase shifts in some
channels can be better fitted by by adjusting the power $\displaystyle n$ in
the contact regulators, which are displayed in green-solid lines.

Fig. 11 shows the local projections of our $\displaystyle NN$ potential
corresponding to Fig. 10. For comparison, we also present the local
projections of EMN500, $\displaystyle\mathrm{N^{2}LO}_{\mathrm{opt}}$ [27] and
the local potential $\displaystyle\mathrm{N^{3}LO}_{\mathrm{Local}}$ [49].
Here the local potential has a cutoff $\displaystyle R_{\pi}=1.0$ fm for pion-
exchange terms and a cutoff $\displaystyle R_{\mathrm{ct}}=0.7$ fm for contact
terms [49]. The nuclear potentials can be different although they all
reproduce the phase shifts. As we have discussed, the nonlocal $\displaystyle
NN$ potential exhibits oscillatory behavior at short ranges. It is difficult
to separate the short-range and long-range contributions distinctly. In
contrast, there is a distinct separation between long-range and contact
contributions in local potentials. Nonlocal $\displaystyle NN$ potentials are
generally softer compared to the local potential. Our potential with
$\displaystyle\Lambda=700$ MeV is harder compared to EMN500. It can be seen
that the $\displaystyle\mathrm{N^{2}LO}_{\mathrm{opt}}$ nuclear force is the
softest amongst them.

In Fig. 12, we show ground-state energy of ${}^{3}\mathrm{H}$, which is
obtained in ab initio no core shell model (NCSM) [6]. The NCSM code is taken
from Ref. [61]. In NCSM calculations, results are dependent on the frequency
$\displaystyle\hbar\omega$ of the harmonic oscillator (HO) basis and also the
truncation parameter $\displaystyle N_{\mathrm{max}}$ of the model space. The
NCSM calculations with our potential ($\displaystyle{\rm N^{3}LO}_{700}$ with
$\displaystyle\Lambda$=700 MeV) converge slower than that with EMN500. To
obtain convergence with an uncertainty of 1$\displaystyle\%$,
$\displaystyle{\rm N^{3}LO}_{700}$ needs $\displaystyle N_{\mathrm{max}}$=24,
while EMN500 needs only $\displaystyle N_{\mathrm{max}}$=20. The optimal HO
frequency in calculations with $\displaystyle{\rm N^{3}LO}_{700}$ (around
$\displaystyle\hbar\omega$=48 MeV) is much larger than the that with EMN500
(around $\displaystyle\hbar\omega$=28 MeV). Such differences between
$\displaystyle{\rm N^{3}LO}_{700}$ and EMN500 are because $\displaystyle{\rm
N^{3}LO}_{700}$ has a higher cutoff. The binding energy with our potential is
7.42 MeV, which is less bound compared to EMN500. The experimental binding
energy of ${}^{3}\mathrm{H}$ is 8.48 MeV. In this case, the role of three-body
force is important. For ${}^{3}P_{0}$, ${}^{1}P_{1}$, and ${}^{1}D_{2}$
channels, if the potentials are replaced by the better fitting of phase shifts
in Fig. 10, the resulted binding energy of ${}^{3}\mathrm{H}$ would be 7.50
MeV. To be more specific, the smaller binding energy by $\displaystyle{\rm
N^{3}LO}_{700}$ is mainly due to the strong repulsive potential in
${}^{3}P_{0}$ and ${}^{3}S_{1}-^{3}D_{1}$. This is because the phase shifts
are fitted up to 500 MeV, which leads to a strong repulsive feature. If the
phase shifts are fitted to 300 MeV for ${}^{3}P_{0}$ and
${}^{3}S_{1}-^{3}D_{1}$, the obtained binding energy of ${}^{3}\mathrm{H}$
would be 7.78 MeV. It seems that the better descriptions of phase shifts at
high energies, the binding energy of ${}^{3}\mathrm{H}$ would be
underestimated using only $\displaystyle NN$ potential.

###  III.4 Discussions

It has been discussed that the renormalization of chiral $\displaystyle NN$
interaction has been an issue [31, 32]. The chiral $\displaystyle NN$
interactions has dependence on the cutoff, as also demonstrated in this work.
The short-range properties are very much related to the regulator function.
Regarding experiments on short-range correlation (SRC), it seems that nuclear
interactions with a high cutoff are needed [50, 41]. This presents a challenge
for the conventional framework of chiral $\displaystyle NN$ interaction. For
the non-local $\displaystyle NN$ interaction in the momentum space, there are
strong oscillations in the coordinate space. Such strong oscillations are more
significant as the cutoff increases. The development of local chiral potential
[62, 63, 64, 49] can avoid such oscillation issues. In our calculations, the
descriptions of $\displaystyle P$ partial waves are not perfect as that of
$\displaystyle S$-wave. Some $\displaystyle D$ waves can also be well
reproduced. There are unphysical deep bound states in $\displaystyle P$
channels as the cutoff increases. There may be some physics missing in the
current framework. It has been argued that the worse descriptions of
$\displaystyle P$ waves can result in better descriptions of nuclear
structures [38]. And the inclusion of $\displaystyle\Delta$ in chiral
interaction at N2LO doesn’t lead to significant differences in calculations of
nuclear matter [38]. Besides, it has been proposed the modified power counting
with the perturbative solving LS equation to resolve the renormalization issue
[65, 66]. The internal structures of nucleons also play a role in SRC and
European Muon Collaboration(EMC) effects [50], impacting our understandings of
nuclear forces and nuclear structures.

##  IV Summery

To summarize, we have investigated the properties of chiral $\displaystyle NN$
interactions at N3LO as the cutoff increases. The cutoff
$\displaystyle\Lambda$ is conventionally taken as around 500 MeV. However,
hard $\displaystyle NN$ interactions with a high cutoff are needed to address
recent experiments on SRC. In this work, we employ the local projection method
to study the range-dependent features of $\displaystyle NN$ interactions in an
intuitively manner. With increasing $\displaystyle\Lambda$, non-local
$\displaystyle NN$ interactions display strong oscillations in coordinate
space. The oscillation width decreases and oscillation amplitude increases as
$\displaystyle\Lambda$ increases. The strong oscillations are related to the
appearance of unphysical bound states in $\displaystyle P$ channels. The local
projections of different pion-exchange terms are analyzed in details. The
remaining polynomial terms in 2PE contributions after the spectral function
representation cause significant oscillations. There are no such unphysical
bound states in $\displaystyle S$ and $\displaystyle D$ waves at N3LO with
cutoffs below 1 GeV.

We also fit a chiral $\displaystyle NN$ interaction at N3LO with a cutoff of
700 MeV. The phase shifts up to 500 MeV have been adopted. As a result, the
phase shifts of $\displaystyle S$ and $\displaystyle D$ waves can be well
reproduced, except for ${}^{1}D_{2}$ channel. However, the fitting of
$\displaystyle P$ waves can not describe phase shifts at low and high energies
simultaneously. The fitting can be improved by choosing different power
$\displaystyle n$ in regulators for contact terms. In this way the better
fitting of phase shifts of $\displaystyle P$ waves would lead to unrealistic
strong oscillations. The $\displaystyle NN$ interaction with 700 MeV cutoff
has been testified by reproducing the properties deuterium with a slightly
larger component of $\displaystyle D$ waves compared to ENM500. For
${}^{3}\mathrm{H}$, the binding energy is calculated to be 7.42 MeV, which is
related to the repulsive potential in ${}^{3}P_{0}$ and
${}^{3}S_{1}-^{3}D_{1}$. Actually the binding energy of ${}^{3}\mathrm{H}$
becomes to be 7.78 MeV if these two channels are fitted to 300 MeV. We also
discussed relevant issues in renormalization and $\displaystyle P$ waves in
the framework of chiral $\displaystyle NN$ interactions.

###### Acknowledgements.

We thank for many helpful discussions with R. Machleidt and D. R. Entem, and
also for helpful comments by Shuang Zhang. This work was supported by the
National Key R$\displaystyle\&$D Program of China (Grant No.
2023YFE0101500,2023YFA1606403) and National Natural Science Foundation of
China under Grants No. 12335007. We also acknowledge the funding support from
the State Key Laboratory of Nuclear Physics and Technology, Peking University
(No. NPT2023ZX01).

##  Appendix A Pion-exchange component

The NN potential is governed by the customary symmetries, hence it manifests
in the general forms for a two-nucleon potential [67]

$\displaystyle V\left({\vec{p}}~{}^{\prime},\vec{p}\right)=$ $\displaystyle
V_{C}+\boldsymbol{\tau}_{1}\cdot\boldsymbol{\tau}_{2}W_{C}+\left[V_{S}+\boldsymbol{\tau}_{1}\cdot\boldsymbol{\tau}_{2}W_{S}\right]\vec{\sigma}_{1}\cdot\vec{\sigma}_{2}$
(17)
$\displaystyle+\left[V_{LS}+\boldsymbol{\tau}_{1}\cdot\boldsymbol{\tau}_{2}W_{LS}\right][-i\vec{S}\cdot(\vec{q}\times\vec{k})]$
$\displaystyle+\left[V_{T}+\boldsymbol{\tau}_{1}\cdot\boldsymbol{\tau}_{2}W_{T}\right]\vec{\sigma}_{1}\cdot\vec{q}\vec{\sigma}_{2}\cdot\vec{q}$
$\displaystyle+\left[V_{\sigma
L}+\boldsymbol{\tau}_{1}\cdot\boldsymbol{\tau}_{2}W_{\sigma
L}\right]\vec{\sigma}_{1}\cdot(\vec{q}\times\vec{k})\vec{\sigma}_{2}\cdot(\vec{q}\times\vec{k}),$

where $\displaystyle\vec{S}=(\vec{\sigma}_{1}+\vec{\sigma}_{2})/2$ is the
total spin, with $\displaystyle\vec{\sigma}_{1}$,
$\displaystyle\vec{\sigma}_{2}$ and
$\displaystyle\boldsymbol{\tau}_{1}\cdot\boldsymbol{\tau}_{2}$ being the spin
and isospin operators.

###  A.1 Leading order

At LO, there is only the 1PE contribution. The charge-independent 1PE is

$V_{1\pi}^{(\mathrm{CI})}\left({\vec{p}}~{}^{\prime},\vec{p}\right)=-\frac{g_{A}^{2}}{4f_{\pi}^{2}}\boldsymbol{\tau}_{1}\cdot\boldsymbol{\tau}_{2}\frac{\vec{\sigma}_{1}\cdot\vec{q}\vec{\sigma}_{2}\cdot\vec{q}}{q^{2}+m_{\pi}^{2}}$
(18)

Note that, for on-shell and in the CMS, there are no relativistic corrections
at any order. As discussed in Ref. [1], by adjusting the value of
$\displaystyle g_{A}$, the familiar expression for 1PE, Eq. (18) is
appropriate at least to fourth order. Therefore, in this paper, 1PE is
represented solely by this term.

###  A.2 Next-to-leading order

The $\displaystyle NN$ diagrams that occur at NLO contribute in the following
way

$\displaystyle W_{C}=$
$\displaystyle\frac{L(\tilde{\Lambda};q)}{384\pi^{2}f_{\pi}^{4}}\left[4m_{\pi}^{2}\left(1+4g_{A}^{2}-5g_{A}^{4}\right)\right.$
(19)
$\displaystyle\left.+q^{2}\left(1+10g_{A}^{2}-23g_{A}^{4}\right)-\frac{48g_{A}^{4}m_{\pi}^{4}}{w^{2}}\right],$
$V_{T}=-\frac{1}{q^{2}}V_{S}=-\frac{3g_{A}^{4}}{64\pi^{2}f_{\pi}^{4}}L(\tilde{\Lambda};q).$
(20)

The loop function $\displaystyle L(\tilde{\Lambda};q)$ is regularized by
spectral-function cutoff $\displaystyle\tilde{\Lambda}=650$ MeV [47, 48]

$L(\tilde{\Lambda};q)=\frac{w}{2q}\ln\frac{\tilde{\Lambda}^{2}\omega^{2}+q^{2}s^{2}+2\tilde{\Lambda}qws}{4m_{\pi}^{2}\left(\tilde{\Lambda}^{2}+q^{2}\right)}$
(21)

with $\displaystyle w=\sqrt{4m_{\pi}^{2}+q^{2}}$ and $\displaystyle
s=\sqrt{\tilde{\Lambda}^{2}-4m_{\pi}^{2}}$.

###  A.3 Next-to-next-to-leading order

The N2LO contribution is given by

$V_{C}=\frac{3g_{A}^{2}}{16\pi
f_{\pi}^{4}}\left[2m_{\pi}^{2}\left(c_{3}-2c_{1}\right)+c_{3}q^{2}\right]\left(2m_{\pi}^{2}+q^{2}\right)A(\tilde{\Lambda};q),$
(22) $W_{T}=-\frac{1}{q^{2}}W_{S}=-\frac{g_{A}^{2}}{32\pi
f_{\pi}^{4}}c_{4}w^{2}A(\tilde{\Lambda};q).$ (23)

The loop function $\displaystyle A(\tilde{\Lambda};q)$ is regularized in the
same way as $\displaystyle L(\tilde{\Lambda};q)$

$A(\tilde{\Lambda};q)=\frac{1}{2q}\arctan\frac{q\left(\tilde{\Lambda}-2m_{\pi}\right)}{q^{2}+2\tilde{\Lambda}m_{\pi}}$
(24)

###  A.4 Next-to-next-to-next-to-leading order

####  A.4.1 Football diagram at N3LO

The football diagram at N3LO generates

$V_{C}=\frac{3}{16\pi^{2}f_{\pi}^{4}}\left\\{\left[\frac{c_{2}}{6}w^{2}+c_{3}\left(2m_{\pi}^{2}+q^{2}\right)-4c_{1}m_{\pi}^{2}\right]^{2}+\frac{c_{2}^{2}}{45}w^{4}\right\\}L(\tilde{\Lambda};q),$
(25)
$W_{T}=-\frac{1}{q^{2}}W_{S}=\frac{c_{4}^{2}}{96\pi^{2}f_{\pi}^{4}}w^{2}L(\tilde{\Lambda};q).$
(26)

####  A.4.2 $\displaystyle c_{i}/M_{N}$ contributions

This contribution consists of diagrams with one vertex proportional to
$\displaystyle c_{i}$ and one $\displaystyle 1/M_{N}$ correction.

$V_{C}=-\frac{g_{A}^{2}L(\tilde{\Lambda},q)}{32\pi^{2}M_{N}f_{\pi}^{4}}\left[\left(c_{2}-6c_{3}\right)q^{4}+4\left(6c_{1}+c_{2}-3c_{3}\right)q^{2}m_{\pi}^{2}+6\left(c_{2}-2c_{3}\right)m_{\pi}^{4}+24\left(2c_{1}+c_{3}\right)m_{\pi}^{6}w^{-2}\right],$
(27)
$W_{C}=-\frac{c_{4}q^{2}L(\tilde{\Lambda},q)}{192\pi^{2}M_{N}f_{\pi}^{4}}\left[g_{A}^{2}\left(8m_{\pi}^{2}+5q^{2}\right)+w^{2}\right],$
(28)
$W_{T}=-\frac{1}{q^{2}}W_{S}=-\frac{c_{4}L(\tilde{\Lambda},q)}{192\pi^{2}M_{N}f_{\pi}^{4}}\left[g_{A}^{2}\left(16m_{\pi}^{2}+7q^{2}\right)-w^{2}\right],$
(29)
$V_{LS}=\frac{c_{2}g_{A}^{2}}{8\pi^{2}M_{N}f_{\pi}^{4}}w^{2}L(\tilde{\Lambda},q),$
(30)
$W_{LS}=-\frac{c_{4}L(\tilde{\Lambda},q)}{48\pi^{2}M_{N}f_{\pi}^{4}}\left[g_{A}^{2}\left(8m_{\pi}^{2}+5q^{2}\right)+w^{2}\right].$
(31)

####  A.4.3 Leading relativistic corrections

The relativistic corrections of the NLO diagrams, count as N3LO and are given
by

$V_{C}=\frac{3g_{A}^{4}}{128\pi
f_{\pi}^{4}M_{N}}\left[\frac{m_{\pi}^{5}}{2w^{2}}+\left(2m_{\pi}^{2}+q^{2}\right)\left(q^{2}-m_{\pi}^{2}\right)A(\tilde{\Lambda};q)\right],$
(32) $W_{C}=\frac{g_{A}^{2}}{64\pi
f_{\pi}^{4}M_{N}}\left\\{\frac{3g_{A}^{2}m_{\pi}^{5}}{2\omega^{2}}+\left[g_{A}^{2}\left(3m_{\pi}^{2}+2q^{2}\right)-2m_{\pi}^{2}-q^{2}\right]\left(2m_{\pi}^{2}+q^{2}\right)A(\tilde{\Lambda};q)\right\\},$
(33) $V_{T}=-\frac{1}{q^{2}}V_{S}=\frac{3g_{A}^{4}}{256\pi
f_{\pi}^{4}M_{N}}\left(5m_{\pi}^{2}+2q^{2}\right)A(\tilde{\Lambda};q),$ (34)
$W_{T}=-\frac{1}{q^{2}}W_{S}=\frac{g_{A}^{2}}{128\pi
f_{\pi}^{4}M_{N}}\left[g_{A}^{2}\left(3m_{\pi}^{2}+q^{2}\right)-w^{2}\right]A(\tilde{\Lambda};q),$
(35) $V_{LS}=\frac{3g_{A}^{4}}{32\pi
f_{\pi}^{4}M_{N}}\left(2m_{\pi}^{2}+q^{2}\right)A(\tilde{\Lambda};q),$ (36)
$W_{LS}=\frac{g_{A}^{2}\left(1-g_{A}^{2}\right)}{32\pi
f_{\pi}^{4}M_{N}}w^{2}A(\tilde{\Lambda};q).$ (37)

####  A.4.4 Leading two-loop contributions

The leading-order 2PE two-loop diagrams are of N3LO. In terms of spectral
functions, the results are

$\operatorname{Im}V_{C}=\frac{3g_{A}^{4}\left(2m_{\pi}^{2}-\mu^{2}\right)}{\pi\mu\left(4f_{\pi}\right)^{6}}\left[\left(m_{\pi}^{2}-2\mu^{2}\right)\left(2m_{\pi}+\frac{2m_{\pi}^{2}-\mu^{2}}{2\mu}\ln\frac{\mu+2m_{\pi}}{\mu-2m_{\pi}}\right)+4g_{A}^{2}m_{\pi}\left(2m_{\pi}^{2}-\mu^{2}\right)\right],$
(38) $\displaystyle\operatorname{Im}W_{C}=$
$\displaystyle\frac{2\kappa}{3\mu\left(8\pi
f_{\pi}^{2}\right)^{3}}\int_{0}^{1}dx\left[g_{A}^{2}\left(\mu^{2}-2m_{\pi}^{2}\right)+2\left(1-g_{A}^{2}\right)\kappa^{2}x^{2}\right]$
(39)
$\displaystyle\left\\{96\pi^{2}f_{\pi}^{2}\left[\left(2m_{\pi}^{2}-\mu^{2}\right)\left(\bar{d}_{1}+\bar{d}_{2}\right)-2\kappa^{2}x^{2}\bar{d}_{3}+4m_{\pi}^{2}\bar{d}_{5}\right]\right.$
$\displaystyle+\left[4m_{\pi}^{2}\left(1+2g_{A}^{2}\right)-\mu^{2}\left(1+5g_{A}^{2}\right)\right]\frac{\kappa}{\mu}\ln\frac{\mu+2\kappa}{2m_{\pi}}+\frac{\mu^{2}}{12}\left(5+13g_{A}^{2}\right)-2m_{\pi}^{2}\left(1+2g_{A}^{2}\right)$
$\displaystyle-3\kappa^{2}x^{2}+6\kappa
x\sqrt{m_{\pi}^{2}+\kappa^{2}x^{2}}\ln\frac{\kappa
x+\sqrt{m_{\pi}^{2}+\kappa^{2}x^{2}}}{m_{\pi}}$
$\displaystyle\left.+g_{A}^{4}\left(\mu^{2}-2\kappa^{2}x^{2}-2m_{\pi}^{2}\right)\left[\frac{5}{6}+\frac{m_{\pi}^{2}}{\kappa^{2}x^{2}}-\left(1+\frac{m_{\pi}^{2}}{\kappa^{2}x^{2}}\right)^{3/2}\ln\frac{\kappa
x+\sqrt{m_{\pi}^{2}+\kappa^{2}x^{2}}}{m_{\pi}}\right]\right\\},$
$\displaystyle\operatorname{Im}V_{S}$
$\displaystyle=\mu^{2}\operatorname{Im}V_{T}=\frac{g_{A}^{2}\mu\kappa^{3}}{8\pi
f_{\pi}^{4}}\left(\bar{d}_{15}-\bar{d}_{14}\right)$ (40)
$\displaystyle+\frac{2g_{A}^{6}\mu\kappa^{3}}{\left(8\pi
f_{\pi}^{2}\right)^{3}}\int_{0}^{1}dx\left(1-x^{2}\right)\left[\frac{1}{6}-\frac{m_{\pi}^{2}}{\kappa^{2}x^{2}}+\left(1+\frac{m_{\pi}^{2}}{\kappa^{2}x^{2}}\right)^{3/2}\ln\frac{\kappa
x+\sqrt{m_{\pi}^{2}+\kappa^{2}x^{2}}}{m_{\pi}}\right],$
$\operatorname{Im}W_{S}=\mu^{2}\operatorname{Im}W_{T}(i\mu)=\frac{g_{A}^{4}\left(4m_{\pi}^{2}-\mu^{2}\right)}{\pi\left(4f_{\pi}\right)^{6}}\left[\left(m_{\pi}^{2}-\frac{\mu^{2}}{4}\right)\ln\frac{\mu+2m_{\pi}}{\mu-2m_{\pi}}+\left(1+2g_{A}^{2}\right)\mu
m_{\pi}\right],$ (41)

where $\displaystyle\kappa=\sqrt{\mu^{2}/4-m_{\pi}^{2}}$.

The potentials in (subtracted) spectral representation are given in momentum
space by

$V_{C,S}(q)=-\frac{2q^{6}}{\pi}\int_{2m_{\pi}}^{\tilde{\Lambda}}d\mu\frac{\operatorname{Im}V_{C,S}(i\mu)}{\mu^{5}\left(\mu^{2}+q^{2}\right)},$
(42)
$V_{T,LS}(q)=\frac{2q^{4}}{\pi}\int_{2m_{\pi}}^{\tilde{\Lambda}}d\mu\frac{\operatorname{Im}V_{T,LS}(i\mu)}{\mu^{3}\left(\mu^{2}+q^{2}\right)}.$
(43)

##  Appendix B Contact terms

The contact terms are given in the general form Eq. (17), with $\displaystyle
V_{\alpha}$ and $\displaystyle W_{\alpha}$
$\displaystyle(\alpha=C,S,LS,T,\sigma L)$ adopting polynomials in terms of
$\displaystyle\vec{k}^{2}$, $\displaystyle\vec{q}^{2}$, and
$\displaystyle\left(\vec{k}\times\vec{q}\right)^{2}$. Terms that include a
factor $\displaystyle\boldsymbol{\tau}_{1}\cdot\boldsymbol{\tau}_{2}$ can be
left out due to Fierz ambiguity.

To be specific, the LO contact potential is given by

$V_{\mathrm{ct}}^{(0)}\left({\vec{p}}~{}^{\prime},\vec{p}\right)=C_{S}+C_{T}\vec{\sigma}_{1}\cdot\vec{\sigma}_{2}$
(44)

and, in terms of partial waves,

$V_{\mathrm{ct}}^{(0)}\left({}^{1}S_{0}\right)=\widetilde{C}_{{}^{1}S_{0}}=4\pi\left(C_{S}-3C_{T}\right),$
(45)
$V_{\mathrm{ct}}^{(0)}\left({}^{3}S_{1}\right)=\widetilde{C}_{{}^{3}S_{1}}=4\pi\left(C_{S}+C_{T}\right).$
(46)

At NLO, we have

$\displaystyle
V_{\mathrm{ct}}^{(2)}\left({\vec{p}}~{}^{\prime},\vec{p}\right)=$
$\displaystyle
C_{1}q^{2}+C_{2}k^{2}+\left(C_{3}q^{2}+C_{4}k^{2}\right)\vec{\sigma}_{1}\cdot\vec{\sigma}_{2}$
(47)
$\displaystyle+C_{5}[-i\vec{S}\cdot(\vec{q}\times\vec{k})]+C_{6}\left(\vec{\sigma}_{1}\cdot\vec{q}\right)\left(\vec{\sigma}_{2}\cdot\vec{q}\right)$
$\displaystyle+C_{7}\left(\vec{\sigma}_{1}\cdot\vec{k}\right)\left(\vec{\sigma}_{2}\cdot\vec{k}\right),$

and partial-wave decomposition yields

$\displaystyle V_{\mathrm{ct}}^{(2)}\left({}^{1}S_{0}\right)$
$\displaystyle=C_{{}^{1}S_{0}}\left(p^{2}+p^{\prime 2}\right),$ (48)
$\displaystyle V_{\mathrm{ct}}^{(2)}\left({}^{3}P_{0}\right)$
$\displaystyle=C_{{}^{3}P_{0}}pp^{\prime},$ $\displaystyle
V_{\mathrm{ct}}^{(2)}\left({}^{1}P_{1}\right)$
$\displaystyle=C_{{}^{1}P_{1}}pp^{\prime},$ $\displaystyle
V_{\mathrm{ct}}^{(2)}\left({}^{3}P_{1}\right)$
$\displaystyle=C_{{}^{3}P_{1}}pp^{\prime},$ $\displaystyle
V_{\mathrm{ct}}^{(2)}\left({}^{3}S_{1}\right)$
$\displaystyle=C_{{}^{3}S_{1}}\left(p^{2}+p^{\prime 2}\right),$ $\displaystyle
V_{\mathrm{ct}}^{(2)}\left({}^{3}S_{1}-{}^{3}D_{1}\right)$
$\displaystyle=C_{{}^{3}S_{1}-{}^{3}D_{1}}p^{2},$ $\displaystyle
V_{\mathrm{ct}}^{(2)}\left({}^{3}D_{1}-{}^{3}S_{1}\right)$
$\displaystyle=C_{{}^{3}S_{1}-{}^{3}D_{1}}p^{\prime 2},$ $\displaystyle
V_{\mathrm{ct}}^{(2)}\left({}^{3}P_{2}\right)$
$\displaystyle=C_{{}^{3}P_{2}}pp^{\prime}.$

The N3LO contacts are

$\displaystyle V_{\mathrm{ct}}^{(4)}$
$\displaystyle\left({\vec{p}}~{}^{\prime},\vec{p}\right)$ (49)
$\displaystyle=$ $\displaystyle
D_{1}q^{4}+D_{2}k^{4}+D_{3}q^{2}k^{2}+D_{4}(\vec{q}\times\vec{k})^{2}$
$\displaystyle+\left[D_{5}q^{4}+D_{6}k^{4}+D_{7}q^{2}k^{2}+D_{8}(\vec{q}\times\vec{k})^{2}\right]\vec{\sigma}_{1}\cdot\vec{\sigma}_{2}$
$\displaystyle+\left(D_{9}q^{2}+D_{10}k^{2}\right)[-i\vec{S}\cdot(\vec{q}\times\vec{k})]$
$\displaystyle+\left(D_{11}q^{2}+D_{12}k^{2}\right)\left(\vec{\sigma}_{1}\cdot\vec{q}\right)\left(\vec{\sigma}_{2}\cdot\vec{q}\right)$
$\displaystyle+\left(D_{13}q^{2}+D_{14}k^{2}\right)\left(\vec{\sigma}_{1}\cdot\vec{k}\right)\left(\vec{\sigma}_{2}\cdot\vec{k}\right)$
$\displaystyle+D_{15}\left[\vec{\sigma}_{1}\cdot(\vec{q}\times\vec{k})\vec{\sigma}_{2}\cdot(\vec{q}\times\vec{k})\right],$

with contributions by partial waves,

$\displaystyle
V_{\mathrm{ct}}^{(4)}\left({}^{1}S_{0}\right)=\widehat{D}_{{}^{1}S_{0}}\left(p^{\prime
4}+p^{4}\right)+D_{{}^{1}S_{0}}p^{\prime 2}p^{2},$ (50) $\displaystyle
V_{\mathrm{ct}}^{(4)}\left({}^{3}P_{0}\right)=D_{{}^{3}P_{0}}\left(p^{\prime
3}p+p^{\prime}p^{3}\right)\text{, }$ $\displaystyle
V_{\mathrm{ct}}^{(4)}\left({}^{1}P_{1}\right)=D_{{}^{1}P_{1}}\left(p^{\prime
3}p+p^{\prime}p^{3}\right),$ $\displaystyle
V_{\mathrm{ct}}^{(4)}\left({}^{3}P_{1}\right)=D_{{}^{3}P_{1}}\left(p^{\prime
3}p+p^{\prime}p^{3}\right),$ $\displaystyle
V_{\mathrm{ct}}^{(4)}\left({}^{3}S_{1}\right)=\widehat{D}_{{}^{3}S_{1}}\left(p^{\prime
4}+p^{4}\right)+D_{{}^{3}S_{1}}p^{\prime 2}p^{2}\text{, }$ $\displaystyle
V_{\mathrm{ct}}^{(4)}\left({}^{3}D_{1}\right)=D_{{}^{3}D_{1}}p^{\prime
2}p^{2},$ $\displaystyle
V_{\mathrm{ct}}^{(4)}\left({}^{3}S_{1}-{}^{3}D_{1}\right)=\widehat{D}_{{}^{3}S_{1}-{}^{3}D_{1}}p^{4}+D_{{}^{3}S_{1}-{}^{3}D_{1}}p^{\prime
2}p^{2}\text{, }$ $\displaystyle
V_{\mathrm{ct}}^{(4)}\left({}^{3}D_{1}-{}^{3}S_{1}\right)=\widehat{D}_{{}^{3}S_{1}-{}^{3}D_{1}}p^{\prime
4}+D_{{}^{3}S_{1}-{}^{3}D_{1}}p^{\prime 2}p^{2}\text{, }$ $\displaystyle
V_{\mathrm{ct}}^{(4)}\left({}^{1}D_{2}\right)=D_{{}^{1}D_{2}}p^{\prime
2}p^{2},$ $\displaystyle
V_{\mathrm{ct}}^{(4)}\left({}^{3}D_{2}\right)=D_{{}^{3}D_{2}}p^{\prime
2}p^{2}\text{, }$ $\displaystyle
V_{\mathrm{ct}}^{(4)}\left({}^{3}P_{2}\right)=D_{{}^{3}P_{2}}\left(p^{\prime
3}p+p^{\prime}p^{3}\right),$ $\displaystyle
V_{\mathrm{ct}}^{(4)}\left({}^{3}P_{2}-{}^{3}F_{2}\right)=D_{{}^{3}P_{2}-{}^{3}F_{2}}p^{\prime}p^{3},$
$\displaystyle
V_{\mathrm{ct}}^{(4)}\left({}^{3}F_{2}-{}^{3}P_{2}\right)=D_{{}^{3}P_{2}-{}^{3}F_{2}}p^{\prime
3}p,$ $\displaystyle
V_{\mathrm{ct}}^{(4)}\left({}^{3}D_{3}\right)=D_{{}^{3}D_{3}}p^{\prime
2}p^{2}.$

The coefficients $\displaystyle C_{{}^{(2S+1)}L_{J}}$ and $\displaystyle
C_{i}$, as well as $\displaystyle D_{{}^{(2S+1)}L_{J}}$ and $\displaystyle
D_{i}$, can be transformed into each other through invertible matrixs. The
specific matrices for these transformations can be found in Ref. [1].

##  Appendix C Spectral-function representation

We denote the tensor term of the pion exchange contribution at NLO (Eq. (20))
as $\displaystyle V^{\mathrm{NLO}}_{T}$. Here we demonstrate the process of
obtaining its spectral-function representation. Our starting point is the
potential obtained through dimensional regularization

$V_{T}(q)=-\frac{3g_{A}^{4}}{64\pi^{2}f_{\pi}^{4}}L(q),$ (51)

where the loop function is

$L(q)=\frac{w}{q}\ln\frac{w+q}{2m_{\pi}}.$ (52)

The corresponding mass spectrum $\displaystyle\rho_{T}(\mu)$ (spectral
function) is

$\rho_{T}(\mu)=\operatorname{Im}V_{T}(i\mu)=\frac{3g_{A}^{4}}{128\pi
f_{\pi}^{4}}\frac{\sqrt{\mu^{2}-4m_{\pi}^{2}}}{\mu}.$ (53)

Then, the spectral-function representation of this term is obtained through a
continuous superposition of Yukawa functions

$V_{T}^{SF}=\frac{2}{\pi}\int_{2M_{\pi}}^{\tilde{\Lambda}}\mathrm{d}\mu\mu\frac{\rho_{T}(\mu)}{\mu^{2}+q^{2}}$
(54)

We do not use the subtracted spectral representation [68] as in Eq. (43),
because the Fourier transform of Eq. (54) directly results in the coordinate
space potential

$\displaystyle\widetilde{V}_{T}(r)=$
$\displaystyle-\frac{1}{6\pi^{2}r^{3}}\int_{2m_{\pi}}^{\tilde{\Lambda}}d\mu\mu
e^{-\mu r}\left(3+3\mu r+\mu^{2}r^{2}\right)$ (55)
$\displaystyle\times\operatorname{Im}V_{T}(i\mu).$

Thus this part captures the finite range part of the potential in coordinate
space

$V_{T}^{SF}(q)=-\frac{3g_{A}^{4}}{64\pi^{2}f_{\pi}^{4}}\left[L(\tilde{\Lambda};q)-\ln\frac{1}{2}\left(\frac{\tilde{\Lambda}+s}{m_{\pi}}\right)\right].$
(56)

Then $\displaystyle V^{\mathrm{NLO}}_{T}$ can be divided as

$V^{\mathrm{NLO}}_{T}=V^{\mathrm{NLO}}_{T,\mathrm{SF}}+V^{\mathrm{NLO}}_{T,\mathrm{pol}},$
(57)

where the remaining polynomial term $\displaystyle
V^{\mathrm{NLO}}_{T,\mathrm{pol}}$ is given by
$\displaystyle-\frac{3g_{A}^{4}}{64\pi^{2}f_{\pi}^{4}}\ln\frac{1}{2}\left(\frac{\tilde{\Lambda}+s}{m_{\pi}}\right)$.

## References

  * Machleidt and Entem [2011] R. Machleidt and D. R. Entem, Phys. Rep. 503, 1 (2011).
  * Epelbaum  _et al._ [2009] E. Epelbaum, H.-W. Hammer, and U.-G. Meißner, Rev. Mod. Phys. 81, 1773 (2009).
  * Hammer  _et al._ [2020] H.-W. Hammer, S. König, and U. Van Kolck, Rev. Mod. Phys. 92, 025004 (2020).
  * Hebeler [2020] K. Hebeler, Phys. Rep. 890, 1 (2020).
  * Tichai  _et al._ [2020] A. Tichai, R. Roth, and T. Duguet, Front. Phys. 8, 164 (2020).
  * Barrett  _et al._ [2013] B. R. Barrett, P. Navrátil, and J. P. Vary, Prog. Part. Nucl. Phys. 69, 131 (2013).
  * Hagen  _et al._ [2014a] G. Hagen, T. Papenbrock, M. Hjorth-Jensen, and D. J. Dean, Rep. Prog. Phys. 77, 096302 (2014a).
  * Hergert  _et al._ [2016] H. Hergert, S. K. Bogner, T. D. Morris, A. Schwenk, and K. Tsukiyama, Phys. Rep. 621, 165 (2016).
  * Soma [2020] V. Soma, Front. Phys. 8, 340 (2020).
  * Ma  _et al._ [2020] Y. Z. Ma, F. R. Xu, N. Michel, S. Zhang, J. G. Li, B. S. Hu, L. Coraggio, N. Itaco, and A. Gargano, Phys. Lett. B 808, 135673 (2020).
  * Hu  _et al._ [2019] B. S. Hu, Q. Wu, Z. H. Sun, and F. R. Xu, Phys. Rev. C 99, 061302 (2019).
  * Li  _et al._ [2019] J. G. Li, N. Michel, B. S. Hu, W. Zuo, and F. R. Xu, Phys. Rev. C 100, 054313 (2019).
  * Michel  _et al._ [2023] N. Michel, W. Nazarewicz, and M. Płoszajczak, Phys. Rev. Lett. 131, 242502 (2023).
  * Viviani  _et al._ [2013] M. Viviani, L. Girlanda, A. Kievsky, and L. E. Marcucci, Phys. Rev. Lett. 111, 172302 (2013).
  * Epelbaum  _et al._ [2002] E. Epelbaum, A. Nogga, W. Glöckle, H. Kamada, U.-G. Meißner, and H. Witała, Phys. Rev. C 66, 064001 (2002).
  * Navrátil  _et al._ [2010] P. Navrátil, R. Roth, and S. Quaglioni, Phys. Rev. C 82, 034609 (2010).
  * Golak  _et al._ [2014] J. Golak, R. Skibiński, K. Topolnicki, H. Witała, E. Epelbaum, H. Krebs, H. Kamada, U.-G. Meißner, V. Bernard, P. Maris, J. Vary, S. Binder, A. Calci, K. Hebeler, J. Langhammer, R. Roth, A. Nogga, S. Liebig, and D. Minossi, Eur. Phys. J. A 50, 177 (2014).
  * Drischler  _et al._ [2021] C. Drischler, J. Holt, and C. Wellenhofer, Annu. Rev. Nucl. Part. Sci. 71, 403 (2021).
  * Drischler  _et al._ [2019] C. Drischler, K. Hebeler, and A. Schwenk, Phys. Rev. Lett. 122, 042501 (2019).
  * Hagen  _et al._ [2014b] G. Hagen, T. Papenbrock, A. Ekström, K. A. Wendt, G. Baardsen, S. Gandolfi, M. Hjorth-Jensen, and C. J. Horowitz, Phys. Rev. C 89, 014319 (2014b).
  * Rios [2020] A. Rios, Front. Phys. 8, 387 (2020).
  * Weinberg [1990] S. Weinberg, Phys. Lett. B 251, 288 (1990).
  * Weinberg [1991] S. Weinberg, Nucl. Phys. B 363, 3 (1991).
  * Entem  _et al._ [2017] D. R. Entem, R. Machleidt, and Y. Nosyk, Phys. Rev. C 96, 024004 (2017).
  * Reinert  _et al._ [2018] P. Reinert, H. Krebs, and E. Epelbaum, Eur. Phys. J. A 54, 86 (2018).
  * Ekström  _et al._ [2015] A. Ekström, G. R. Jansen, K. A. Wendt, G. Hagen, T. Papenbrock, B. D. Carlsson, C. Forssén, M. Hjorth-Jensen, P. Navrátil, and W. Nazarewicz, Phys. Rev. C 91, 051301 (2015).
  * Ekström  _et al._ [2013] A. Ekström, G. Baardsen, C. Forssén, G. Hagen, M. Hjorth-Jensen, G. R. Jansen, R. Machleidt, W. Nazarewicz, T. Papenbrock, J. Sarich, and S. M. Wild, Phys. Rev. Lett. 110, 192502 (2013).
  * Entem and Machleidt [2003] D. R. Entem and R. Machleidt, Phys. Rev. C 68, 041001 (2003).
  * Epelbaum  _et al._ [2015] E. Epelbaum, H. Krebs, and U. G. Meißner, Eur. Phys. J. A 51, 53 (2015).
  * Nogga  _et al._ [2005] A. Nogga, R. G. E. Timmermans, and U. van Kolck, Phys. Rev. C 72, 054006 (2005).
  * van Kolck [2020] U. van Kolck, Front. Phys. 8, 79 (2020).
  * Entem  _et al._ [2008] D. R. Entem, E. R. Arriola, M. P. Valderrama, and R. Machleidt, Phys. Rev. C 77, 044006 (2008).
  * Beane  _et al._ [2001] S. R. Beane, P. F. Bedaque, L. Childress, A. Kryjevski, J. McGuire, and U. Van Kolck, Phys. Rev. A 64, 042103 (2001).
  * Epelbaum  _et al._ [2004a] E. Epelbaum, W. Glöckle, and U.-G. Meißner, Eur. Phys. J. A 19, 401 (2004a).
  * Wendt  _et al._ [2012] K. A. Wendt, R. J. Furnstahl, and S. Ramanan, Phys. Rev. C 86, 014003 (2012).
  * Arndt  _et al._ [2007] R. A. Arndt, W. J. Briscoe, I. I. Strakovsky, and R. L. Workman, Phys. Rev. C 76, 025209 (2007).
  * Piarulli  _et al._ [2018] M. Piarulli, A. Baroni, L. Girlanda, A. Kievsky, A. Lovato, E. Lusk, L. E. Marcucci, S. C. Pieper, R. Schiavilla, M. Viviani, and R. B. Wiringa, Phys. Rev. Lett. 120, 052503 (2018).
  * Nosyk  _et al._ [2021] Y. Nosyk, D. R. Entem, and R. Machleidt, Phys. Rev. C 104, 054001 (2021).
  * Wiringa  _et al._ [1995] R. B. Wiringa, V. G. J. Stoks, and R. Schiavilla, Phys. Rev. C 51, 38 (1995).
  * Arrington  _et al._ [2022] J. Arrington, N. Fomin, and A. Schmidt, Annual Review of Nuclear and Particle Science 72, 307 (2022).
  * Schmidt  _et al._ [2020] A. Schmidt _et al._ , Nature 578, 540 (2020), (CLAS Collaboration).
  * Entem  _et al._ [2015] D. R. Entem, N. Kaiser, R. Machleidt, and Y. Nosyk, Phys. Rev. C 91, 014002 (2015).
  * Blankenbecler and Sugar [1966] R. Blankenbecler and R. Sugar, Phys. Rev. 142, 1051 (1966).
  * Machleidt [1993] R. Machleidt, One-boson-exchange potentials and nucleon-nucleon scattering, in  _Computational Nuclear Physics 2: Nuclear Reactions_ , edited by K. Langanke, J. A. Maruhn, and S. E. Koonin (Springer, New York, 1993) Chap. 1, p. 1.
  * Machleidt [2001] R. Machleidt, Phys. Rev. C 63, 024001 (2001).
  * Note [1] Compared to Ref. [35], a factor of $\displaystyle i^{L-L^{\prime}}$ is absorbed in the definiation of matrix elements in the momentumm space.
  * Kaiser  _et al._ [1997] N. Kaiser, R. Brockmann, and W. Weise, Nucl. Phys. A 625, 758 (1997).
  * Epelbaum  _et al._ [2004b] E. Epelbaum, W. Gloeckle, and U.-G. Meißner, Eur. Phys. J. A 19, 125 (2004b).
  * Saha  _et al._ [2023] S. K. Saha, D. R. Entem, R. Machleidt, and Y. Nosyk, Phys. Rev. C 107, 034002 (2023).
  * Hen  _et al._ [2017] O. Hen, G. A. Miller, E. Piasetzky, and L. B. Weinstein, Rev. Mod. Phys. 89, 045002 (2017).
  * Abbott  _et al._ [2017] B. P. Abbott _et al._ , Phys. Rev. Lett. 119, 161101 (2017).
  * Abbott  _et al._ [2018] B. P. Abbott _et al._ , Phys. Rev. Lett. 121, 161101 (2018).
  * Stoks  _et al._ [1993] V. G. J. Stoks, R. A. M. Klomp, M. C. M. Rentmeester, and J. J. De Swart, Phys. Rev. C 48, 792 (1993).
  * Charbonneau [2002] P. Charbonneau, An introduction to genetic algorithms for numerical optimization, http://cobweb.cs.uga.edu/potter/CompIntell/no_tutorial.pdf (2002).
  * van Kolck  _et al._ [1998] U. van Kolck, M. C. M. Rentmeester, J. L. Friar, T. Goldman, and J. J. de Swart, Phys. Rev. Lett. 80, 4386 (1998).
  * Vincent and Phatak [1974] C. M. Vincent and S. C. Phatak, Phys. Rev. C 10, 391 (1974).
  * Bergervoet  _et al._ [1988] J. R. Bergervoet, P. C. Van Campen, W. A. Van Der Sanden, and J. J. De Swart, Phys. Rev. C 38, 15 (1988).
  * Chen  _et al._ [2008] Q. Chen, C. R. Howell, T. S. Carman, W. R. Gibbs, B. F. Gibson, A. Hussein, M. R. Kiser, G. Mertens, C. F. Moore, C. Morris, A. Obst, E. Pasyuk, C. D. Roper, F. Salinas, H. R. Setze, I. Slaus, S. Sterbenz, W. Tornow, R. L. Walter, C. R. Whiteley, and M. Whitton, Phys. Rev. C 77, 054002 (2008).
  * Miller  _et al._ [1990] G. A. Miller, B. M. K. Nefkens, and I. Šlaus, Phys. Rep. 194, 1 (1990).
  * Jentschura  _et al._ [2011] U. Jentschura, A. Matveev, C. Parthey, J. Alnis, R. Pohl, Th. Udem, N. Kolachevsky, and T. Hänsch, Phys. Rev. A 83, 042505 (2011).
  * Michel and Płoszajczak [2021] N. Michel and M. Płoszajczak,  _Gamow Shell Model: The Unified Theory of Nuclear Structure and Reactions_, Vol. 983 (Springer, 2021).
  * Gezerlis  _et al._ [2014] A. Gezerlis, I. Tews, E. Epelbaum, M. Freunek, S. Gandolfi, K. Hebeler, A. Nogga, and A. Schwenk, Phys. Rev. C 90, 054323 (2014).
  * Piarulli  _et al._ [2015] M. Piarulli, L. Girlanda, R. Schiavilla, R. N. Pérez, J. E. Amaro, and E. R. Arriola, Phys. Rev. C 91, 024003 (2015).
  * Piarulli  _et al._ [2016] M. Piarulli, L. Girlanda, R. Schiavilla, A. Kievsky, A. Lovato, L. E. Marcucci, S. C. Pieper, M. Viviani, and R. B. Wiringa, Phys. Rev. C 94, 054007 (2016).
  * Long and Yang [2012a] B. Long and C.-J. Yang, Phys. Rev. C 85, 034002 (2012a).
  * Long and Yang [2012b] B. Long and C.-J. Yang, Phys. Rev. C 86, 024001 (2012b).
  * Goto and Machida [1961] J. Goto and S. Machida, Prog. Theor. Phys. 25, 64 (1961).
  * Donoghue [1995] J. F. Donoghue, Dispersion relations and effective field theory (1995), arXiv: hep-ph/9607351.