Effective field theory in the harmonic oscillator basis

# Effective field theory in the harmonic oscillator basis

S. Binder Department of Physics and Astronomy, University of Tennessee, Knoxville, TN 37996, USA Physics Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA    A. Ekström Department of Physics and Astronomy, University of Tennessee, Knoxville, TN 37996, USA Physics Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA    G. Hagen Physics Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA Department of Physics and Astronomy, University of Tennessee, Knoxville, TN 37996, USA    T. Papenbrock Department of Physics and Astronomy, University of Tennessee, Knoxville, TN 37996, USA Physics Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA    K. A. Wendt Department of Physics and Astronomy, University of Tennessee, Knoxville, TN 37996, USA Physics Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA
###### Abstract

We develop interactions from chiral effective field theory (EFT) that are tailored to the harmonic oscillator basis. As a consequence, ultraviolet convergence with respect to the model space is implemented by construction and infrared convergence can be achieved by enlarging the model space for the kinetic energy. In oscillator EFT, matrix elements of EFTs formulated for continuous momenta are evaluated at the discrete momenta that stem from the diagonalization of the kinetic energy in the finite oscillator space. By fitting to realistic phase shifts and deuteron data we construct an effective interaction from chiral EFT at next-to-leading order. Many-body coupled-cluster calculations of nuclei up to Sn converge fast for the ground-state energies and radii in feasible model spaces.

###### pacs:
21.30.-x, 21.30.Fe, 21.10.Dr, 21.60.-n

## I Introduction

The harmonic oscillator basis is advantageous in nuclear-structure theory because it retains all symmetries of the atomic nucleus and provides an approximate mean-field related to the nuclear shell model. However, interactions from chiral effective field theory (EFT) Epelbaum et al. (2009); Machleidt and Entem (2011) are typically formulated in momentum space, while the oscillator basis treats momenta and coordinates on an equal footing, thereby mixing long- and short-ranged physics. This incommensurability between the two bases is not only of academic concern but also makes oscillator-based ab initio calculations numerically expensive. Indeed, the oscillator basis must be large enough to accommodate the nucleus in position space as well as to contain the high-momentum contributions of the employed interaction. Furthermore, one needs to perform computations at different values of the oscillator spacing to gauge model-space independence of the computed results Maris et al. (2009); Hagen et al. (2010); Jurgenson et al. (2013); Roth et al. (2014). Several methods have been proposed to alleviate these problems. Renormalization group transformations, for instance, are routinely used to “soften” interactions Bogner et al. (2003, 2007); Roth et al. (2014), and many insights have been gained through these transformations Jurgenson et al. (2009). However, such transformations of the Hamiltonian and observables Lisetskiy et al. (2009); Schuster et al. (2014) add one layer of complexity to computations of nuclei.

One can contrast the effort of computations in the oscillator basis to, for instance, computations in nuclear lattice EFT Lee (2009). Here, the effective interaction is tailored to the lattice spacing and, thus, to the ultraviolet (UV) cutoff, and well-known extrapolation formulas Lüscher (1986) can be used to estimate corrections due to finite lattice sizes. The lattice spacing is fixed once, reducing the computational expenses. This motivates us to seek a similarly efficient approach for the oscillator basis, i.e., to formulate an EFT for nuclear interactions directly in the oscillator basis.

In recent years, realistic ab initio nuclear computations pushed the frontier from light -shell nuclei Navrátil et al. (2009); Barrett et al. (2013) to the medium-mass regime Hagen et al. (2012); Holt et al. (2012); Wienholtz et al. (2013); Somà et al. (2014); Lähde et al. (2014); Hagen et al. (2014); Hergert et al. (2014); Hagen et al. (2016). At present, the precision of computational methods considerably exceeds the accuracy of available interactions Binder et al. (2014), and this is the main limitation in pushing the frontier of ab initio computations to heavy nuclei. To address this situation, several efforts, ranging from new optimization protocols for chiral interactions Ekström et al. (2013, 2015); Carlsson et al. (2016) to the inclusion of higher orders Entem et al. (2015); Epelbaum et al. (2015) to the development of interactions with novel regulators Gezerlis et al. (2014) are under way. To facilitate the computation of heavy nuclei we propose to tailor interactions from chiral EFT to the oscillator basis.

There exist several proposals to formulate EFTs in the oscillator basis. Haxton and coworkers proposed the oscillator based effective theory (HOBET) Haxton and Song (2000); Haxton and Luu (2002); Haxton (2007, 2008). They focused on decoupling low- and high-energy modes in the oscillator basis via the Bloch-Horowitz formalism, and on the resummation of the kinetic energy to improve the asymptotics of bound-state wave functions in configuration space. The HOBET interaction is based on a contact-gradient expansion, and the matrix elements are computed in the oscillator basis. The resulting interaction exhibits a weak energy dependence. The Arizona group Stetcu et al. (2007, 2010); Rotureau et al. (2012) posed and studied questions related to the UV and infrared (IR) cutoffs imposed by the oscillator basis, and developed a pion-less EFT in the oscillator basis. This EFT was also applied to harmonically trapped atoms Rotureau et al. (2010). In this approach, the interaction matrix elements are also based on a contact-gradient expansion and computed in the oscillator basis. Running coupling constants depend on the UV cutoff of the employed oscillator basis. In a sequence of other papers, Tölle et al. studied harmonically trapped few-boson systems in an effective field theory based on contact interactions with running coupling constants Tölle et al. (2011, 2013).

Our oscillator EFT differs from these approaches. For the interaction we choose an oscillator space with a fixed oscillator frequency and a fixed maximum energy . The matrix elements of the interaction are taken from an EFT formulated in momentum space and evaluated at the discrete momentum eigenvalues of the kinetic energy in this fixed oscillator space. This reformulation, or projection, of a momentum-space EFT onto a finite oscillator model space requires us to re-fit the low-energy coefficients (LECs) of interactions at a given order of the EFT. We determine these by an optimization to scattering phase shifts (computed in the finite oscillator basis via the -matrix approach Heller and Yamani (1974); Shirokov et al. (2004)) and from deuteron properties. The power counting of the oscillator EFT is based on that of the underlying momentum-space EFT. The finite oscillator space introduces IR and UV cutoffs Stetcu et al. (2007); Coon et al. (2012); Furnstahl et al. (2012); More et al. (2013); Furnstahl et al. (2014); König et al. (2014); Furnstahl et al. (2015a); Wendt et al. (2015), and these are thus fixed for the interaction.

In practical many-body calculations we will keep the oscillator frequency and the interaction fixed at and , but employ the kinetic energy in larger model spaces. This increase of the model space increases (decreases) its UV (IR) cutoff but does not change any interaction matrix elements and thus leaves the IR and UV cutoff of the interaction unchanged. As UV convergence of the many-body calculations depends on the matrix elements of an interaction König et al. (2014), oscillator EFT guarantees this UV convergence by construction because no new potential matrix elements enter beyond . We stress that our notion of UV convergence relates to the convergence of the many-body calculations and should not be confused with the expectation that observables are independent of the regulator or cutoff. Infrared convergence builds up the exponential tail of bound-state wave functions in position space, as the effective IR cutoff of a finite nucleus is set by its radius. Thus, the increase of the model space for the kinetic energy achieves IR convergence. Regarding IR convergence, oscillator EFT is similar to the HOBET of Ref. Haxton (2007). In practice, IR-converged values for bound-state energies and radii can be obtained applying “Lüscher-like” formulas for the oscillator basis Furnstahl et al. (2012).

We view oscillator EFT similar to lattice EFT Lee (2009). The latter constructs an interaction on a lattice in position space while the former builds an interaction on a discrete (but non-equidistant) mesh in momentum space. In both EFTs the UV cutoff of the interaction is fixed once and for all, and LECs are adjusted to scattering data and bound states. The increase in lattice sites achieves IR convergence in lattice EFT, while the increase of the number of oscillator shells achieves IR convergence in oscillator EFT.

As we will see, the resulting EFT interaction in the oscillator basis exhibits a fast convergence, similar to the phenomenological JISP interaction Shirokov et al. (2007). From a practical point of view, our approach to oscillator EFT allows us to employ all of the existing infrastructure developed for nuclear calculations.

The discrete basis employed in this paper is actually the basis set of a discrete variable representation (DVR) Harris et al. (1965); Dickinson and Certain (1968); Light et al. (1985); Baye and Heenen (1986); Littlejohn et al. (2002); Light and Carrington (2007); Bulgac and McNeil Forbes (2013) in momentum space. While coordinate-space DVRs are particularly useful and popular in combination with local potentials, the results of this paper suggest that DVRs are also useful in momentum-space-based EFTs, because they facilitate the evaluation of matrix elements.

This paper is organized as follows. In Section II we analyze the momentum-space structure of a finite oscillator basis. In Section III we validate our approach by reproducing an interaction from chiral EFT at next-to-leading order (NLO). In Section IV we construct a NLO interaction from realistic phase shifts and employ this interaction in many-body calculations, demonstrating that converged binding energies and radii can be obtained for nuclei in the mass-100 region in model spaces with  10 to 14 without any further renormalization. We finally present our summary in Section V.

## Ii Theoretical considerations

In this Section we present the theoretical foundation of an EFT in the oscillator basis. We derive analytical expressions for the momentum eigenstates and eigenvalues in finite oscillator spaces and present useful formulas for interaction matrix elements.

### ii.1 Momentum states in finite oscillator spaces

The radial wave functions of oscillator basis states can be represented in terms of generalized Laguerre polynomials as

 ψn,l(r)= (1) (−1)n√2n!Γ(n+l+\nicefrac32)b3(rb)le−12(rb)2Ll+\nicefrac12n(r2b2).

Here, is the oscillator length expressed in terms of the nucleon mass and oscillator frequency . In what follows, we use units in which .

In free space, the spherical eigenstates of the momentum operator are denoted as with being continuous. The corresponding wave functions

 ⟨r,l|k,l⟩=√2πjl(kr) (2)

are spherical Bessel functions up to a normalization factor. These continuum states are normalized as

 ⟨k′,l|k,l⟩=∞∫0drr2⟨k′,l|r,l⟩⟨r,l|k,l⟩=δ(k−k′)kk′. (3)

Introducing the momentum-space representation of the radial oscillator wave functions via the Fourier-Bessel transform

 ~ψn,l(k)≡∞∫0drr2⟨k,l|r,l⟩ψn,l(r) =√2n!b3Γ(n+l+\nicefrac32)(kb)le−12k2b2Ll+\nicefrac12n(k2b2),

enables us to expand the continuous momentum states (2) in terms of the oscillator wave functions

 ⟨r,l|k,l⟩=∞∑n=0~ψn,l(k)ψn,l(r). (5)

As we want to develop an EFT, it is most important to understand the squared momentum operator . An immediate consequence of Eq. (2) is of course that the spherical Bessel functions are also eigenfunctions of with corresponding eigenvalues , i.e., . In the oscillator basis, the matrix representation of is tri-diagonal, with elements

 ⟨n′,l|^p2|n,l⟩ = b−2[(2n+l+\nicefrac32)δn′n (6) −√n(n+l+\nicefrac12)δn′+1n −√(n+1)(n+l+\nicefrac32)δn′−1n].

We want to solve the eigenvalue problem of the operator in a finite oscillator basis truncated at an energy . For partial waves with angular momentum the basis consists of wave functions (1) with , i.e., the sum in Eq. (5) is truncated at

 N≡[Nmax−l2]. (7)

Here denotes the integer part of . While clearly depends on and , we will suppress this dependence in what follows. Motivated by Eq. (5) we act with the matrix of on the component vector . The well-known three-term recurrence relation for Laguerre polynomials (see, e.g., Eq. 8.971(4) of Ref. Gradshteyn and Ryzhik (2000)) implies

 0 = (2n+l+\nicefrac32−b2k2)~ψn,l(k) (8) −√n(n+l+\nicefrac12)~ψn−1,l(k) −√(n+1)(n+l+\nicefrac32)~ψn+1,l(k),

and for our eigenvalue problem we arrive at

 ^p2⎛⎜ ⎜ ⎜ ⎜ ⎜ ⎜⎝~ψ0,l(k)⋮~ψN−1,l(k)~ψN,l(k)⎞⎟ ⎟ ⎟ ⎟ ⎟ ⎟⎠ = (9)

For such that , the second term on the right-hand side of Eq. (9) vanishes, and we obtain an eigenstate of the momentum operator (6) in the finite oscillator space.

Thus, momenta (with ) such that is a root of the the Laguerre polynomial solve the eigenvalue problem of the operator in the finite oscillator space. We recall that has roots. Thus, in a finite model space consisting of oscillator functions with the eigenvalues of the squared momentum operator are the roots of the Laguerre polynomial . We note that depends on the angular momentum as well as . To avoid a proliferation of indices, we suppress the latter dependence in what follows. By construction, the basis built on discrete momentum eigenstates is a DVR Baye and Heenen (1986); Littlejohn et al. (2002); Light and Carrington (2007).

Previous studies showed that the finite oscillator basis is equivalent to a spherical cavity at low energies Furnstahl et al. (2012), and the radius of this cavity is related to the wavelength of the discrete momentum eigenstate with lowest momentum. References More et al. (2013); Furnstahl et al. (2014) give analytical results for the lowest momentum eigenvalue in the limit of . The exact determination of the eigenvalues of the momentum operator in the present work allows us to give exact values for the radius of the cavity corresponding to the finite oscillator basis. This radius is relevant because it enters IR extrapolations Coon et al. (2012); More et al. (2013); Furnstahl et al. (2014). Let be the smallest root of the Laguerre polynomial , and let denote the smallest root of the spherical Bessel function . Then

 z0,lL(N,l)=k0,l (10)

defines the effective radius we seek.

The radial momentum eigenfunction corresponding to the eigenvalue in the partial wave with angular momentum has an expansion of the form

 ϕμ,l(r)≡cμ,lN∑n=0~ψn,l(kμ,l)ψn,l(r). (11)

This wave function is the projection of the spherical Bessel function (5) onto the finite oscillator space. It is also an eigenfunction of the momentum operator projected onto the finite oscillator space because the specific values of decouple this wave function from the excluded space. In Eq. (11) is a normalization constant that we need to determine. In order to do so we consider the overlap

 ⟨ϕμ,l|ϕν,l⟩=∞∫0drr2ϕμ,l(r)ϕν,l(r) =cμ,lcν,lN∑n=0~ψn,l(kμ,l)~ψn,l(kν,l) =cμ,lcν,l√(N+1)(N+l+\nicefrac32) ×~ψN,l(kμ,l)~ψN+1,l(kν,l)−~ψN+1,l(kμ,l)~ψN,l(kν,l)(k2μ,l−k2ν,l)b2.

Here, we used the Christoffel-Darboux formula for orthogonal polynomials, see, e.g., Eq. 8.974(1) of Ref. Gradshteyn and Ryzhik (2000). As and are roots of , we confirm orthogonality. For we use the rule by l’Hospital and find (with help of Eq. 8.974(2) of Ref. Gradshteyn and Ryzhik (2000))

 c−1μ,l=√(N+1)(N+l+\nicefrac32)kμ,lb~ψN,l(kμ,l). (13)

It is also useful to compute the overlap

 ⟨kμ,l,l|ϕν,l⟩ = √π2∞∫0drr2jl(kμ,lr)ϕν,l(r) (14) = δνμc−1μ,l.

This overlap vanishes for , thus, the eigenstates of the operator in finite oscillator spaces are orthogonal to the continuous momentum eigenstates when the latter are evaluated at the discrete momenta. This exact result is very useful for the computation of matrix elements of a potential operator .

For arbitrary continuous momenta we obtain from Eq. (II.1)

 ~ϕν,l(k)≡⟨k,l|ϕν,l⟩ = kν,l/bk2ν,l−k2~ψN+1,l(k). (15)

The wave function (15) is the Fourier-Bessel transform of the discrete radial momentum wave function . We note that Eq. (11) relates the discrete momentum eigenfunctions to the oscillator eigenstates via an orthogonal transformation, implying

 N∑μ=0c2μ,l~ψn,l(kμ,l)~ψn′,l(kμ,l)=δnn′, (16)

and

 N∑n=0~ψn,l(kμ,l)~ψn,l(kν,l)=δνμc−2μ,l. (17)

We remind the reader that the discrete set of momenta is fixed once a particular is chosen. Equation (16) can be used to relate oscillator basis functions to the discrete momentum eigenfunctions (11). Thus,

 ψn,l(r)=N∑μ=0cμ,l~ψn,l(kμ,l)ϕμ,l(r). (18)

### ii.2 Matrix elements of interactions from EFT

Nucleon-nucleon () interactions from EFT are typically available for continuous momenta in a partial-wave basis in form of the matrix elements . Numerical integration techniques are used to transform these matrix elements into the oscillator basis. However, there is a very simple approximative relationship between the matrix elements with continuous momenta and the matrix elements in the discrete momentum basis. This relationship is motivated by EFT arguments and we use it in our applications of oscillator EFT. We consider the matrix element

 ⟨k′,l′|^V|ϕμ,l⟩=∞∫0dkk2⟨k′,l′|^V|k,l⟩⟨k,l|ϕμ,l⟩ (19) = 12b3∞∫0dxxl+\nicefrac12e−x ×(⟨k′,l′|^V|x\nicefrac12b−1,l⟩⟨x\nicefrac12b−1,l|xle−x|ϕμ,l⟩).

Here, we introduced the dimensionless integration variable and factored out a weight function from the integrand (given in brackets) in preparation for the next step. We evaluate the integral using ()-point generalized Gauss-Laguerre quadrature based on the selected weight function. Thus, the matrix element (19) becomes

 ⟨k′,l′|^V|ϕμ,l⟩= 12b3N∑ν=0wν,l⟨k′,l′|^V|x\nicefrac12ν,lb−1,l⟩⟨x\nicefrac12ν,lb−1,l|ϕμ,l⟩xlν,le−xν,l +ΔN+1.

Here, are the roots of the Laguerre polynomial , the weights are

 wν,l ≡ Γ(N+l+\nicefrac52)xν,l(N+1)![(N+2)Ll+\nicefrac12N+2(xν,l)]2 (21) = Γ(N+l+\nicefrac32)xν,l(N+1)!(N+l+\nicefrac32)[Ll+\nicefrac12N(xν,l)]2,

and the error term is

 ΔN+1≡(N+1)!Γ(N+l+\nicefrac52)(2N+1)!f(2N+2)(ξ), (22)

see, e.g., Ref. Concus et al. (1963). For the weights, we also used Eq. 8.971(6) of Ref. Gradshteyn and Ryzhik (2000). In the error term, denotes the -th derivative of the integrand (given in round brackets) of Eq. (19), evaluated at which is somewhere in the integration domain.

We want to estimate the order of the error term when using EFT interactions. For this purpose, we write the potential as a sum of separable potentials

 V(k′,k)=∑avaga(k′)ga(k), (23)

and write

 ga(k)=(k/Λ)le−12k2b2~ga(k/Λ). (24)

Here, is a high-momentum cutoff. The function is an even function of its arguments (see, e.g., Ref. Machleidt and Entem (2011)), and can be expanded in a Taylor series

 ~ga(k)=∞∑n=0~g(n)a(0)n!(xΛ2b2)n. (25)

Here, we again used . With this expansion in mind, and noting that the wave function can be expanded in terms of oscillator wave functions, the integrand (given in round brackets) of Eq. (19) is a product of and a sum of Laguerre polynomials (from the wave function ). The (+1)-point Gauss-Laguerre integration is exact for monomials up to . As the wave function contains monomials up to , the Gauss-Laguerre integration becomes inexact for terms starting at in the Taylor series (25). Thus, in the error term (22) scales as . In the oscillator EFT, typical momenta scale as , and the error term scales as

 ΔN+1=O((k/Λ)2N+l+4). (26)

Therefore,

 ⟨k′,l′|^V|ϕμ,l⟩ = ⟨k′,l′|^V|kμ,l,l⟩cμ,l (27) +O((k/Λ)2N+l+4).

Repeating the calculation for the bra side yields the final result

 ⟨ϕν,l′|^V|ϕμ,l⟩ = cν,l′cμ,l⟨kν,l′,l′|^V|kμ,l,l⟩ (28) +O((k/Λ)2N+l+4).

In oscillator EFT, we will omit the correction term and set

 ⟨ϕν,l′|^V|ϕμ,l⟩ = cν,l′cμ,l⟨kν,l′,l′|^V|kμ,l,l⟩. (29)

We note that this assignment seems to be very natural for an EFT built on a finite number of discrete momentum states. For sufficiently large , the difference to the matrix element obtained from an exact integration can be view as a correction that is beyond the order of the power counting of the EFT we build upon. In Eq. (29) the matrix elements between the discrete and continuous momentum states simply differ by normalization factors because of the different normalization (14) of discrete and continuous momentum eigenstates. Thus, partial-wave decomposed matrix elements of any momentum-space operator can readily be used to compute the corresponding oscillator matrix elements. In the Appendix A we present an alternative motivation for the usage of Eq. (29).

In the remainder of this Subsection, we give useful formulas that relate the matrix elements and wave functions of the discrete momentum basis and the oscillator basis. We note that the oscillator basis states are related to the discrete momentum states via Eq. (18). Thus, we can also give an useful formula that transforms momentum-space matrix elements to the oscillator basis according to

 ⟨ψn′,l′|^V|ψn,l⟩= (31) N∑ν,μ=0c2ν,l′~ψn′,l′(kν,l′)⟨kν,l′,l′|^V|kμ,l,l⟩c2μ,l~ψn,l(kμ,l) +O(k2N+2).

This formula also reflects the well known fact that the oscillator basis mixes low- and high-momentum physics.

The relation between matrix elements in the oscillator basis and the discrete momentum basis is given by

 ⟨ϕν,l′|^V|ϕμ,l⟩= (32) cν,l′cμ,lN∑n,n′=0~ψn′,l′(kν,l′)⟨ψn′,l′|^V|ψn,l⟩~ψn,l(kμ,l).

Finally, we discuss the inversion of Eq. (31), e.g., for situations where scattering processes in the continuum have to be considered for interactions based on oscillator spaces. We obtain

 ⟨kν,l′,l′|^V|kμ,l,l⟩= N∑n,n′=0~ψn′,l′(kν,l′)⟨ψn′,l′|^V|ψn,l⟩~ψn,l(kμ,l) +O(k2N+l+4),

because of Eq. (17).

For arbitrary momenta, one needs to use the overlaps (15) in the evaluation of the matrix elements, and finds the generalization of Eq. (II.2) as

 ⟨k′,l′|^Π^V^Π|k,l⟩=~ψN+1,l′(k′)~ψN+1,l(k) (34) ×N∑ν,μ=0kν,l′/b(k′)2−k2ν,l′⟨ϕν,l′|^V|ϕμ,l⟩kμ,l/bk2−k2μ,l.

Here, we introduced the projection operator onto the finite oscillator space

 ^Π≡N∑ν=0Nmax−2ν∑l=0|ϕν,l⟩⟨ϕν,l|. (35)

We note that the projection operator acts as a UV (and IR) regulator. It is nonlocal, and can be written in many ways. Examples are

 ⟨k′,l|^Π|k,l⟩=N∑n=0~ψn,l(k′)~ψn,l(k) = √(N+1)(N+l+\nicefrac32)b2 × ~ψN,l(k)~ψN+1,l(k′)−~ψN+1,l(k)~ψN,l(k′)k2−(k′)2 (38) = ~ψN+1,l(k′)~ψN+1,l(k)b2 ×N∑ν=0k2ν,l[(k′)2−k2ν,l][k2−k2ν,l].

Here, the first identity comes directly from the definition of the projector in terms of the oscillator eigenfunctions. The second identity follows from the calculation displayed in Eq. (II.1), while the third identity follows from Eq. (15).

This presents us with an alternative motivation (but not derivation) of Eq. (29). We evaluate the projected matrix elements (34) at discrete momenta and find with Eq. (15) that

 ⟨kν,l′,l′|^Π^V^Π|kμ,l,l⟩=c−1ν,l′c−1μ,l⟨ϕν,l′|^V|ϕμ,l⟩. (39)

We note that

 ⟨kν,l′,l′|^Π^V^Π|kμ,l,l⟩=⟨kν,l′,l′|^V|kμ,l,l⟩ (40)

approximately holds for the discrete momenta in the finite oscillator space (cf. Eq. (28)).

## Iii Chiral interactions in finite oscillator bases

In this Section we present a proof-of-principle construction of a chiral interaction in the framework of the oscillator EFT. First, we study the effects that the truncation to a finite oscillator basis has on phase shifts of existing interactions. Second, we demonstrate that a momentum-space chiral interaction at NLO can be equivalently constructed in oscillator EFT.

We consider the chiral interactions NLO Entem and Machleidt (2003) and NLO Carlsson et al. (2016). Both interactions employ regulators of the form

 f(q)=exp⎡⎣−(qΛχ)2n⎤⎦. (41)

Here, is a relative momentum, is an integer, and is the high-momentum cutoff, specifically  MeV. We use in what follows. This cutoff needs to be distinguished from the (hard) UV cutoff König et al. (2014)

 ΛUV≈√2(Nmax+\nicefrac72)/b (42)

of the oscillator-EFT interaction.

Let us comment on using the projector (35) in combination with the regulator (41). In momentum space, the projector (35) is approximately the identity operator for momenta between the IR and UV cutoffs and of the oscillator basis. For momenta the projector (35) falls off as a Gaussian. For the regulator (41) we choose , which introduces a super-Gaussian falloff for momenta . As an example, let us consider  MeV and  MeV. Then, at the point where the super-Gaussian falloff goes over into a Gaussian falloff. Assuming a ratio that is typical for the power counting in chiral EFT, , and the asymptotic Gaussian falloff is not expected to introduce significant contributions to contact interactions at NLO. Eventually, one might want to consider removing the regulator from an oscillator-based EFT. At this moment however, it is also useful in taming the oscillations discussed below and shown in Fig. 1.

### iii.1 Effects of finite oscillator spaces on phase shifts

It is instructive to study the effects that a projection onto a finite oscillator basis has on phase shifts. For this purpose, we employ the well known interaction NLO Entem and Machleidt (2003). First, we transform its matrix elements to a finite oscillator space using numerically exact quadrature, and subsequently compute the phase shifts using the -matrix approach of Ref. Shirokov et al. (2004). The two parameter combinations ,  MeV, and ,  MeV respectively yield a UV cutoff  MeV, see Eq. (42). This considerably exceeds the high-momentum cutoff of the interaction. Figure 1 shows the resulting phase shifts in selected partial waves.

The projection onto finite oscillator bases introduces oscillations in the phase shifts. On the one hand, many ab initio calculations of atomic nuclei yield practically converged results for bound-state observables in model spaces with and oscillator frequencies around  MeV. On the other hand, oscillator spaces consisting of only 10 to 20 shells are clearly too small to capture all the information contained in the original interaction. We note (i) that the NLO interaction is formulated for arbitrary continuous momenta whereas oscillator EFT limits the evaluation to a few discrete momenta, and (ii) that the present cuts off the high-momentum tails of the NLO interaction.

For we expect to arrive at the original phase shifts, and this is supported in Fig. 1 by the reduced oscillations in the phase shifts compared to their counterparts. The period of oscillations is approximately given by the IR cutoff. In the bottom panel of Fig. 1, the oscillator spacing is increased to 80 MeV (for ) and 46 MeV (for ), respectively, yielding a UV cutoff MeV in the oscillator basis. The phase shift oscillations are significantly reduced for such large values of .

We note that the eigenvalues of the Hamiltonian matrix in oscillator shells play a special role for the phase shifts in uncoupled partial wave channels. At these energies, the eigenfunctions are standing waves with a Dirichlet boundary condition at the (energy-dependent) radius of the spherical cavity that is equivalent to the finite oscillator basis, and one can alternatively use this information in the computation of the phase shifts Luu et al. (2010); More et al. (2013). The filled circles in Fig. 1 indicate the values of the phase shifts at these energies, which are close to those of the original interaction.

Subsection II.2 discusses two ways of projecting momentum-space interactions onto a finite oscillator space. The first approach involves the determination of the matrix element (19) via an exact numerical integration over continuous momenta. The second approach, Eq. (29), uses -point Gauss-Laguerre quadrature to compute the integral in Eq. (19). This only requires us to evaluate the interaction at those momenta that are physically realized in the finite oscillator basis, which is more in the spirit of an EFT. Because we want to follow the oscillator EFT approach in later sections, we study the effect that the error term (22) associated with the -point Gauss-Laguerre quadrature has on the projected phase shifts. Figure 2 shows a comparison of projected NLO phase shifts obtained from the two projection approaches to the original NLO phase shifts. Overall, both versions yield very similar phase shifts. For the phase shifts associated with the -point Gauss-Laguerre quadrature, the oscillations seem to be somewhat reduced for small energies. Also, we find a notably improved agreement between the -point Gauss-Laguerre phase shifts and the original ones at the discrete eigenenergies of the scattering channel Hamiltonians, as indicated by the full circles. From now on we exclusively use the -point Gauss-Laguerre integration to compute matrix elements in oscillator EFT.

### iii.2 Reproduction of the phase shifts of a NLO interaction

In what follows we employ the oscillator EFT at NLO. All matrix elements in oscillator EFT are based on Eq. (29), i.e., the matrix elements from continuum momentum space are evaluated at the discrete momenta of the finite oscillator basis. At this order the chiral interactions depend on 11 LECs and exhibit sufficient complexity to qualitatively describe nuclear properties. Specifically, in this Section we aim at reproducing the phase shifts and selected deuteron properties of the chiral interaction NLO Carlsson et al. (2016) by optimizing these 11 LECs. Throughout this work, the fits evaluate the phase shifts at 20 equidistant energies in the laboratory energy range up to 350 MeV, with weights . We note that more sophisticated weights (including also the pion mass or the oscillatory patterns) would be needed for a quantification of uncertainties, see, e.g. Refs. Stump et al. (2001); Carlsson et al. (2016); Furnstahl et al. (2015b). In this work we only investigate the feasibility of an oscillator EFT with regard to the computation of heavy nuclei.

The main goal of oscillator EFT is to enable the computation of heavy nuclei. Therefore, we set for the interaction. This allows us to perform IR extrapolations based on calculations in spaces with , 12, and 14 for the kinetic energy. To determine we study phase shifts in selected partial waves in Fig. 3. The first value,  MeV, again corresponds to a UV cutoff of  MeV and yields the familiar oscillations in the phase shifts. The second value,  MeV, corresponds to a UV cutoff  MeV that significantly exceeds the chiral cutoff  MeV. In this case, the oscillations are drastically reduced. Consequently, we use MeV in the following.

In Fig. 4 we compare further phase shifts from the oscillator EFT to the phase shifts. We note that the channel is a prediction because there is no corresponding LEC at this chiral order.

For completeness, and to assess the large- behavior, we also show the results of an optimization in oscillator EFT with and MeV ( MeV) as green crosses. In general, oscillator EFT reproduces the phase shifts well over the whole energy range up to the pion-production threshold.

In the deuteron channel we fit not only to phase shifts but also to its binding energy , point-proton radius , and quadrupole moment . For , we reproduce the deuteron properties well, as can be seen in Table 1. For the effective interaction, it becomes more difficult to simultaneously reproduce all data. Therefore, we relax the requirement to reproduce the quadrupole moment in favor of the other deuteron properties and the phase shifts. We converge the deuteron calculations by employing large model spaces of 100 oscillator shells where only the kinetic energy acts beyond the used to define the effective interaction.

In Table 1 we also compare the LECs of the effective interactions for = 10 and 80 to the LECs of NLO. For = 80, our approach quite accurately reproduces the NLO LECs, with being the only exception. Because this LEC is associated with the deuteron channel, its value is affected by the weighting of the deuteron properties in the fit, and we expect that a better reproduction can be achieved by assigning different weights to the deuteron properties. More importantly, with the exception of , the values of the LECs for the = 10 interaction are very similar to the LECs of the NLO interaction. Also, all of the values are of natural size.

Let us also address the sensitivity of our results to the UV cutoff of the employed oscillator space. For this purpose we keep  MeV fixed and optimize two more interactions defined in model spaces with and , respectively. We recall that the UV cutoff of the model space increases with increasing , and that corresponds to UV cutoff 650, 700, 750 MeV, respectively. Resulting few-body observables from these interactions are shown in Table 2, and also compared to the infinite-space interaction . We see that results converge slowly toward those of the interaction as is increased. We also note that the few-body observables from the different interactions exhibit differences that are consistent with uncertainty expectations at NLO Epelbaum et al. (2015); Lynn et al. (2016); Carlsson et al. (2016).

Clearly, the NLO interaction from oscillator EFT differs from NLO through the complicated projection that introduces IR and UV cutoffs and is highly nonlocal, see Eq. (35). In the EFT sense the difference between these interactions should be beyond the order at which we are currently operating. While we cannot prove this equivalence, the numerical results of this Section encourage us to pursue the construction of a chiral NLO interaction within oscillator EFT by optimization to the phase shifts from a high-precision potential in the following Section.

## Iv NLO interactions in oscillator EFT and many-body results

In what follows, we set for the interaction. On the one hand, lower oscillator frequencies correspond to larger oscillator lengths and lead to a rapid IR convergence. On the other hand, lower oscillator frequencies also correspond to lower UV cutoffs. In what follows we choose  MeV which results in a UV cutoff  MeV and an IR length  fm according to Eq. (10). Considering the tail of the regulator function (41) we set  MeV, which ensures that significantly exceeds .

Practical calculations are performed in the laboratory system, and use the interaction together with the intrinsic kinetic energy in oscillator spaces with . Here and are the radial and angular momentum quantum numbers of the harmonic oscillator in the laboratory frame. Results for are feasible and will allow us to perform IR extrapolations for bound-state energies and radii.

In this Section, we construct a chiral interaction at NLO from realistic phase shifts, and subsequently utilize it in coupled-cluster calculations of He, O, Ca, Zr, and Sn. Our main objectives are (i) to present a proof-of-principle optimization of a realistic interaction within the framework of the oscillator EFT, and (ii) to demonstrate that such an interaction converges fast even in heavy nuclei. We compute the matrix elements in oscillator EFT, i.e., based on Eq. (29) and omit the high-order correction terms.

In the optimization, the low-energy coefficients are obtained from a fit to realistic scattering data (represented here by phase shifts from the CD-Bonn potential Machleidt (2001)) and deuteron properties. The fitting procedure is identical to the one described in Section III. Figure 5 presents the resulting phase shifts for a selection of scattering channels.

We reproduce the phase shifts over a large energy range for several partial waves. The channel is an obvious exception, as it deviates more clearly from the CD-Bonn phase shifts at higher energies, but we note that at NLO there is no LEC to adjust in this channel. For , the deviations are also considerable, however, it is of NLO quality, see Fig. 4. For the deuteron, we obtain a good reproduction of the binding energy and radius, and a reasonably well result for the quadrupole moment, as shown in Table 3. The LECs, also shown in Table 3, are natural in size and similar to NLO; the largest deviation is for the LEC which is only about of its NLO counterpart.