{'fraction_non_alphanumeric': 0.04167391820267384, 'fraction_numerical': 0.023669013527410805, 'mean_word_length': 4.027921406411582, 'pattern_counts': {'":': 0, '<': 8, '<?xml version=': 0, '>': 20, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 9, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Subtractive renormalization of the NN interaction in chiral effective theory up to next-to-next-to-leading order: S waves

We extend our subtractive-renormalization method in order to evaluate the 1S0 and 3S1-3D1 NN scattering phase shifts up to next-to-next-to-leading order (NNLO) in chiral effective theory. We show that, if energy-dependent contact terms are employed in the NN potential, the 1S0 phase shift can be obtained by carrying out two subtractions on the Lippmann-Schwinger equation. These subtractions use knowledge of the the scattering length and the 1S0 phase shift at a specific energy to eliminate the low-energy constants in the contact interaction from the scattering equation. For the J=1 coupled channel, a similar renormalization can be achieved by three subtractions that employ knowledge of the 3S1 scattering length, the 3S1 phase shift at a specific energy and the 3S1-3D1 generalized scattering length. In both channels a similar method can be applied to a potential with momentum-dependent contact terms, except that in that case one of the subtractions must be replaced by a fit to one piece of experimental data. This method allows the use of arbitrarily high cutoffs in the Lippmann-Schwinger equation. We examine the NNLO S-wave phase shifts for cutoffs as large as 5 GeV and show that the presence of linear energy dependence in the NN potential creates spurious poles in the scattering amplitude. In consequence the results are in conflict with empirical data over appreciable portions of the considered cutoff range. We also identify problems with the use of cutoffs greater than 1 GeV when momentum-dependent contact interactions are employed. These problems are ameliorated, but not eliminated, by the use of spectral-function regularization for the two-pion exchange part of the NN potential


I. INTRODUCTION
Chiral perturbation theory (χPT) has now been applied to the problem of nucleon-nucleon (NN) interactions for almost two decades. Since Refs. [1,2], chiral potentials have been deduced from χPT Lagrangians [3,4,5,6,7,8] through the computation of all NN-irreducible diagrams that occur up to a given order in χPT. The resulting potential includes a variety of contact interactions, which parameterize our ignorance of physics coming from energies higher than the chiral-symmetry-breaking scale, Λ χ . Some of these contact interactions occur in the πN Lagrangian, and the low-energy constants (LECs) that multiply these operators can be determined independently from the pion-nucleon scattering data [9,10]. More of a challenge-both conceptual and practical-is posed by the string of operators that represent the short-distance physics in the NN sector. These operators renormalize the divergent loops found when calculating the χPT potential, and so encode the contribution to observables from high-energy NN states. Each such contact interaction is associated with an unknown LEC. χPT power counting applied to the NN potential mandates that at a fixed order in the chiral expansion only a finite number of LECs enter the potential. For instance, at O(P 3 ) in the NN S-waves (which will be where the bulk of our attention is in this work) there are five contact operators that must be considered: two act in the 1 S 0 channel, two in the 3 S 1 channel, and one causes coupling between the 3 S 1 and 3 D 1 partial waves.
For low partial waves the χPT potential is non-perturbative, and the standard approach involves use of the Lippmann-Schwinger equation (LSE) to reconstruct the full NN amplitude from the NN-irreducible part [11,12,13,14]. (C.f. Refs. [15,16,17,18,19,20,21,22,23,24] where parts of the potential are treated in perturbation theory.) Since the NN interactions obtained in χPT do not fall off as the nucleon momenta p go to infinity, a cutoff, denoted here by Λ, must be imposed on the intermediate states in the LSE. The LECs are then fitted to data for a variety of cutoffs, and should absorb any strong (power-law with power > 0 or logarithmic) dependence of observables on Λ. A corollary is that the predictions of the effective theory should not depend on the particular quantity chosen to fix the LECs, as long as the kinematic point is within the domain of validity of the chiral effective theory. A potential which does not have these desirable properties (approximate cutoff and renormalization-point independence) is not "properly renormalized". If the set of contact operators employed cannot achieve this we conclude that this "chiral effective theory" (χET) is unable to give reliable results, since its predictions are too sensitive to the treatment of the unknown short-distance physics at scale Λ. In Refs. [3,5,6,7,8] the potential V was computed to a fixed order, and then the NN LECs that appear in V were fitted to NN data for a range of cutoffs between 500 and 800 MeV. The resulting predictions-especially the ones obtained with the O(P 4 ) potential derived in Refs. [7,8]-are approximately Λ−independent and describe NN data with considerable accuracy.
However, the range of Λ considered in these papers is quite narrow. Does the fact that these analyses only consider Λ ≤ 800 MeV represent an intrinsic limitation on the χET approach?
In the case of leading-order potential (which is O(P 0 ) and whose long-range part consists only of one-pion exchange) the answer has been shown to be "No" in the 1 S 0 and and 3 S 1 − 3 D 1 channels. The problem can be properly renormalized for arbitrarily large Λ's with only constant contact interactions present in V in each channel-in accordance with the χPT power counting [24,25,26,27,28,29,30,31,32] 1 . However, once Λ > 600 MeV the χPT counting for the NN potential fails in waves with L > 0 where the tensor force from one-pion exchange is attractive [27,30,31,34,35]. That counting predicts no contact interaction in these channels, but a contact interaction is necessary to stabilize the LSE predictions at these cutoffs. Furthermore, in these partial waves it has been found that even when higher-order χPT potentials are considered and the phase shifts obtained from the LSE are approximately cutoff independent, the results can still show significant renormalization-point dependence [35].
At next-to-leading order (NLO, O(P 2 )) and next-to-next-to-leading order (NNLO, O(P 3 )), the χPT potential includes diagrams with two-pion-exchange (TPE). Since the loop integrals there diverge as the square or cube of the momentum transfer |q| ≡ |p ′ − p|, the TPE must be associated with contact terms up to order P 2 -just as one would have predicted using naive-dimensional analysis in powers of NN momentum.
In Refs. [5,36,37,45] the role of various different O(P 2 ) contact terms in the LSE was examined and the renormalization of the NN S-waves was discussed. Refs. [5,36,37] found that momentum-dependent contact terms had difficulty absorbing the cutoff dependence at larger cutoffs, and also that the contribution from the (supposedly higher-order) TPE was larger than one-pion exchange. (See Ref. [38] for a similar conclusion in the 1 S 0 channel using a co-ordinate-space analysis.) Ref. [36] also considered an energy-dependent O(P 2 ) contact interaction in the 1 S 0 channel and showed that the use of a contact interaction linear in energy (which formally should be equivalent to the momentum-dependent interaction at the orders in the chiral expansion being considered here) leads to the appearance of resonances there for Λ → ∞. In either case, the results of these studies suggest that for cutoffs larger than about 1 GeV the contribution from contact terms dominates over that from the long-range part of the potential, raising the question of whether χPT is acting as a phenomenological fit form, rather than an effective theory. As argued in Refs. [13,35,38], the chiral-symmetry-breaking scale could be responsible for this situation.
This suggests that there may be a critical cutoff ∼ Λ χ above which it does not make sense to iterate the TPE inside the LSE 2 . It is the purpose of this paper to critically examine the behavior of the NN phase shifts predicted by χET at NLO and NNLO in the 1 S 0 and 3 S 1 -3 D 1 channels over a range of cutoffs from Λ = 500 MeV to Λ > 2 GeV. In this way we hope to identify this critical cutoff, if it exists, and discuss the signatures of, and mechanisms responsible for, its appearance.
In order to address these issues we extend our previously developed subtractive-renormalization method [32,35] to evaluate the 1 S 0 and 3 S 1 − 3 D 1 phase shifts with the higher-order χET potentials. If energy-dependent contact interactions are employed in these potentials we show how to relate the χET phase shifts to on-shell quantities. For momentum-dependent contact interactions we formulate the scattering equations in terms of two (one) on-shell quantities in the triplet (singlet) channel and then fit the single remaining unknown parameter to experimental data. This allows us to easily go to high cutoffs in the LSE and provides clean information about the effect of the renormalization point. Our method works for any long-range potential. Here we consider long-range parts evaluated using both dimensional regularization (DR) and spectral-function-regularization (SFR) [37] at O(P 2 ) and O(P 3 ). The forms of these potentials were given in Ref. [35]. Since issues like differing extractions of the πN LECs c 1 , c 3 and c 4 could produce sizeable changes in the phase shifts we are not overly focused on the quality of our fit. Instead, our purpose is to see whether these potentials behave self-consistently with respect to cutoffs in the LSE. Thus, our strategy is to evaluate in detail the cutoff-dependence of the S-wave phase shifts, and examine the renormalizationpoint dependence. These features depend much less on the precise choice made for c 1 , c 3 and c 4 .
The structure of our work is as follows. First, in Sec. II, we introduce our subtractive-renormalization method and explain how to deal with both momentum-and energy-dependent contact terms in the singlet channel. We take as input to the method the 1 S 0 scattering length, a s , and the phase shift at a particular energy, E * . Then, in Sec. III, we extend our subtractive-renormalization method to the coupled-channels case and deal with constant, energy-dependent and momentum-dependent contact terms. While the constant contact term is, as in Ref. [32], solved by one subtraction with a t as input, for the energy-dependent contact term we pin down the three unknown constants with information with a t , the generalized scattering length α 20 that governs the low-energy behavior of ǫ 1 , and the phase shift at a particular energy as inputs. In Sec. IV, we discuss our results in the 1 S 0 channel and their implications. In Sec. V, we discuss our results in the 3 S 1 -3 D 1 channel and their implications. We summarize our findings in Sec. VI.

II. SUBTRACTIVE RENORMALIZATION IN THE SINGLET CHANNEL
The short-range part of the potential, i.e., the contact terms, is a parameterization of unknown high-energy physics. Before introducing our subtractive renormalization scheme, we need to decide what types of contact term we should consider. Throughout this work we will consider chiral two-pion-exchange [3,4,5,39] as the long-range part of the NN potential. This is associated with contact terms up to O(P 2 ). However, since for non-relativistic particles we have p 2 /M ∼ E there is the possibility to consider energy-dependent, instead of momentum-dependent contact interactions. In particular, this may have certain advantages in terms of evading theorems that limit the impact of short-distance potentials on phase shifts [41,42]. Thus, in the 1 S 0 channel we will consider the following three types of contact terms: where λ, γ, C 2 are unknown constants. (Of course, the numerical value of λ is different for each case.) The overall potential is then given by with the subscripts SR and LR standing for the short-and long-range parts of the potential. Specifically, in the 1 S 0 channel we have with expressions for the S-wave projected central, spin, and tensor isoscalar and isovector potentials taken from Ref. [4,37] and collected in Ref. [35].

A. The constant contact interaction
The constant contact term can be evaluated by one subtraction with the 1 S 0 scattering length a s as input [32]. This is possible because knowledge of the on-shell amplitude plus the long-range potential yields the fully off-shell amplitude, as long as only a constant contact interaction is present [32,43,44]. The step from the on-shell amplitude at an arbitrary energy E * , t(p * 0 , p * 0 ; E * ) (with p * 0 ≡ M E * ), to the half-off-shell amplitude at the same energy is explained in Appendix A. The symmetry property of the amplitude t(p, p * 0 ; E * ) = t(p * 0 , p; E * ) then means that the same sequence of manipulations can be used to construct t(p ′ , p; E * ). With the fully-off-shell amplitude at a single energy in hand, it is well known how to construct the t-matrix at any arbitrary energy [26]. In operator form we have where g 0 is the free resolvent of the Lippmann-Schwinger equation, (Here M is the nucleon mass and E + = E + iε, with ε a positive infinitesimal.)

B. Energy-dependent contact interaction
In this case there are two unknown constants. The strategy is to use one subtraction (as in Ref. [32]) to first eliminate the constant part of the contact term. For the energy-dependent contact term we then perform a second subtraction to eliminate the remaining unknown.
We start from the partial-wave LS equation, which is written explicitly in this channel as The energy-dependent 1 S 0 potential is given by Here v LR represents the long range part of the χPT potential, which is computed up to NLO or NNLO. However, the derivation presented below holds for any energy-independent long-range potential which is a function of p and p ′ and satisfies v LR (p ′ , p) = v LR (p, p ′ ). From now on we drop the partial-wave index, since this section only describes the singlet channel 1 S 0 .
To further simplify the presentation, we adopt the following operator notation of the LS equation Setting E = 0 in Eq. (2.7) leads to Equation (2.8) contains only one unknown, λ. Therefore, the matrix element t(p ′ , p; 0) can be obtained from one experimental point, i.e., the NN scattering length, a s , using the procedure explained in Subsection A, and given in detail in Ref. [32]. Using properties of the LSE, Eq. (2.8) can be rewritten as This allows to reconstruct λ if we wish to do so: Leaving E finite, Eq. (2.7) can be written as leading to a formal expression for the potential term Subtracting Eq.(2.10) from Eq.(2.12) isolates the term γE, Proceeding in a similar fashion as in Ref. [32], we can obtain t(E * ) at a fixed energy E * from the value of the NN phase shift at E * (for details see Appendix A). We find 3 (2.14) 3 Note that some care is required if the standard subtraction method is used to deal with the singularity when the integration momentum becomes equal to p 0 . In that case it is useful to introduce an additional term in the integral equation that involves the potential evaluated at the point p ′ = p 0 , p = p 0 . However, if other techniques are used (e.g. contour rotation) this additional term is not necessary and so we do not include it in our derivation.
Taking the difference of Eqs. (2.14) and (2.12) leads to To eliminate γ, we insert Eq. (2.13) into Eq. (2.15), Rearranging Eq. (2.16) and multiplying both sides with 1 + g 0 (E)t(E) from the right and 1 + t(0)g 0 (0) from the left, we arrive at can be calculated numerically. Since we already obtained the fully off-shell t(0) and t(E * ), Eq. (2.17) is an integral equation for t(E), which we can solve by standard methods.
To summarize, we perform two subtractions to the LS equation to eliminate the two unknown constants λ and γ. The resulting equation requires as input the scattering length a s and the phase shift at one specific energy E * . The only restriction on E * is that it must be within the domain of validity of our theory. Hence one can test the consistency of the theory by examining the extent to which results depend upon the choice of E * .

C. Momentum-dependent contact interaction
The coefficient of the momentum-dependent contact interaction, denoted here by C 2 , is not straightforwardly related to any S-matrix element. Thus we cannot apply our subtraction procedure in the case of a momentum-dependent contact interaction. Hence we instead adopt a "mixed" procedure, which involves a single subtraction plus fitting of C 2 . 4 This "mixed" procedure is carried out as follows: First we guess a value for the constant C 2 , and then insert it into the half-off-shell and on-shell LSE: Taking the difference of the two equations cancels the constant λ and leads to . (2.21) The above equation gives the value of λ that is consistent with the experimental value of the scattering length and our choice of C 2 . It thus defines a relationship λ = λ(C 2 ; a s ). (Note that-in spite of the form of Eq. (2.21)-the relationship is not linear, since t(p ′ , 0; 0) is also affected by the choice of C 2 .) Therefore, when trying to determine λ and C 2 we only need to guess C 2 and can then obtain λ from Eqs. (2.20) and (2.21). These two constants are then entered into the on-shell LS equation which is solved for the phase shifts. Finally C 2 is adjusted to fit the desired observable. In Sec. IV we will examine the results obtained when C 2 is adjusted to reproduce the 1 S 0 effective range, r 0 , and those found when we enforce the requirement that the theory correctly predict the phase shift at a particular energy E * .

III. SUBTRACTIVE RENORMALIZATION IN THE J=1 TRIPLET CHANNEL
Since the NN interaction is non-central, the S = 1 waves constitute a coupled-channel problem. Here we consider the 3 S 1 -3 D 1 coupled partial waves, and we again consider three different contact terms: Here we write the contact terms explicitly in their matrix form, where the diagonal represents the direct channels 3 S 1 -3 S 1 and 3 D 1 -3 D 1 and the off-diagonal the channels 3 S 1 -3 D 1 and 3 D 1 -3 S 1 . The unknown constants are λ, γ (or C 2 ) and λ t , and their value is different in each case. Case A is the leading-order contact interaction discussed in Refs. [24,32]. Cases B and C include the structures that appear at NLO [O(P 2 )] in the standard chiral counting for short-distance operators.
At this order in the chiral expansion the J=1 coupled-channel problem acquires two additional contact interactions. One is of the form σ 1 ·q σ 2 ·q [5], and hence gives a non-zero matrix element for the 3 S 1 -3 D 1 transition. The other can be written as either an energy-dependent or momentum-dependent piece of the 3 S 1 -3 S 1 potential, although until now most works on χET have considered only the latter [3,5,7,8], but see Ref. [36]. For the energy-dependent contact interaction we will show below that three subtractions can be performed to eliminate the unknown LECs. This leaves us with an integral equation for the t-matrix at an arbitrary energy that takes as input three experimental quantities. For case B, due to the momentum-dependence in the contact term, we use a double-subtraction-plus-one-fit scheme to solve the LSE. The procedures for solving case B and case C are quite similar, as they both employ subtractions to eliminate λ and λ t .

A. Constant contact interaction
This case can be solved with a single subtraction, as described in Ref. [32]. Once again, the key fact is that once a t , the triplet scattering length, is known, the value of t(0, 0; 0) is fixed, and that, together with knowledge of the long-range potential, is sufficient to obtain the two non-zero elements of t(p, 0; 0). The procedure by which we obtain these two elements, t 00 and t 20 , can be reconstructed from Appendix B by setting C 2 = 0 in the first part of the derivation presented there. Hermiticity of v then implies that t 00 (p, 0; 0) = t 00 (0, p; 0) and t 20 (p, 0; 0) = t 02 (0, p; 0), and the same set of manipulations can then be used to obtain t l ′ l (p ′ , p; 0). The operator equation (2.3), with E * = 0, then allows us to reconstruct t(E), albeit this time as a two-by-two matrix of operators, rather than the single operator t that is needed to describe the 1 S 0 channel.

B. Momentum-dependent central part
In this case we will eliminate two constants from the integral equation: λ, and the coefficient of the tensor shortdistance interaction λ t . We will then adjust the coefficient C 2 to reproduce one piece of 'experimental' information.
As usual, we begin by writing the on-shell and half-off-shell partial-wave LS equations: indicates the angular-momentum quantum number of the channels. To simplify notation, we write the 2 × 2 matrix as t. Subtracting Eq. (3.2) from Eq. (3.1) cancels λ, Here we used the threshold behavior of the partial-wave potential, v l ′ l LR (p, k) ∼ p l ′ k l , to infer that, e.g. v 20 LR (0, p) = 0. Equation (3.4) shows that this feature of v LR has as consequence that t 02 (p, 0; 0) = t 20 (0, p; 0) = 0 and t 22 (p, 0; 0) = t 22 (0, p; 0) = 0. Using Eq. (3.4) we can thus obtain the first part of the t-matrix equation that needs to be solved, namely where we have used that t 00 (0, 0; 0) = at M . We observe that λ t has been eliminated from this equation, but that we do not have a closed set of equations, since t 20 (p ′ , 0; 0) appears on the right-hand side.
In order to use the threshold behavior of t 20 we need to divide Eq. (3.3) by p 2 − p 2 0 and Eq. (3.6) by p 2 0 and then add both equations. A quick examination shows that γ and λ t both drop out of the final result-but the final result is quite lengthy. Hence, to simplify the notation arising from these manipulations we define the following map where the operator O can represent either the t-matrix or a potential function. With this notation we obtain The second row of the kernel of this equation simplifies considerably since , and extracting the S → D element of the resulting equation, we can also employ (again because of threshold behavior) together with the analogous result for v 22 LR , in order to obtain This is the key result of our derivation, since it provides a second integral equation that, together with Eq. (3.5), is a closed set for the half-off-shell t-matrix at zero energy. These equations require as input the long-range potential v LR , a choice of C 2 , and the experimental observables a t and the value of (3.14) This limit, and all the limits of v LR in Eq. (3.13), are well defined due to the threshold behavior of the potential and the t-matrix. The limit (3.14) is a specific case of a "generalized scattering length", which is defined as Since α l ′ l is related to an on-shell element of the t-matrix it can be obtained from an experimental phase-shift analysis (see, e.g. Ref. [52]).
Equations (3.5) and (3.13), together with an initial guess for the value of C 2 , now yield t 00 (p, 0; 0) and t 20 (p, 0; 0) by standard methods for the solution of integral equations. The other elements of the zero-energy t, t 02 (p, 0; 0) and t 22 (p, 0; 0), are equal to zero, as mentioned above.
It is formally possible to eliminate C 2 from Eq. (3.8) by a similar set of manipulations to the ones that led to Eq. (3.13). This yields a pair of coupled integral equations for t 00 (p, 0; 0) and t 20 (p, 0; 0) that is NN LEC free. However, these equations contain the quantity lim p0→0 F t 00 (0, p 0 ; E) as part of the driving term. This is a limit of a half-off-shell matrix element, and so cannot be obtained in a model-independent way from experimental data.
Therefore we now proceed similarly to the case of the momentum-dependent interaction in the singlet channel. We guess a value for C 2 , and solve Eqs. (3.5) and (3.13) for the half-off-shell t at zero energy. We then insert these quantities into the consistency condition for the LEC λ t : Once we have λ t , the consistency condition for λ reads The derivation of these two constraints for the LECs is given in Appendix B. Eqs. (3.15) and (3.16) define the full potential for the coupled-channels problem for a given trial value of C 2 . This potential is then inserted into the LSE and phase shifts computed. C 2 is adjusted until the desired experimental datum is reproduced with sufficient accuracy. This method has the advantage that Eqs. (3.16) and (3.15) define functional relationships: λ = λ(C 2 ; a t , α 20 ), and λ t = λ t (C 2 ; a t , α 20 ). (The generalized scattering length α 20 enters implicitly through its effect on the half-off-shell t-matrices appearing in the formulae for λ and λ t .) Thus, we need to perform only one-parameter fits in order to obtain all three of the LECs that are pertinent in the S-wave J=1 channel.

C. Energy-dependent central part
Here we could follow the method of Section II B and eliminate all three constants from the integral equation: λ, the coefficient of the tensor short-distance interaction, λ t , and the coefficient of the energy-dependent part, γ. However, it is simpler to mimic the steps of the previous Subsection, and reconstruct the underlying contact interaction from knowledge of the phase shifts. The manner in which this is done is described in detail in Appendix B but can be summarized in the following steps: 1. Use knowledge of a t , the triplet scattering length, to eliminate λ from the integral equation. Hence obtain t 00 (p, 0; 0).
2. Substitute the result into Eq. (3.4), with γ = 0 (since the energy-dependent piece of the contact interaction has no impact on t at E = 0) and rearrange the equation so as to obtain λ t from the on-shell quantity lim p0→0 t 20 (p 0 , p 0 ; E)/p 2 0 ≡ α 20 /M . 3. The value of λ can then be obtained from the requirement that t 00 (0, 0; 0) = a t /M . 4. Use the 3 S 1 and 3 D 1 phase shifts and the mixing parameter ǫ 1 at a specific value of the energy, say E * , to reconstruct t(p * 0 , p * 0 ; E * ). Then, from a suitably modified (3.3), solve for t(p, p * 0 ; E). 5. With this half-off-shell t at a finite energy in hand, γ can be obtained from the version of Eq. (3.2) that is suitable for the energy-dependent contact interaction. The final formula for γ is given in Eq. (B10).
It may appear that in this procedure one needs five pieces of input data to fix the constants (a t , the generalized scattering length α 20 , and the three phase shifts at energy E * ). However, only one of the four matrix elements of t is used to extract γ. Hence, in reality only three experimental inputs: t 00 (0, 0; 0), t 00 (p * 0 , p * 0 ; E * ), and the derivative of t 20 (p 0 , p 0 ; E) with respect to E at a specific energy, are needed to solve this problem.
We then obtain phase shifts for arbitrary E by by substituting the obtained values of λ, λ t and γ into the potential that appears in the LSE and solving that equation for the on-shell t-matrix.

IV. RESULTS: THE SINGLET S-WAVE
In this section we present the results obtained in the 1 S 0 channel when we employ the three different types of contact terms introduced in Sect. II. First, we only adopt the leading-order contact term, i.e. a constant. Then, we consider the energy-dependent contact term, i.e., v SR,0 = λ + γE. Finally, we use the more standard momentum-dependent contact term v SR,0 = λ + C 2 (p 2 + p ′2 ). For each of these cases, we examine the results found with the following three different forms of the long-range potential v LR : 1. TPE computed using dimensional regularization up to O(P 2 ) (denoted as DR NLO). 2. TPE computed using dimensional regularization up to O(P 3 ) (denoted as DR NNLO). 3. TPE computed using spectral-function regularization up to O(P 3 ) (denoted as SFR NNLO).
For the case of spectral-function regularization, we adoptΛ = 800 MeV as intrinsic cutoff. The effect that varying the intrinsic cutoff from 600 − 800 MeV has on this potential is discussed in reference [48]. When combined with the momentum-dependent contact interaction, the O(P 2 ) TPE yields the usual NLO calculation in χET, while the two O(P 3 ) forms of v LR both give calculations of NNLO accuracy.
A. The lowest-order contact term A proper renormalization of these TPE potentials entails the presence of contact terms up to second order in |q|. However, it is still of interest to investigate what happens if only the lowest-order contact term is included in the potential (see also Refs. [36,46]). This enables us to see how the phase shift converges as the LSE cutoff Λ is increased. Then, when higher-derivative contact interactions are added to the potential, we can see how the cutoff dependence changes.
The renormalization for the leading-order contact term is performed with one subtraction, with the scattering length a s = −23.7 fm [49] used as renormalization condition. This is the same as procedure as described for the 1 S 0 channel in Ref. [32], except that here we are using a long-range potential that is higher-order in the chiral power counting.
The phase shifts 1 S 0 calculated with the DR NLO, DR NNLO and SFR NNLO potential are shown for different cutoffs Λ in Fig. 1. The results show that a single, constant, contact interaction is sufficient to stabilize the 1 S 0 phase shift with respect to the cutoff once Λ > 800 MeV for the SFR NNLO TPE and Λ > 1000 MeV for the DR NLO TPE. The DR NNLO potential requires Λ > 1200 MeV to become stable.
The DR NLO potential exhibits a resonance-like behavior when Λ ≤ 700 MeV, but the phase shift quickly converges after Λ exceeds 800 MeV. This indicates that-at least for this potential-a cutoff smaller than 700 MeV removes essential physics, and brings the phase shift rather close to that obtained with v LR = v OP E in Ref. [32]. On the other hand, it is amusing to note that after δ( 1 S 0 ) converges with respect to the cutoff, the DR NLO result is actually closer to the PWA93 phase shifts [47,49] than is the result from the other two potentials.
This indicates that the pure (i.e. un-compensated by contact interaction) NNLO TPE produces too large a correction, although the general trend of providing more repulsion is correct. We note that the phase shift calculated with the SFR NNLO potential stabilizes as a lower value of Λ, and that this Λ → ∞ result is closer to the PWA93 result, compared to that obtained with the DR NNLO potential.
In the upper panel of Fig. 2 we plot the values of the effective range r 0 obtained from the above potentials. Since a s is fixed to the same value in all calculations, this figure shows the influence of Λ on the low-energy phase shifts. Comparing our results for the three different potentials to the "experimental" value r 0 = 2.7 fm [49], we see that the SFR NNLO potential gives a value of r 0 that is closest to the experimental value. The results with the DR NLO and SFR NNLO TPE are stable with respect to Λ once Λ ≥ 700 MeV. The DR NNLO TPE value eventually also stabilizes within small variations once Λ > 1200 MeV (not shown in Fig. 2). This cutoff is rather large, reflecting the too-strong energy dependence that was already discussed in the previous paragraph.

B. Second order: energy-dependent contact term
Next we use a second-order contact term with the simplest energy-dependent form, namely v SR,0 = λ + γE. This form arises because the contact operator N † i∂ t N N † N can absorb the divergences in the two-pion exchange potential, with the difference absorbed into terms that are proportional to the free-nucleon equation of motion. (See, e.g. Ref. [4], where the amplitude is calculated in Born approximation, and so the contact terms required are only a function of the nucleon on-shell momentum p 0 .) These additional terms can be transformed into higher-order multi-body operators. Thus, up to O(P 3 ) in χPT, the energy-dependent contact term is just as valid a representation of the short-distance physics as the more commonly chosen momentum-dependent one. Of course, the use of energy-dependent potentials engenders additional complications-especially when one wishes to embed the two-nucleon Hilbert space in a larger Fock space. Nevertheless, we can choose this form for the short-distance operators to describe NN scattering, and it has the advantage that the two unknown constants can be solved by our double-subtraction method discussed in Sect. II.
As input this method requires the value of a s and the phase shift at a specific energy E * . In principle, this energy can be chosen arbitrarily. In our calculation we choose T lab = 2E * = 2.8 MeV, i.e. at the energy where the 1 S 0 phase shift has its turning point. We tested other choices E * < 10 MeV and found that they have only a minor effect on the results. We also tested a relatively high energy, T lab > 100 MeV. Such a choice generates results which fit the phase shift at low energies (due to a s as input) and at the chosen higher energy (due to the use of δ(E * ) as input). However, we did not adopt this procedure, since in it the ultimate dependence of results on the choice of renormalization point (i.e. E * ) reflects the the cutoff-dependence in the phase shift in a manner that we shall discuss below.
In Fig. 3 the 1 S 0 phase shift calculated with the TPE potentials DR NLO, DR NNLO and SFR NNLO with the energy-dependent contact term is displayed for cutoffs Λ = 600 − 2000 MeV. A general feature is that the phase shift oscillates as a function of Λ. There are also strongly diverging phase shifts at Λ = 800 and 2000 MeV for the DR NLO potential because a resonance state is created by the contact term. (A similar behavior for the DR NLO potential was observed in Ref. [36] for a momentum-dependent contact term with Λ = 2000 MeV.) For the DR NNLO (SFR NNLO) TPE, at Λ = 1000(2000) MeV the phase shift diverges violently from the previous converged value.
In order to further understand the effect of the energy-dependent contact term, we calculate the DR NNLO TPE phase shift at T lab = 10 and 100 MeV for Λ up to 19 GeV and display the results in Fig. 4. The phase shifts have the same oscillatory behavior, with the period becoming longer as Λ becomes larger. Comparison with the results obtained by adopting only a constant contact term (second panel of Fig. 1, where this phase shift converges with respect to Λ for Λ > 1200 MeV), suggests that the quasi-periodic behavior of δ(T lab = 10 & 100 MeV) is caused by the energy-dependent contact term.
The oscillatory behavior of the phase shift at this energy is associated with resonances at specific Λ. An argand plot confirms that there is a resonance for the cutoff at which this phase shift diverges. Tracking the energy of the resonance as a function of Λ shows that it moves in a regular way, with the real and imaginary parts of E R decreasing in size as Λ increases, until eventually the resonance becomes a bound state and the pole in t moves towards an energy of −∞. In order to understand this better, we perform a two-potential analysis of the t-matrix resulting from this particular v SR . The details are given in Appendix C. We find that the appearance of poles is a consequence of the specific form of the contact interaction. When these poles are in the domain of validity of the theory the resulting phase shifts are in poor agreement with the Nijmegen PSA. Only at those cutoffs where the artificial poles correspond to deeply bound states or high-energy resonances do the phase shifts match the Nijmegen analysis well.
The effective range r 0 obtained from all three different potentials is about r 0 = 2.68fm and varies by < 2% over the range of cutoffs considered, with the variation being only ∼ 0.5% if the long-range part is computed to O(P 3 ). This is a consequence of our choice of a low-energy phase shift as second input. Comparing with the upper panel of Fig. 2, one can clearly see that r 0 receives significant correction from the O(P 2 ) contact term, especially for the DR NLO potential.
Thus, the energy-dependent contact term improves the overall agreement of the 1 S 0 phase shift, and in particular the effective range r 0 , with respect to the average value r 0 = 2.7 fm extracted from the different Nijmegen potentials. Once one tries to increase the cutoff Λ beyond 1000 MeV, this form of short-distance physics leads to problems if the DR TPE is used as long-range potential. However, adopting TPE computed using spectral-function regularization means a wider range of cutoffs can be employed before the same issues occur as for the DR TPE.
It is worthwhile mentioning that there exists a particular regularization of two-pion exchange which avoids resonant behavior in the phase shifts. This involves using dimensional regularization for the NLO part of TPE, and spectralfunction regularization for the NNLO part 5 . Then, even if an energy-dependent contact term is chosen for the short-range potential, only minor oscillations of the phase shifts with respect to Λ occur for Λ > 1 GeV. However, strictly speaking, a potential which is obtained using two different regularization schemes is not self-consistent. Thus, in the later discussion we will only present the results obtained from the three TPE potentials mentioned at the beginning of this section.

C. Momentum-dependent contact term
A momentum-dependent contact term associated with the TPE has the form λ + C 2 (p 2 + p ′2 ), and can be solved by one subtraction plus the use of a one-parameter fit, as explained in Sect. II. The input chosen for the subtraction is again a s . In what follows we will try to perform two different fits to fix the value of C 2 . One uses only low-energy information, namely the effective range, r 0 = 2.7 fm. The other fit fixes C 2 by attempting to reproduce the phase shift at T lab = 200 MeV.
The resulting 1 S 0 phase shift is shown in Fig. 5 for the DR NLO TPE. In this case we cannot reproduce r 0 , or the phase shift at T lab = 200 MeV, once Λ > 600 MeV. This pathology of the fit at NLO was discussed in Ref. [46], and is related to the Wigner bound that limits the impact of short-distance physics on the effective range-or more generally on the energy dependence of phase shifts [42,50]. Therefore, for the DR NLO TPE we can only perform an overall best fit to the phase shift. Fig. 5 shows that, within this limitation, the phase shift quickly becomes cutoff independent once Λ > 900 MeV.
The results for the DR NNLO TPE are shown in Fig. 6, where we see that the two different fit procedures generate different results for the same Λ. This is especially visible at Λ = 500 and 1000 MeV, where a resonance-like behavior is present in the latter case when C 2 is fitted to r 0 . For values of Λ not close to these problematic cutoffs the phase shift is almost independent of the renormalization point.
This independence of renormalization point is even more apparent when the SFR NNLO TPE is adopted. As can be seen in Fig. 7, for Λ between 700 − 1800 MeV, the two different fitting procedures lead to almost the same 1 S 0 phase shift. The same figure also suggests that Λ < 600 MeV is too small and cuts off too much of the potential to allow a good fit to the Nijmegen analysis. Finally, the agreement between the two fits breaks down at Λ = 2000 MeV. Thus, by adopting the SFR NNLO TPE, we achieve renormalization-point independence for a wider range of cutoffs, but this renormalization-point independence breaks down for higher Λ.
Finally, we show in the lower panel of Fig. 2 the extracted effective range r 0 as function of Λ for the 1 S 0 phase shift fitted to the Nijmegen value at T lab = 200 MeV. This allows to assess the sensitivity of the results to renormalization point (E * ≈ 0 or E * ≈ 200 MeV lab. energy). The DR NLO TPE has the least variation of r 0 with Λ over the range of Λ considered here, but it needs to be remembered that the experimental value cannot be reproduced with a real value of C 2 [46,51]. At order NNLO, the SFR TPE gives values of r 0 which vary less with Λ than do the ones from the DR TPE once Λ > 800 MeV, and are stable until Λ = 2000 MeV. This can already be inferred from Fig. 7.
In summary, the results obtained from the NNLO DR TPE with a momentum-dependent contact term have a peculiar behavior at Λ = 1000 MeV. Adopting the SFR potential increases the validity range of Λ as in the energydependent case.

V. RESULTS IN THE J=1 TRIPLET CHANNEL
In the J=1 triplet channel we again adopt the three different potentials, i.e., the DR TPE up to NLO, the DR TPE up to NNLO and the SFR TPE up to NNLO. To renormalize the potential in this channel we employ three different contact terms (cases A-C in Sec. III). For case B (energy-dependent) and C (momentum-dependent), the value of the generalized scattering length is required. We first adopt the value given in [52] for α 20 . Then we adjust α 20 to find the best fit. Since α 20 is not given from experiment, we list the values extracted from several so called "high-precision" potentials in Table I.

A. The lowest-order contact term
We first perform the renormalization with the constant contact term (A) in the 3 S 1 − 3 S 1 channel. The results obtained with this short-distance potential, and long-distance potentials DR NLO, DR NNLO and SFR NNLO are shown in Figs. 8-10, respectively. In the first two cases the 3 S 1 phase shifts reach cutoff independence as Λ approaches 1000 MeV, and then diverge once more at higher Λ. The 3 D 1 phase shifts exhibit a similar feature, although the pattern of convergence is not as obvious as in 3 S 1 . For the mixing angle ǫ 1 , there is no obvious convergence with respect to Λ for the DR NLO and DR NNLO TPE. This can be explained by the fact that we only adopt a contact term in the 3 S 1 − 3 S 1 channel and fit to the scattering length a t = 5.428 fm. This was sufficient to stabilize the J=1 triplet channel in the case of the one-pion-exchange potential [24,28,32]. However, the NNLO TPE is singular and attractive in both eigen-channels and thus needs two contact interactions (or, equivalently, two subtractions) in order to stabilize its predictions with respect to Λ [46]. In contrast, the NLO TPE has a repulsive behavior as r → 0 in both J=1 eigen-channel. This also leads to predictions that are unstable after one subtraction, albeit for the opposite reason than at NNLO. The SFR NNLO potentials do not have the same singular attractive behavior at short distances, and so it is no coincidence that its phase shifts show the best overall convergence in 3 S 1 , 3 D 1 , and ǫ 1 . The SFR strongly reduces the short-range attraction in the NNLO pieces of the TPE, and a lower-order contact term is then able to remove the Λ dependence of the phase shifts. Fig. 11 shows that the phase shifts at T lab = 10 and 50 MeV obtained with the SFR NNLO TPE and a constant contact term become stable for Λ > 2500 MeV.

B. Energy-dependent contact term
The results of using an energy-dependent contact term to renormalize the three different long-range potentials are shown in Figs. 12-15. In Fig. 12 we adopt a t = 5.428 fm, α 20 = 2.28 × 10 −10 MeV −3 [52], and the Nijmegen PWA values of δ( 3 S 1 ), δ( 3 D 1 ) and ǫ 1 at T lab = 10 MeV as our input to renormalize the DR NNLO TPE. This choice does not produce an optimal description of ǫ 1 with respect to the Nijmegen analysis, although it is pretty close.
Actually, ǫ 1 is a notoriously delicate observable, and the choice made for α 20 affects it appreciably. Therefore, in the case of DR NNLO we further adjust α 20 to 2.25 × 10 −10 MeV −3 , which yields the best fit of ǫ 1 with respect to the Nijmegen values for Λ = 700 − 1200 MeV. This result is shown in Fig. 13. Comparison with Fig. 12 shows that this 1.5% shift in α 20 leads to noticeable improvement in ǫ 1 over this cutoff range, with few other significant changes. The same value, α 20 = 2.25 × 10 −10 MeV −3 , also gives the best fit for DR NLO, and this α 20 produces results for SFR NNLO which coincide quite well with the Nijmegen data. In addition, α 20 = 2.25 × 10 −10 MeV −3 is still within a few per cent of the value obtained for this quantity from the potentials constructed in Ref. [53]. Such a difference is certainly within the range of expected N 3 LO corrections to our results.
The DR NNLO results present several interesting features. First, there exists a range of cutoffs Λ = 800 − 900 MeV where all three phase shifts agree with the Nijmegen analysis remarkably well 6 . However, with higher cutoffs the phase shifts gradually move away from the Nijmegen analysis. Second, for Λ = 2000(1400) MeV, δ( 3 S 1 )(δ( 3 D 1 )) diverge and show a resonance-like behavior.
In order to analyze this feature, we calculate the J=1 triplet phase shifts for a fixed energy of T lab = 10 and 50 MeV and a wide range of cutoffs Λ = 500 − 5500 MeV, and plot them in Figs. 16 and 17. Both figures show oscillatory behavior of the phase shifts with respect to Λ. The S and D-waves have different points of divergence in this cycle, and (unsurprisingly) affect one another when either gets large. The oscillation pattern of the mixing angle ǫ 1 is a combination of that appearing in the two phase shifts δ( 3 S 1 ) and δ( 3 D 1 ).
In Fig. 15 we show the results for the SFR NNLO potential with α 20 = 2.25 × 10 −10 MeV −3 . The phase shifts only show minor cutoff dependence for Λ < 1000 MeV and converge nicely for higher Λ. As far as the fit to the Nijmegen PWA is concerned, a decent agreement is observed for all Λ for δ( 3 S 1 ) and δ( 3 D 1 ). For ǫ 1 , the converged curves do not fit as well as those for Λ = 600 − 800 MeV. This is not a surprise, since we adopted a value of α 20 which fits ǫ 1 in this region of cutoff. We found that at higher cutoffs, i.e., Λ ∼ 2300 MeV, the SFR NNLO TPE presents a divergent pattern in the phase shifts. Similar to the singlet channel, SFR NNLO can defer, but cannot avoid, the resonance caused by the energy-dependent contact term.
The effective range r 0 obtained at NLO is r 0 = 1.75 ± 0.01 fm for 900 MeV ≤ Λ ≤ 1800 MeV, with the error representing the extent of the cutoff dependence in that range. For Λ = 2000 MeV, r 0 = 1.70fm, which is due to the resonance-like behavior in the phase shifts created by the energy-dependent contact term. For DR NNLO, r 0 = 1.75 ± 0.01 fm except for Λ = 1400 and 2000 MeV, which are just the cutoffs where the phase shifts also diverge, as shown in Figs. 16 and 17. The same behavior in r 0 is present in the SFR NNLO, but for higher cutoffs.
In summary, in the triplet channel the energy-dependent contact term creates a similar pattern as in the singlet channel. For the DR TPE (SFR NNLO), there is a highest cutoff Λ = 1200(2300) MeV with which the potential can be iterated in the LSE before resonances are generated at energies within the domain of validity of χET.

C. Momentum-dependent contact term
The results of using a momentum-dependent contact term to renormalize the DR NNLO TPE are shown in Fig. 18. Here we adopt a t = 5.428 fm and α 20 = 2.25 × 10 −10 MeV −3 as input and then perform a fit to the phase shift δ( 3 S 1 ) at T lab = 200 MeV. This pins down the three unknown LECs in the potential at this order.
For phase shifts obtained by our two-subtractions-plus-one-fitting method, once the trial value of C 2 is imposed, the other two constants λ and λ t are fixed by the two pieces of experimental information that are our input: α 20 and a t . Thus, the results shown in Fig. 18 satisfy the constraints at T lab = 0 associated with this information and reproduce the 3 S 1 phase shift at T lab = 200 MeV. However, they do not provide the best possible fit to the phase shifts. In contradistinction to the previous section, the fit is not appreciably improved by variations of α 20 . Compared to the results obtained using the same contact terms in Refs. [5,7,8], the overall fit to the Nijmegen phase shifts is noticeably worse. Furthermore, the results obtained by fixing δ( 3 S 1 ) at T lab = 200 MeV do not show stability for δ( 3 D 1 ) and ǫ 1 with respect to Λ. This suggests that the TPE associated with a momentum-dependent contact term depends strongly on the choice of renormalization point once Λ > 1 GeV.
We wish to see up to what value of Λ it is possible to ameliorate these difficulties by a different choice of input. Thus we remove the two constraints imposed by a t = 5.428 fm and α 20 = 2.25 × 10 −10 MeV −3 and allow all three NN LECs to vary freely, so that we can obtain the best overall fit to the phase shifts. The results for DR NLO, DR NNLO and SFR NNLO are listed in Figs. 19 and 20, respectively. The fits are a noticeable improvement over Fig. 18. For all three potentials, starting from Λ = 700 MeV, the 3 S 1 and 3 D 1 phase shifts fit the Nijmegen analysis quite well and have only minor variation with respect to Λ up to 1000 MeV. A stronger cutoff dependence shows up only in the mixing parameter ǫ 1 . There we observe a clear decrease in the Λ sensitivity of the result as the long-range potential is improved from DR NLO to DR NNLO and then again to SFR NNLO. We plot r 0 extracted from those best fits in Fig. 21. Again, we see better agreement with the Nijmegen value r 0 = 1.833 fm [52] for DR NNLO than for DR NLO. For 700 ≤ Λ ≤ 1000 MeV, r 0 varies by less than 3% for all three potentials.

VI. CONCLUSION
In this paper we developed subtractive renormalization schemes which can be applied to NN scattering in χET at NNLO and maximize the use of information from experiment. We showed how to manipulate LSEs containing the contact interactions of χET so as to eliminate constant, tensor, and energy-dependent contact operators from the equation in favor of pieces of experimental data. A by-product of our analysis is the result that the momentumdependent contact operator cannot be eliminated in this way, and thus is not straightforwardly related to any NN S-matrix element. We can still, however, analyze the χET potential containing that operator, by using a mixture of subtractive renormalization and more standard fitting methods.
The subtractive-renormalization method shows particular advantages in cases where the LECs undergo limit-cycle or limit-cycle-like behavior, and thus adjusting several of them to reproduce data represents a difficult fitting problem. It is also very useful when carrying out calculations at cutoffs Λ > 1 GeV, since fine-tuning of the unknown NN LECs occurs in that regime. Finally, the direct input of experimental information allows a straightforward determination of the impact of the renormalization point-i.e. the choice of which experimental information we input-on the resulting phase shifts.
Looking first at the case of a constant contact interaction, in the 1 S 0 channel our results confirm the finding of Ref. [36]: when a constant contact term is adopted, the resulting 1 S 0 phase shift converges and becomes independent of Λ in the limit Λ → ∞. This occurs once Λ > 1000(1200) for the DR NLO and SFR NNLO (DR NNLO) potentials. However, we found that in this calculation the contribution from the NNLO(Q 3 ) part of the TPE is larger than the NLO part.
In contrast, in the coupled 3 S 1 − 3 D 1 channels the dimensionally regularized interactions cannot be renormalized with a single, constant, contact interaction. However, when the two-pion exchange is computed using spectral-function regularization Λ-independent results are again obtained at large Λ: presumably because the short-distance behavior of this potential is the same as that of one-pion exchange, which is well-known to be renormalized by a single constant contact term in the 3 S 1 − 3 S 1 part of the potential [24,28,32].
We also considered the contact terms that can occur up to O(Q 2 ) in χET, examining both energy-dependent and momentum-dependent contact terms. For both the 1 S 0 and 3 S 1 − 3 D 1 channels, the phase shifts obtained from dimensionally regularized two-pion exchange show oscillatory behavior as a function of Λ when energy-dependent contact terms are used. This is associated with the appearance, and subsequent movement with energy, of poles in the NN amplitude. The first cutoff at which such a pole appears in the domain of validity of χET is around Λ = 1000 MeV.
We also found that cutoffs less than about 600 MeV sometimes produced resonant features in the phase shifts. We attribute this failure to the propensity of these cutoffs to eliminate too much of the central attraction in the TPE that is responsible for reproducing the energy dependence of the S-wave phase shifts extracted from experiment-especially in the 1 S 0 .
In between these two limits, i.e. for 600 MeV < Λ < 1 GeV, we found that the 1993 Nijmegen PWA phase shifts, in both the 1 S 0 and 3 S 1 − 3 D 1 channels, could be quite well reproduced using a short-distance potential that includes an energy-dependent contact term, together with a long-distance potential calculated to O(Q 3 ) in χPT. With an energy-dependent contact term good fits to the Nijmegen phase shifts can be obtained at cutoffs above 1 GeV, but only in fairly narrow windows of Λ, where the energy of the spurious bound state produced by the energy-dependent contact term is large compared to the energy scales of interest.
Finally, if the spectral-function regularized TPE is adopted, then the range of the cutoff one can use before the resonance-like behavior occurs is greatly increased.
Turning to the case of the momentum-dependent contact interaction: we found that in the case of this type of contact interaction there is a significant improvement with respect to the Nijmegen phase shift when one improves the long-range part of the potential from NLO to NNLO. Indeed, if the χET NLO potential is employed there are certain observables (e.g. the effective range) that cannot be reproduced if we restrict ourself to real values of the LECs in the NN potential. If we wish our NLO calculation to reproduce the effective range in the 1 S 0 we are restricted to Λ < 600 MeV. This problem is absent at NNLO.
We examined the renormalization-point dependence of our results, and found that in the 3 S 1 − 3 D 1 channel at Λ = 800 MeV the results from both the DR and SFR TPE up to NNLO are almost independent of the renormalization point. However, once Λ = 1000 MeV is reached, the results with the DR NNLO potential show appreciable dependence on the choice of experimental input. Once again, the SFR TPE is better behaved in this regard: there is only minor renormalization-point dependence for the SFR NNLO potential from Λ = 600 MeV until Λ = 1800 MeV, and the renormalization-point independence breaks down at Λ = 2000 MeV.
The use of a momentum-dependent potential and the expressions for the NNLO TPE should allow us to confirm previous results for the triplet channel (e.g. those of Ref. [5]). However, we found that if we impose strict relations between the three unknown LECs in order to reproduce particular pieces of experimental data (e.g. a t and α 20 ), the resulting fit to the Nijmegen phase shifts is in general worse than that obtained by direct fitting of data over a finite range of energies [5]. Only when all the three constants are allowed to vary independently do we get best fits, which are comparable in quality to the results in [5,46]. However, we want be emphasize that this fit can only be obtained for cutoffs Λ = 700-900 MeV for the DR NNLO TPE. (For the SFR NNLO TPE, it is possible to obtain a reasonable fit up to Λ = 2000 MeV.) There is thus an appreciable difference between using two subtractions that fix the zero-energy behavior of the amplitude and performing a best fit involving all three J = 1, S = 1 NN LECs. That difference is an indication that the triplet S-wave phase shifts obtained from chiral TPE depend on the renormalization point rather strongly, despite the fact that they are quite cutoff independent in this window of Λ.
In summary, we have probed the renormalization-point and cutoff dependence of chiral potentials which include TPE when these potentials are iterated using the LSE. We find that energy-dependent contact terms-which might be considered useful as they allow one to evade certain theorems regarding short-distance potentials [41,42]-tend to introduce unphysical resonances in the results once Λ > 1 GeV. These resonances are generically absent when a momentum-dependent contact term is employed. Since the χET power counting suggests that energy-and momentum dependent contact interactions are equivalent up to the order to which we work this implies that the power counting is only a useful guide for Λ < 1 GeV.
Although the momentum-dependent contact interaction does not produce the dramatic failures to agree with experiment that occur in the energy-dependent case, it does have some difficulties. Chief among these is the failure of certain aspects of the NLO calculation. Serious doubt must exist regarding the convergence of an expansion where (at least for the 1 S 0 channel) LO is in poor agreement with experiment, NLO cannot reproduce the effective range, and it is not until NNLO that anything like a reasonable picture of the phase shifts emerges.
The less singular interactions obtained with spectral-function regularization avoid many of the difficulties inherent in the use of their dimensionally regularized counterparts. The cutoff dependence is almost always weaker, with the cutoff dependence of phase shifts found with the SFR NNLO potential and a constant contact term disappearing altogether for Λ > 1 GeV. However, ǫ 1 cannot be reproduced in this range of cutoffs, and the issue with the effective range at NLO persists.
We conclude that it makes little sense to iterate the two-pion-exchange potentials obtained from chiral perturbation theory using the Lippmann-Schwinger equation unless one restricts the momentum in the LSE integral to lie below ∼ 1 GeV [11,13,14]. In this domain, and for the case of a momentum-dependent contact term, our subtractiverenormalization method produces results which agree with previous analyses [5,34,36,37,46].

Acknowledgments
This work was performed in part under the auspices of the U. S. Department of Energy, Office of Nuclear Physics, under contract No. DE-FG02-93ER40756 with Ohio University. We thank the Ohio Supercomputer Center (OSC) for the use of their facilities under grant PHS206. D.R.P. is grateful for financial support form the Mercator programme of the Deutsche Forschungsgemeinschaft, for the hospitality of the Theoretical Physics group at the University of Manchester during the initial stages of this work, and for useful conversations with Evgeny Epelbaum.
For any given energy E * the on-shell value of the t matrix is related to the phase shift δ(E * ) by where f 00 (E * ) represents the scattering amplitude and p * 0 = √ M E * . Then we write the half-shell and on-shell LS equations as The above equation does not contain any unknown quantities and can solved with standard techniques. With t(p * ′ , p * 0 ; E * ) at hand, we can carry out another subtraction to get t(p * ′ , k * ; E * ), where k * represents any arbitrary momentum.

APPENDIX B: DETAILS OF THE SUBTRACTION METHOD IN THE J=1 TRIPLET CHANNEL
Let us start with the on-shell and half-off-shell LSE in this, two-channel, case. We will consider a momentumdependent short-distance potential, so we have: Subtracting Eq. (B1) from Eq. (B2) gives The only unknown in Eq. (B3) is λ t , thus we can express Similarly, once λ t is calculated, from Eq. (B1) we can express λ as Thus, by using the two scattering length a t and α 20 together with a trial value for C 2 , we can obtain the other two unknown constants λ and λ t . This determines the complete NN potential at this order, and it is then a straightforward matter to compute the phase shifts using the LSE. The LSE with the energy-dependent contact term (case B) satisfies exactly the same relation as Eq. (3.8) and Eq. (3.13), while Eq. (3.5) is replaced by Equation (B6), together with Eq. (3.8), can be used to solve for t(p, 0; 0). Note that here we don't need to input a trial value for C 2 . Once we have t 00 (p, 0; 0) at hand, we can substitute it, together with the t 20 obtained from Eq. (3.13), back into Eq. (3.4) (which is unaffected by the issue of energy-vs. momentum-dependence and so is the same as Eq. (B4)) and obtain λ t . The value of λ then can be obtained from the new version of Eq. (B1), Finally, in order to solve for γ, one needs a third piece of information as input. A convenient choice is the phase shifts at an energy E. From 3 S 1 , 3 D 1 and ǫ 1 at a specific energy E we compute t(p 0 , p 0 ; E). Then, from a modified Eq. (3.3), we solve for t(p, p 0 ; E). Once we have t(p, p 0 ; E) at hand, γ can be obtained from a modified Eq. (3.2): as: Note that here we have chosen to solve γ by matching the l = l ′ = 0 part of the matrix element in both sides of Eq. (B9). One also could have chosen to solve γ by matching the ll ′ = 02 or the ll ′ = 22 part of the equation. This corresponds to the fact that we have three phase shifts in the coupled triplet, namely 3 S 1 , 3 D 1 and ε 1 , however we only need one additional condition to determine the value of the third unknown γ.

APPENDIX C: TWO POTENTIAL ANALYSIS OF THE TRANSITION MATRIX AND ITS RELATION TO BOUND STATES
Here we want to consider the TPE potential with an energy dependent contact term, i.e., where V 2 should contain the energy dependence. The general derivation of solution of the LS equation based on the Gellman-Goldberger relation [55] leads to the following expression for the t-matrix Eq. (C2) can be solved by iteration, and the first order gives the well known DWBA expression. However, we do not want to follow along that line. Instead we will use Eq. (C3) to derive a series representation of t in terms of V 2 and easily calculable integrals. Explicitly the result up to first order in V 2 , which we denote by t (1) , is given as Since V 2 = γE does not depend on the momenta p and p′ it can always be taken out of the integral. Therefore, we define and observe that we may always make the replacement t 1 G 0 → Γ and G 0 t 1 → Γ if the operators in question appear next to a V 2 . Using this notation we can rewrite Eq. (C4) as In the same notation, the sum up to second order in V 2 , t (2) , will be Let us further define and note that, when evaluated between two V 2 's, g 0 t 1 g 0 may be replaced by f (E). Consequently t (1) and t (2) can be simplified as The next order, t (3) reads Continuing, the N th order t (N ) can be written as Letting N → ∞, we obtain the formal sum This tells us that t−matrix diverges when 1 V2 = f (E) + G(E). We have performed such a calculation for the bound states discussed in Sec. IV B and found that Eq. (C12) indeed holds at the energies where they occur.            Here r0 is shown as a function of the cutoff Λ in the LSE, and is extracted from a best fit of the phase shifts to the Nijmegen PWA93. The solid line represents r0 obtained from Ref. [52]. Note that the value of r0 at Λ = 600 MeV for DR NNLO is −11.5 fm, which is not plotted in the figure.