Improving Predictions with Reliable Extrapolation Schemes and Better Understanding of Factorization

New insights into the inter-nucleon interactions, developments in many-body technology, and the surge in computational capabilities has led to phenomenal progress in low-energy nuclear physics in the past few years. Nonetheless, many calculations still lack a robust uncertainty quantification which is essential for making reliable predictions. In this work we investigate two distinct sources of uncertainty and develop ways to account for them.

Harmonic oscillator basis expansions are widely used in ab-initio nuclear structure calculations. Finite computational resources usually require that the basis be truncated before observables are fully converged, necessitating reliable extrapolation schemes. It has been demonstrated recently that errors introduced from basis truncation can be taken into account by focusing on the infrared and ultraviolet cutoffs induced by a truncated basis. We show that a finite oscillator basis effectively imposes a hard-wall boundary condition in coordinate space. We accurately determine the position of the hard-wall as a function of oscillator space parameters, derive infrared extrapolation formulas for the energy and other observables, and discuss the extension of this approach to higher angular momentum and to other localized bases. We exploit the duality of the harmonic oscillator to account for the errors introduced by a finite ultraviolet cutoff.

Nucleon knockout reactions have been widely used to study and understand nuclear properties. Such an analysis implicitly assumes that the effects of the probe can be separated from the physics of the target nucleus. This factorization between nuclear structure and reaction components depends on the renormalization scale and scheme, and has not been well understood. But it is potentially critical for interpreting experiments and for extracting process-independent nuclear properties. We use a class of unitary transformations called the similarity renormalization group (SRG) transformations to systematically study the scale dependence of factorization for the simplest knockout process of deuteron electrodisintegration. We find that the extent of scale dependence depends strongly on kinematics, but in a systematic way. We find a relatively weak scale dependence at the quasi-free kinematics that gets progressively stronger as one moves away from the quasi-free region. Based on examination of the relevant overlap matrix elements, we are able to qualitatively explain and even predict the nature of scale dependence based on the kinematics under consideration.

Improving Predictions with Reliable Extrapolation Schemes and Better Understanding of Factorization

Presented in Partial Fulfillment of the Requirements for

the Degree Doctor of Philosophy in the

Graduate School of The Ohio State University
Sushant N. More, M.S.

Graduate Program in Physics

The Ohio State University

Dissertation Committee:Prof. Richard Furnstahl, AdvisorProf. Robert PerryProf. Michael LisaProf. John Beacom

© Copyright by

Sushant N. More


To the memory of my late grandparents, who valued integrity and education above everything else.

Acknowledgments \@mkbothAcknowledgmentsAcknowledgments

It was an absolute pleasure to work with Prof. Richard Furnstahl. Even after four years of working closely with Dick, I am still in awe of the time and effort he invests in his students. His passion towards physics, dedication towards mentorship, and caring about the scientific community is exemplary. I certainly consider it a privilege to be his padawan and hope to have imbibed at least part of his scientific ethos.

I am also thankful to my committee members, in particular Robert Perry. Their constructive criticism and feedback has helped me develop professionally. I was fortunate to learn from a few exceptional teachers at Ohio State. I thank them for instilling in me a joy for learning advanced physics. I would like to express my gratitude towards Scott Bogner and Thomas Papenbrock; they have been a constant source of support, and I look up to them as mentors. I am grateful to Sabine Jeschonnek for helping with the deuteron disintegration theory. I have also benefited tremendously from the interactions with professors and peers in the wider nuclear physics community.

My time at Ohio State coincided with some of the best post-docs in the group—Heiko, Kai, and Sebastian. They were always ready to help and my PhD experience has been enriched immensely through interactions with them. I have worked closest with Sebastian; I admire his efficiency and his attention to details.

I was also lucky to overlap with excellent fellow graduate students in the group. It was fun growing up with my physics siblings—Kyle, Sarah, and Alex. Special mention goes to Kyle for doing a great job of an elder sibling. It was a pleasure to share office space with Brian, Russell, and Khal. I thank them for being a jolly company and for putting up with my idiosyncrasies. The other physics graduate students that I would like to acknowledge are my gym buddies Hudson and Chuck and my former officemate Omar. I have enjoyed discussion with them over a full gamut of topics and have always turned to them for friendly advice when needed.

My friend Kaushik and I landed in the US together five years ago. We have had a good run as flatmates since then. From teaching me cooking to being my driving instructor, I am indebted to him for many things. One of the best things that happened to me at OSU was joining the student group Sankalpa. I thank Gayatri, Vinayak, and Keerthi for getting me started along that path. Friends formed through Sankalpa have been like a family to me and have helped me grow as a person. I thank Mithila, Dhruv, and Garima for being my safety net; Janani, Pooja, Shubhendu, and Keshav for intellectually engaging discussions; and Vipul, Nandan, Subhasree, Bakul, and Madhura for being a lot of fun.

I would like to take this opportunity to visit my roots and express gratitude towards my undergraduate teachers, friends, and relatives back in India. They trusted in me when I had little to show. I value the decade-long friendship of my pals Sheetal Kumar, Ajit, Shishir, Vivek, and Lakshmi Priya.

Finally, I must put on record the unremitting love of my doting yet disciplinarian parents and my headstrong but sweet younger sister. Their unwavering encouragement through thick and thin keeps me going and motivates me to scale new heights.


May 2011 Dual B.S. M.S., Indian Institute of Science Education and Research (IISER), Pune, India

December 2013 M.S, The Ohio State University, Columbus, Ohio


Deuteron electrodisintegration with unitarily evolved potentials
S. N. More, S. König, R. J. Furnstahl, and K. Hebeler, Phys. Rev. C 92, 064002 (2015)

Ultraviolet extrapolations in finite oscillator bases
S. König, S. K. Bogner, R. J. Furnstahl, S. N. More, and T. Papenbrock, Phys. Rev. C 90, 064007 (2014)

Systematic expansion for infrared oscillator basis extrapolations
R. J. Furnstahl, S. N. More, and T. Papenbrock, Phys. Rev. C 89, 044301 (2014)

Universal properties of infrared oscillator basis extrapolations
S. N. More, R. J. Furnstahl, A. Ekström, G. Hagen, and T. Papenbrock, Phys. Rev. C 87, 044326 (2013)

Non-local linear stability of ion beam eroded surfaces
S. N. More, R. Kree. Applied Surface Science 258 (2012) 4179-4185

Fields of Study

Major Field: Physics

Table of Contents\@mkbothTable of ContentsTable of Contents



List of Figures\@mkbothList of FiguresList of Figures



List of Tables\@mkbothList of TablesList of Tables



Chapter 1 Introduction

1.1 Overview of Nuclear Physics

Nuclear physics deals with the properties of atomic nuclei. The nuclear landscape shown in Fig. 1.1 has been the traditional playground for nuclear physics. The questions historically driving nuclear physics have been: how do protons and neutrons make stable nuclei and rare isotopes? What are the limits of nuclear existence? What are the nuclear binding energies, excitation spectra, radii and so on? We would also like to describe nuclear reactions, make predictions about the shape of the nuclei and understand how the shape dictates the nuclear properties.

Figure 1.1: Nuclear landscape. A total of isotopes (black squares) are stable on the time scale of the solar system. As more protons or neutrons are added to these stable nuclei, we enter the regime of short-lived radioactive nuclei (green squares). ‘Drip lines’ mark the limit of nuclear existence, where the last nucleon is no longer bound. The uncertainties around drip lines (in red) were obtained by averaging the results of different theoretical models. Figure from [1].

By the mid-1970s, it was generally accepted that the nucleons (proton and neutrons) and other hadrons are composed of quarks, and that the quarks are held together through the exchange of gluons [2]. The following decades witnessed rapid development in the fundamental theory of strong interactions describing the interactions between quarks and gluons. This theory goes by the name of Quantum Chromodynamics (QCD) [3]. One of the active areas of investigation is obtaining the hadron structure from QCD. This includes, for example, understanding the origin of proton spin, which is studied experimentally at Jefferson Laboratory [4]. A related focus area is understanding the nature of the quark-gluon plasma (QGP)—the phase in which the universe is believed to exist for up to a few milliseconds after the Big Bang [5]. The energies involved in this subfield (few GeVs) are higher than the energies in ‘traditional’ nuclear physics (few MeVs) introduced in the opening paragraph. Therefore it is conventional to refer to the two subfields as high-energy nuclear physics and low-energy nuclear physics. The work in this thesis will mainly focus on questions in low-energy nuclear physics (LENP).

Apart from the questions at the core of LENP, mentioned in connection to Fig. 1.1, inputs from LENP are immensely important in other areas as well. One such broad area is that of nuclear astrophysics. The majority of the stable and known nuclei shown in Fig. 1.1 were formed in big bang, stellar, or supernova nucleosynthesis [6, 7]. Inputs from LENP are critical in understanding the processes involved in nucleosynthesis and predicting the observed abundances of isotopes. Neutron stars are another fascinating astrophysical objects for low-energy nuclear physicists. We would like to determine the equation of state for neutron stars and understand how and why stars explode [8].

Finally, there are questions about the fundamental symmetries of the universe where nuclear physics hopes to make significant contributions. For instance, why is there more matter than antimatter in the universe? What is the nature of dark matter [9]? What is the nature of the neutrinos (Majorana or Dirac fermions) and how have they shaped the evolution of the universe? In fact, as we will see later, accurate calculations of nuclear matrix elements are critical for the experiments undertaken to understand the nature of neutrinos [10].

In addition to the broad scientific impact of nuclear physics that we have already mentioned, it also has many real-life applications. Our knowledge of nuclei and ability to produce them has led to an increase in the quality of life for humankind. Applications of nuclear physics encompass a diverse domain including but not limited to energy, security, medicine, radioisotope dating, and material sciences.

1.2 Checkered past; promising future

In 1935, Hideki Yukawa proposed the seminal idea of nuclear interactions being mediated by a massive boson [11]. This could explain how protons and neutrons would stay bound in a nucleus, overcoming the Coulomb repulsion between protons. The fact that such a model described scattering data well at low energies (few MeVs) and the eventual discovery of pions in 1947 led to a wide acceptance of this model. Very soon other heavy mesons (, , ) were discovered as well. Scattering experiments also indicated that the strength of the nuclear potential depended on distance and at short distances the potential was repulsive.

By 1950, there emerged an industry for coming up with better nuclear potentials. These boson-exchange models shared some common features. The long-range part of the nucleonic interaction was given by pion exchange, the intermediate range was governed by multiple (mostly two) pion exchange, and short-range repulsion was thought to be because of overlap of nucleons. When heavy mesons were discovered, they were added to the intermediate range sector. The pion, being the lightest meson, has the longest range. These general considerations form the basis for phenomenological potentials used even today as seen in Fig. 1.2.

Figure 1.2: AV18 [12], Reid93 [13], and Bonn [14] potentials for channel as functions of internucleonic distance. These potentials accurately describe neutron-proton scattering up to laboratory energies of 300 MeV. Regions I, II, and III correspond to long-range, intermediate-range, and short-range parts discussed in the text. Figure from [15].

This intense effort is well summarized by Hans Bethe’s quote in his essay ‘What Holds the Nucleus Together?’ in Scientific American (1953): “In the past quarter century physicists have devoted a huge amount of experimentation and mental labor to this problem – probably more man-hours than have been given to any other scientific question in the history of mankind.” The boson models did not have a smooth sailing though. In particular the intermediate range multi-pion sector was beset with problems. The pessimism this resulted in is palpable in Marvin Goldberger’s comment in 1960: “There are few problems in nuclear theoretical physics which have attracted more attention that that of trying to determine the fundamental interaction between two nucleons. It is also true that scarcely ever has the world of physics owed so little to so many…It is hard to believe that many of the authors are talking about the same problem or, in fact, that they know what the problem is.” A running joke was that nuclear physics is really ‘unclear’ physics!

There was relatively slow progress with regards to the development of internucleonic potentials in 1970’s and 80’s. However, this period saw a rapid development of perturbative QCD. It was realized that the nucleons and pions are composed of quarks which are held together by exchange of gluons (cf. Fig. 1.3). This pushed the effort to derive the nuclear potential from the ‘fundamental’ theory of QCD.

(a) Inter-nucleon interaction in 1940
(b) Inter-nucleon interaction in 1980
Figure 1.3: Evolution of the inter-nucleon interaction picture over time. Figures from a talk by Witold Nazarewicz.

However, the effort to replace the hadronic descriptions at ordinary nuclear densities with a quark description as in Fig. 1.2(b) was not very fruitful.

Figure 1.4: Summary of measurements of the QCD coupling as a function of energy scale [16].

As seen in Fig. 1.4, the strength of the QCD coupling increases with decreasing energies. This makes QCD non-perturbative in the low-energy regime of nuclear physics, limiting the success of analytical calculations.

Another aspect that makes low-energy nuclear physics difficult is that it is a many-body problem. It exhibits some emergent phenomena that are difficult to capture in a reductionist approach. This issue has been well-summarized by the famous article ‘More is different’ by Phillip Anderson [17] (albeit with a focus on many-body problem in condensed matter).

Despite all these challenges, great strides have been made in LENP in the last few decades. As indicated in Fig. 1.5

Figure 1.5: LENP version of Moore’s law and its violation. axis is the mass number of nuclei that can be calculated from ab-initio calculations. In past few years, it has been possible to push the ab-initio frontier to heavier nuclei. Figure from a talk by Gaute Hagen.

particularly the last few years have seen an explosion in the capabilities of low-energy nuclear theory. This phenomenal progress has been possible due to the combination of a few factors—new insights about the nuclear force, developments in many-body technology, and a surge in computational capabilities. In the following sections, we look briefly at each of these developments which will lead us to how the author’s PhD work fits into the bigger picture.

1.3 Understanding the Force

We saw that non-pertubativeness of QCD at low energies (cf. Fig. 1.4) motivated phenomenological descriptions of nuclear forces. The two popular categories of phenomenological interactions are the meson exchange models (which we touched upon during discussion of Fig. 1.2) and local111The potential is local if . phenomenological potentials. Examples of these include the Bonn potentials [14, 18, 19] and the Argonne potential [12].

The Argonne potential is probably the most widely used phenemenological potential. One of the reasons being that until recently, it was the only precision interaction usable for Quantum Monte Carlo calculations. The Argonne interactions are built by writing down all the operators that satisfy the required symmetries—translational and Galilean invariance, rotational invariance in space and spin, rotational invariance in isospin, time reversal, and spatial reflection. These operators are given below.


is the tensor force.

There are a total of 14 operators in Eq. 1.1. The AV18 potential has four more operators which are the charge-dependent and charge-symmetry breaking terms; they are small but needed to get for np, nn, and pp scattering. There are also AV8, AV6 potentials which use a limited set of operators.

The AV18 potential is written as


where is the inter-nucleon separation and . The coefficients in the potential are fit to nucleon scattering up to MeV, and to deuteron bound state properties such as binding energy, radii, and quadrupole moment. A similar exercise has been done for 3N (3-body) interactions. However, the large number of possible three-body operators makes it difficult to get rid of model dependence.

Meson exchange models are formulated in terms of exchange of mesons taking into account the quantum nature of mesons (scalar, vector, pseudo-scalar, so on). The masses are those of real mesons, but couplings are fit parameters. In the simplest form, the interaction is a sum of Yukawa potentials,


Meson exchange potentials and the AV18 potential share the common shortcoming that there is no scope for systematic improvements. It is also unclear how to seek model independence and do robust uncertainty quantification. Finally, these models offer limited guidance on the strength and relevance of three- and higher-body forces.

1.3.1 Chiral EFT

Figure 1.6: Hierarchy of degrees of freedom and associated energy scales in nuclear physics [20].

An intriguing aspect of the world we live in, is that there are interesting phenomena at virtually all energy and length scales we can probe. From TeV energies at the Large Hadron Collider (LHC) to the life-defining process of respiration which has the energy scale of only few meV, there are physical processes of interest at each step. Nuclear physics spans a wide range of energy and length scales (cf. Fig. 1.6); a wider range than most subfields. This hierarchy provides both challenges and opportunities.

Figure 1.6 indicates the relevant degrees of freedom for the given energy scales. Even though degrees of freedom are a matter of choice, in practice, appropriate degrees of freedom often dictate the success of a theory. To quote Steven Weinberg [21]: “You can use any degrees of freedom you want, but if you use the wrong ones, you’ll be sorry.” Weinberg in his seminal paper [22] applied the concept of effective field theory (EFT) to low-energy QCD. This effort proceeds by writing down the most general Lagrangian consistent with the (approximate) symmetries.

QCD has the expected symmetries of translational, Galilean, and rotational invariance, and spatial reflection and time reversal. Along with that, in the limit of vanishing quark masses, the QCD Lagrangian also possess an exact chiral symmetry [23]. If the chiral symmetry holds, “left”- and “right”-handed fields do not mix. As with any other continuous symmetry, spontaneous breaking of chiral symmetry leads to massless Goldstone boson(s). The masses of the up and down quark (quarks relevant in LENP) are both small ( MeV [16]), but non-zero. Therefore the chiral symmetry of the QCD Lagrangian is only approximate and the resulting Goldstone boson—pion—is light (compared to mass of nucleon), but not massless.

Chiral EFT (-EFT) uses nucleons and pions as degrees of freedom 222Delta-full chiral EFTs also include the —a resonant state of the nucleon—as a degree of freedom. . Heavy mesons are “integrated out”. The crucial difference that distinguishes -EFTs from meson theories of the 1950s is that they are constrained by chiral symmetry. Broken chiral symmetry serves as a connection with the underlying theory of QCD. A major advantage of -EFT is that it permits systematic improvements and allows the possibility of having reliable uncertainty quantification.

Figure 1.7: Diagrams for a chiral Lagrangian at each order. Solid lines are the nucleons and dashed lines pions. Figure from [24].

Diagrams in a -EFT are shown in Fig. 1.7. The three- and higher-body forces appear naturally in -EFT with an expected hierarchy (higher-body forces are suppressed successively). The coupling constants at the vertices in Fig. 1.7 are called low-energy coupling constants (LECs); they encode the QCD physics and cannot be calculated in -EFT. LECs are fit to experimental data. Lattice QCD provides a promising method for extracting them in the near future (cf. Subsec. 1.4.1). Despite some concerns over power-counting (see Ref. [25] and references therein), it is fair to say that -EFT has been a major breakthrough in low-energy nuclear theory. For more details and applications of -EFT, please see Refs. [24] and [26].

1.4 Many-body methods

We learned in kindergarten quantum mechanics that a two-body problem can be solved by reducing it to a one-body problem. The many-body problem is far more relevant to the nuclear physics, and is more challenging. For the nuclear properties such as the bound state energies, the equation we need to solve is the time-independent Schrödinger equation


The Hamiltonian is given by


where is the kinetic energy and are two- and three-body nuclear potentials (e.g., from -EFT or phenomenological potentials).

Over the years various methods have been developed to tackle this problem. We will list some of the broad categories below. We only provide a brief explanation for each of them and refer the reader to the cited references for details.

  • Direct diagonalization: This category involves expanding the many-body wave function in an appropriate complete basis. Very often this basis is chosen to be the Slater determinant of harmonic oscillator (HO) wave functions. HO wave functions form a complete basis with discrete energy levels that depend on only one scale (the oscillator frequency). Moreover, HO wave functions are simple analytic functions that go to zero at large distance just as nuclear bound state wave functions.

    Finite computational power forces us to truncate the infinite sum of Slater determinants at some point. One of the common truncation scheme is the truncation. It keeps all particle states such that


    With truncation, solving the many-body Schrödinger equation (Eq. 1.4) becomes a matrix diagonalization problem, and can be solved using standard algorithms. The truncation in HO basis also allows separation of center-of-mass motion from relative motion, ensuring that only intrinsic properties are being calculated.

    (a) He energies as a function of .
    (b) Li energies as a function of .
    Figure 1.8: Convergence of energies as a function of the truncation parameter [27].

    Figure 1.8 shows convergence plots for ground state energies for He and Li. As is increased the energies approach asymptotic values333When the asymptotic value doesn’t match the experimental value, it points to the shortcomings in the nuclear Hamiltonian..

    Figure 1.9: Matrix dimension grows factorially with the number of nucleons. Figure courtesy of Pieter Maris.

    A limitation of the diagonalization method is that the size of the Hamiltonian matrix that we need to diagonalize grows factorially as we go to higher and higher (cf. Fig. 1.9). This often forces us to truncate the basis before convergence is reached (cf. Fig. 1.7(b)), and necessitates development of reliable extrapolation schemes. Chapter 2 describes the author’s original work devoted to the development of extrapolation schemes relevant to this problem.

    Many different methods fall under the broad umbrella of diagonalization methods. This includes the traditional shell model [28], the No Core Shell Model (NCSM) [29], the No Core Full Configuration (NCFC) [30], and the Importance Truncated No Core Shell Model (IT-NCSM) [31]. The last three methods primarily differ in the nature of the truncation they employ.

  • Monte Carlo methods: These methods use the imaginary time evolution of the Schrödinger equation


    This equation can be solved by making an ansatz for . can be expanded in a complete basis of eigenvectors of as


    As long as the trial wave function is not orthogonal to the actual ground state, i.e, in Eq. 1.8, it can be shown that the imaginary time evolution (Eq. 1.7) projects out the ground state in the limit. There are many different ways to do the stochastic time evolution of Eq. 1.7 to project out the ground state. The two most popular in nuclear physics are the Green’s function Monte Carlo (GFMC) [32, 33] and the Auxiliary Field Diffusion Monte Carlo (AFDMC) [34, 35]. Monte Carlo methods work best with a local potential. The AV18 potential has been the potential of choice for these methods. Recently, it has been possible to derive the low-order -EFT in local form, making its use possible for Monte Carlo methods [35].

  • Coupled Cluster: The energy of the many-body state is given by


    The coupled cluster (CC) method tries to build the state from the reference state (which for instance can be a Hartree-Fock state) using the transformation


    Thus, Eq. 1.9 becomes


    The cluster operator in Eq. 1.10 is defined with respect to the reference state.


    generates -particles--holes excitations. In practice, Eq. 1.12 is truncated at (more recently corrections from are included as well). CC scales much better than diagonalization or Monte Carlo methods discussed above and is therefore possible to use for medium-mass nuclei [36].

  • Density Functional Theory: Density Functional Theory (DFT) is based on the principle that the many-body ground state can be written as a functional of the density , i.e, . Consequently, the energy (or any other observable) can be written as a functional of density


    DFT proceeds by writing down an energy density functional (EDF) guided by intuition and general theoretical arguments [37, 38]. Inputs from experiments and exact calculations for simple few-body systems are also used for constructing the EDF. The properties of the physical system are then found by the two step minimization of the EDF—first minimization is at a fixed density and the second minimization is over . Once the EDF is decided upon, DFT does not scale prohibitively with the number of nucleons , and is therefore the most popular for calculating properties of heavy-mass nuclei.

  • IM-SRG: In-Medium Similarity Renormalization Group (IM-SRG) is based on the SRG technique we will look at in Subsec. 1.5.1. It uses a series of continuous unitary transformations to decouple the ground state of the many-body Hamiltonian from the excitations. IM-SRG has made it possible to apply the ab-initio (starting with 2N and 3N forces) methods to medium-mass nuclei and beyond. Please see the Ref. [39] for a recent review on progress acheived by IM-SRG.

1.4.1 “It from the bit”444The phrase “It from the bit” was originally used by John Wheeler while elucidating his ideas on digital physics. —lattice theories

Lattice QCD

We saw through Fig. 1.4 that the largeness of the QCD coupling at low-energies makes it unamenable to analytical calculations. A well-established non-perturbative approach in this regime is lattice QCD [40]. In lattice QCD, one discretizes space-time; fields representing quarks are defined at lattice sites and the gluon fields are defined on links connecting neighboring sites.

Ideally we would like the lattice size to be as large as possible and the lattice spacing to be as small as possible. However, lattice calculations are computationally extremely intensive, thereby severely constraining the lattice size and spacing. Moreover, the computational cost of simulations scales with the quark mass roughly as [41]. The simulations are therefore often done at quark masses larger than the physical quark masses. From the Gell-Mann-Oakes-Renner relation , and therefore the pion mass in the lattice calculations is larger than its physical value as well. The results are then extrapolated down to physical quark (or pion) mass. In the nuclear case, for accurate extrapolation one must also take into account the crucial non-analytic structure associated with chiral symmetry breaking [42].

Recent calculations have been able to use the physical pion mass though only for very light hadrons [43]. In near future, we hope to extract the low-energy constants in Fig. 1.7 from lattice QCD.

Lattice EFT

Quarks and gluons have many degrees of freedom in terms of spin, color charge, and flavor. This along with strong nonlinearity and non-perturbativeness of the problem makes getting nuclear physics from lattice QCD computationally difficult. An alternative approach is to have nucleons on the lattice site rather than the quarks.

(a) Lattice QCD            
(b) Lattice EFT            
Figure 1.10: Lattice QCD vs. lattice EFT. Figures from a talk by Dean Lee.

The difference between the two methods is illustrated in Fig. 1.10. Lattice EFT combines the framework of effective field theory and computational lattice methods and is a promising tool for studying light nuclear systems.

Lattice theories inherently come with the associated graininess and the finite size. Thus, we have a cutoff for both the maximum length and the maximum momentum scale that we can have on a lattice. This is qualitatively similar to the harmonic oscillator basis truncation problem that we will look extensively at in Chapter 2. Various methods have been used to obtain the continuum limit from the lattice. We will touch upon these in Subsec. 2.3.1, where we look at similarities and differences between the lattice methods and our work with oscillator basis truncation.

1.5 RG techniques

We saw in Fig. 1.8 that the convergence in many-body calculations is slow. To understand why this is the case, recall from Fig. 1.2 that the nuclear potentials have a strong short-range repulsion. This hard core leads to high-momentum components in the potentials.

Figure 1.11: The AV18 potential in momentum space. Figure from [44].

This can be seen in Fig. 1.11 where we plot the AV18 potential in the channel in momentum space. We note that for the AV18 potential, only for . However, the Fermi momentum of the nucleon in a heavy nucleus like Pb is only about . The Fermi momentum sets the momentum scale of low-energy nuclear processes we wish to study. Thus, we have a mismatch of resolution scales; the processes we wish to describe are (), whereas the momentum scale in the potential is much higher. To use the analogy due to Tom Banks, ‘it is like trying to understand the properties of waves in the ocean in terms of Feynmann diagrams’. Though in principle this can be done, it makes calculations intractably complicated.

This bring us back to Fig. 1.6. The progression from top to bottom in Fig. 1.6 can be thought of as reduction in resolution. This can be established theoretically using renormalization group (RG) methods. As mentioned before, the focus in LENP is the intermediate region, where nucleons are the degrees of freedom. But even within this limited region, the concept of changing resolution by RG methods has been extremely advantageous [44].

(a) Nucleus under a high-resolution probe.
(b) Fine details (nucleon substructure) not resolved when probed at low energies.
(c) A painting by Georges Seurat, which is an example of pointillism. Small, distinct dots of color are applied to form a pattern. The pattern would be lost under a high-resolution probe. Image from Wikimedia Commons.
Figure 1.12: Physics interpretation often changes with resolution.

We would like to stress that contrary to the popular notion, high resolution is not always the best thing, especially when the processes we are looking at are low-momentum. Also, though the value of the calculated observable is independent of the resolution, the physical interpretation often changes with the resolution. This is illustrated in Fig. 1.12.

(a) An illustration of a low-pass filter.
(b) phase shifts for the AV18 potential and for the AV18 potential after the low-pass filter which sets [45].
Figure 1.13: Low-pass filter on nuclear potential—illustration and effect on phase shifts.

One of the methods to get rid of the high-momentum components is to apply a ‘low-pass filter’ on the potential (cf. Fig. 1.12(a)). This is routinely done, for example, in image processing. Compression of a digital photograph is achieved by Fourier transforming it, setting the high-momentum modes in the Fourier transform equal to zero, and then transforming back. However, as seen in Fig. 1.12(a), the straightforward application of a low-pass filter fails to reproduce nuclear phase shifts even at low energies.

This failure of a low-pass filter can be understood by recalling from Fig. 1.11 that the high and low-momentum modes are coupled. For instance, consider (schematically) the expression for the tangent of phase shift in perturbation theory


The second term in Eq. 1.14 involves a sum over off-diagonal matrix elements of . Therefore even phase shifts for small will have significant contributions from high if the coupling matrix elements are large.

A correct way to get rid of the high-momentum components is to use the RG evolution and lower the momentum cutoff in small steps.

(a) The running in .
(b) The SRG running in .
Figure 1.14: Schematic illustration of two types of RG evolution. For given or the matrix elements outside the corresponding lines are negligible. This decouples the high-momentum modes from the low-momentum ones. Figure from [44].

Two common choices for RG running are shown in Fig. 1.14. The RG running shown in Fig. 1.13(a) is referred to as “” was historically developed first for LENP. It attempts to get the potential in the low-pass-filter form through successive RG transformations [46, 47]. Though successful for two-nucleon forces, it has been difficult to systematically treat many-body forces in the framework.

A more recent approach through the similarity renormalization group (SRG) is illustrated in Fig. 1.13(b). SRG running drives the potential to a band-diagonal form. It is possible to systematically account for the many-body forces in the SRG framework. We will be using the SRG framework in our work presented in Chapter 3. Next we present a brief introduction to the SRG technique.

1.5.1 Srg

SRG for LENP was developed at Ohio State by S. Bogner (then a post-doc), R. Furnstahl, and R. Perry [48]. Their approach was inspired by the RG flows equations developed by Wegner for condensed matter applications [49]. The basic idea in SRG is to apply a series of unitary transformations to transform the Hamiltonian into a band-diagonal form shown in Fig. 1.13(b):


where is the RG flow parameter. . In practice, instead of using Eq. 1.15, the SRG evolution is done through the flow equation


is the operator which generates the flow. , the relative kinetic energy operator, is the most popular choice for , though other choices for have been explored [50].

Note that since we are just doing unitary transformations, eigenvalues such as energies are unchanged under SRG


It is beneficial to change the flow parameter to , where . has the units of momentum (). With this change of variable and setting , for a given partial-wave channel, the Eq. 1.16 in momentum basis becomes


The first term on the right side of Eq. 1.18 drives the potential to the band diagonal form shown in Fig. 1.13(b) and the second term on right side of Eq. 1.18 makes sure that the unitarity is maintained.

Figure 1.15: SRG evolution of the chiral potential [51] in the channel. Figure from [46].

Figure 1.15 shows Eq. 1.18 in action. We see that as we go to lower , we achieve the decoupling between high- and low-momentum components. Once this decoupling has been achieved, one can apply the low-pass filter and work with smaller matrices if desired.

Effects of SRG evolution
Figure 1.16: probability density in deuteron for the AV18 potential and the AV18 potential evolved to three SRG ’s. Figure from [45].

The conventional nuclear potentials have a strong short-range repulsion. As a result, the wave functions at short distance are suppressed. This suppression is seen for the unevolved potential in Fig. 1.16. This suppression is called the short-range correlation (SRC). The RG evolution gets rid of the short-range repulsion in the potential (it shifts the strength from high-momentum modes into the low-momentum modes). The wave functions from the evolved potential therefore do not exhibit the SRC as seen in Fig. 1.16. It is also noteworthy that even though the deuteron wave functions change at short distances, the long-distance part is independent of the SRG . The long-distance part of the wave function is related to the asymptotic normalization coefficient (ANC), which is an observable. We will see later in Chapter 2 that the ANC being an observable plays an important role in establishing the universality of the extrapolation formulas derived in Chapter 2. The observables should be independent of the resolution scale, which at the end of the day is a theoretical choice.

(a) Effect of evolution on binding energy, the -state probability, and the asymptotic -state ratio . Figure from [46].
(b) The solid black line indicate the phase shifts for the unevolved potential. The dotted blue line and the dashed red line are the phase shifts when the low-pass filter is applied to the unevolved potential and the evolved potential respectively.
Figure 1.17: Effect of RG evolution on (non-)observables.

In Fig. 1.17, we look at the effect of evolution on a few more (non-)observables. We see in Fig. 1.16(a) that the deuteron binding energy and the ratio of -state ANC to the -state ANC are independent of the resolution scale as expected. However, the -state probability depends on the resolution. This indicates that the -state probability of the deuteron is not an observable. We return to the phase shifts with a low-pass filter problem in Fig. 1.16(b). We see that when the low-pass filter is applied to the SRG evolved potential, it reproduces the unevolved phase shifts up to the cut-off energy of the filter. This discussion also makes it clear that the nuclear potential itself is not an observable and there is no “correct” nuclear potential. In fact, there is more than one “correct” potential. For instance, all the potentials in Fig. 1.15 give the same values for observables such as binding energies, phase shifts, etc.

Figure 1.18: SRG evolution greatly accelerates the convergence in many-body calculations. Figure from [45].

Since the SRG evolved potentials only have low-momentum modes, they are amenable to numerical calculations. As seen in Fig. 1.18, SRG evolution greatly accelerates the convergence in many-body calculations.

We note that even when we start with just two-body forces, SRG evolution introduces three- and higher-body forces. This can be understood from the flow equation (Eq. 1.16)


The second equality in Eq. 1.19 demonstrates how the commutators give rise to three- and higher-body (up to -body) forces. The initial potential in Fig. 1.18 includes both two- and three-body forces; it has been demonstrated that 3-body forces are crucial in getting the correct experimental values from theory [52]. To keep the invariance of energy with respect to the resolution scale (as in Fig. 1.16(a)), it is important to keep also the induced 3-body forces [27]. A major development in the SRG technology has been the ability to consistently evolve three-body forces [53, 54].

A related important development is that of In-Medium SRG (IM-SRG) [39]. IM-SRG uses a reference state which is different from the particle vacuum used in SRG. For example, can be a Slater determinant that is fair approximation to nucleus’ ground state. Just like in SRG, IM-SRG then uses a series of unitary transformations to decouple the reference state from excitations. IM-SRG also maintains the hierarchy of many-body forces, namely .

RG techniques have made possible calculations of medium-mass nuclei starting from inter-nucleonic interactions. However, if we keep pushing the calculations to higher-mass nuclei, we run into the same problem as indicated in Fig. 1.7(b), i.e, we run out of computational power before we reach convergence. So, along with RG techniques we also need reliable extrapolation techniques that will allow us to extrapolate the results from finite to . This problem will form the basis of Chapter 2.

1.6 Path forward for LENP

We have seen that nuclear theory has come a long way from the pion theories of the ’s. The focus these days is on doing precision calculations and making reliable predictions. The tool box of a nuclear theorist includes a wide variety of techniques. In most cases, different techniques have complimentary strengths. In other cases, alternative methods provide a means to cross-check answers. Having made great strides in the evaluation of the bound state properties, the push recently has been on calculating resonant and scattering states, and using the improved understanding of nuclear structure to study nuclear reactions.

To pick a particular example, let’s focus on the neutrinoless double beta decay () example. is of wide interest because of its potential to shed light on the nature of neutrinos (i.e, if neutrino is a Majorana or a Dirac fermion). has not been observed yet and there are experiments around the world looking for this decay [55]. The experiments need guidance from theory to help design the experiment and to interpret the measurements. The theoretical calculation of cross section involves computing a nuclear matrix element for the transition.

Figure 1.19: Nuclear matrix element for for various nuclei calculated using different methods. The legends indicate the different methods used for calculations. Figure from a talk by Petr Vogel.

Figure 1.19 indicates the current level of agreement for calculations from different methods. We see that these results differ by a factor of two. The reaction cross section is proportional to the matrix element squared, and therefore the current level of agreement translates to an uncertainty of factor four in the cost of the experiment.

We thus see that it is important that the theoretical calculations are accompanied by a reliable uncertainty estimate. The sources of uncertainty in a theoretical calculation can broadly be classified into two categories—uncertainty in the input (e.g., shortcomings in the potential, limitations of the assumptions made, etc.) and uncertainty arising from the method (e.g., errors due to various approximations made while solving). As presented in the next section, the work in this thesis focuses on investigating two common sources of uncertainties.

1.7 Thesis organization

The next two chapters in this thesis will present the author’s original work as a PhD student. As demonstrated in Fig. 1.7(b), we often run out of computational power before the convergence is reached, necessitating a need for trustworthy extrapolation schemes. Chapter 2 describes our published work in this regard [56, 57, 58]. The extrapolation schemes we developed were physically motivated as opposed to phenomenological forms previously in use. The author of this thesis was a lead author on Ref. [56]. All authors of Ref. [57] contributed equally. We will therefore look in detail at the results presented in Refs. [56, 57] in this thesis. The author’s contribution to work presented in Ref. [58] was secondary, and hence only the summary from that work will be presented.

Theoretical calculations of nuclear cross sections involve evaluating nuclear structure (which involves description of the initial and final state) and nuclear reaction (which involves a description of the probe). To make accurate predictions, it is important to understand the uncertainty stemming from the renormalization scale and scheme dependence of nuclear structure and reaction components. We address this problem in Chapter 3 by using the SRG to look at the scale dependence of deuteron electrodisintegration. This effort was led by the author and has been published in Ref. [59].

Both these projects were the first of its kind, and are fertile grounds for further development. In fact, the work presented in Chapter 2 has already had a lot of impact as testified by the number of articles citing our publications. We believe that the work in Chapter 3 will also receive wide attention soon. We present the details of our calculations in Appendix, which would enable anyone interested to reproduce our results.

Chapter 2 Extrapolation

As we have seen in the introduction, the harmonic oscillator (HO) basis is routinely used in low-energy nuclear physics (LENP) calculations. We also saw that the size of Hamiltonian matrix that we need to diagonalize grows factorially with the number of nucleons (cf. Fig. 1.9), severely restricting the number of terms that can be kept in the basis expansion. The single particle nuclear wave function with the truncation introduced in Chapter 1 is given by


in Eq. 2.1 are the HO wave functions; is the frequency of the HO 111 In LENP, the oscillator frequency is often denoted by rather than .. For reference, the -wave HO wave function is given by


where is the normalization constant, is the reduced mass, and denotes the generalized Laguerre polynomial.

The energy obtained in the HO basis——is a function of and . This is illustrated in Fig. 2.1.

Figure 2.1: Ground state energy for He as a function of and . Figure taken from [60].

We see that as we go to higher , the curves get flatter with respect to , or in other words the dependence on drops out.

The goal is to extrapolate to from a finite . The most widely used extrapolation scheme employs an exponential in form


where and are determined separately for each (with the option of constraining the fit to get the same asymptotic value). Figure 2.2 shows estimate for the ground state energy for He obtained using the extrapolation form of Eq. 2.3.

Figure 2.2: The estimate for exact He ground state energy using Eq. 2.3. Extrapolated answer from the constrained fit and the experimental binding energy are indicated by horizontal lines. Figure from [61].

The exponential in extrapolation is widely used in literature and seems to work quite well [62, 63, 64, 30, 31]. There are however many open questions about this extrapolation scheme such as the answer for depends on the oscillator frequency and it is not clear which is the best choice for . The terms and in Eq. 2.3 are fit to data. There is no way to extract these terms for one nucleus and use it to predict something else. Moreover, the physical motivation for an exponential in extrapolation is slim at best. It has been claimed that for larger nuclei is a logarithmic measure of the number of states [63]. This would account for the exponential behavior, but there is no demonstration that it follows in general or with a specific logarithmic dependence.

An alternative approach to extrapolations is motivated by effective field theory (EFT) and based instead on explicitly considering the infrared (IR) and ultraviolet (UV) cutoffs imposed by a finite oscillator basis [65]. The truncation in the oscillator basis introduces a maximum length scale (or an IR cutoff) and also a maximum momentum scale (or an UV cutoff). These length and momentum scales can be motivated by the classical turning points denoted by and respectively. We have


is the oscillator length given by . The errors due the finite IR (UV) cutoff are called the IR (UV) errors. To draw a lattice analogy, IR errors stem from the finite box size, and the UV errors are a result of the finite lattice spacing (or the graininess of the lattice). Ideally, we would like the box size to be as large as possible and the lattice spacing to be as small as possible. Because of finite computational power, this is not always possible though and therefore we need reliable extrapolation schemes in both IR as well as UV.

This approach of thinking of the HO truncation in terms of IR and UV cutoffs has led to lot of development in the past three years. Note that appears in the numerator for and in the denominator for , so it is not possible to make both cutoffs large simultaneously. We can choose the oscillator parameters such that one of the cutoffs in Eqs. 2.4 is large, making the errors due to that cutoff small and focus on the errors due to the other cutoff. The first attempt and test for a theoretically motivated IR correction was made in [60]. These corrections were made theoretically sound in [56, 57]. In [58] we looked at the UV correction for the deuteron. Our papers [56, 57, 58] have led to physically motivated extrapolation schemes and will form the basis of the next two sections. Insights from our work have also led to development of extrapolation schemes (both in IR and UV) for the many-body case by other groups. We will touch upon these developments in Subsec. 2.3.1.

2.1 Infrared story 222Based on [56] and [57]

As mentioned in the introduction, there was a lack of well motivated extrapolation schemes in LENP and this is where our work comes in. We started with the two-body case because it is more tractable mathematically. Note that extrapolation is usually not necessary for the two-body problem, because convergence is reached before we run out of computational power. This allows us to test our extrapolation formulas. Once we establish that the approach works for the two-body case, we can hope to extend the approach to few- and many-body case.

As mentioned earlier, IR cutoff effectively puts system in a finite box. We need to find appropriate box length such that


Note that our original problem was to find given at a finite . Once we make the correspondence in Eq. 2.5 and express energy as a function of box length, we can use various techniques (discussed later) to estimate which equals .

2.1.1 Tale of tails

The box size is usually bigger than the range of the potential. Thus imposing the IR cutoff modifies the asymptotic part (or tail) of the bound-state wave function. Our early work focused on trying to estimate the appropriate box size by matching the tails of wave functions in the truncated HO basis to the tails of wave functions in boxes of different lengths.

Our strategy was to use a range of model potentials for which the Schrödinger equation can be solved analytically or to any desired precision numerically to broadly test and illustrate various features, and then turn to the deuteron for a real-world example. In particular we considered:


where for each of the models we work in units with , reduced mass , and , but consider different values for . For the realistic potential we use the Entem-Machleidt 500 MeV chiral EFT NLO potential [66] and unitarily evolve it with the similarity renormalization group (SRG). These potentials provide a diverse set of tests for universal properties. Because we can go to very high  and for the two-particle bound states (and therefore large ), it is possible to always ensure that UV corrections are negligible.

We start with empirical considerations before presenting an analytical understanding. An example of how this correspondence between the HO truncation and a hard wall at specific length plays out is presented in Fig. 2.3.

Figure 2.3: (a) The exact radial wave function (dashed) for a square well Eq. 2.6 with depth (and ) is compared to the wave function obtained from an HO basis truncated at with (solid). The spatial extent of the wave function obtained from the HO basis truncation is dictated by the square of HO wave function for the highest radial quantum number (dot-dashed). (b) The wave functions obtained from imposing a Dirichlet boundary condition (bc) at , and are compared to the wave function in truncated HO basis.

In the top panel, the exact ground-state radial wave function (dashed) for the square well in Eq. 2.6 is compared to the solution in an oscillator basis truncated at determined by diagonalization (solid). The truncated basis cuts off the tail of the exact wave function because the individual basis wave functions have a radial extent that depends on  (from the Gaussian part; cf. Eq. 2.2) and on the largest power of (from the polynomial part). The latter is given by . With and , this means that gives the largest power.

The cutoff will then be determined by the oscillator wave function, , whose square (which is the relevant quantity) is also plotted in the top panel (dot-dashed). It is evident that the tail of the wave function in the truncated basis is fixed by this squared wave function. Our premise is that the HO truncation is well modeled by a hard-wall (Dirichlet) bc at . If so, the question remains how best to quantitatively determine given and . In the bottom panel of Fig. 2.3 we show the wave functions for several possible choices for . corresponds to choosing the classical turning point (i.e. the half-height point of the tail of ); it is manifestly too small. The authors of [60] advocated an improved choice for given by


The length in Eq. 2.10 is obtained by linear extrapolation from the slope at the half-height point. However, choosing


was found to work the best in almost all examples.

The most direct illustration of this conclusion comes from the bound-state energies. In the example in Fig. 2.3, the exact energy (in dimensionless units) is while the result for the basis truncated at is , which is therefore what we hope to reproduce. With , the energy is , with it is , and with it is . While this is only one example of a model problem, we have found that always gives a better energy estimate than (and something like is almost always worse).

Another signature that demonstrates the suitability of is that points from many different and values all lie on the same curve. Figures 2.4 and 2.5 show the energies from a wide range of HO truncations for , and for the Gaussian well and the square well potential, respectively. The energies for different and lie on the same smooth and unbroken curve if we use but not with the other choices. For and , one finds that sets of points with different but same fall on smooth, -dependent curves. For the square well, there are small discontinuities visible even for . At the square well radius, the wave function’s second derivative is not smooth, and this is difficult to approximate with a finite set of oscillator functions. This lack of UV convergence is likely the origin of the very small discontinuities.

Figure 2.4: Ground-state energies versus (top), (middle), and (bottom) for a Gaussian potential well Eq. 2.8 with and . The crosses are the energies from HO basis truncation. The energies obtained by numerically solving the Schrödinger equation with a Dirichlet bc at lie on the solid line. The horizontal dotted lines mark the exact energy .
Figure 2.5: Ground-state energies versus (top), (middle), and (bottom) for a square well potential well Eq. 2.6 with and . The crosses are the energies from HO basis truncation. The energies obtained by numerically solving the Schrödinger equation with a Dirichlet bc at lie on the solid line. The horizontal dotted lines mark the exact energy .

As a further test, we solve the Schrödinger equation with a vanishing Dirichlet bc (solid lines in Figs. 2.4 and 2.5), and compare to the energies obtained from the HO truncations (crosses). The finite oscillator basis energies are well approximated by a Dirichlet bc with a mapping from the oscillator and to an equivalent length given by . Note that for large , the differences between , and may be smaller than other uncertainties involved in nuclear calculations, but for practical calculations one will want to use small results, where these considerations are very relevant.

These results from model calculations are consistent with those from realistic potentials applied to the deuteron. To illustrate this, we use the NLO 500 MeV potential of Entem and Machleidt [66]. We generate results in an HO basis with  ranging from to and from to (in steps of 4 to avoid HO artifacts for the deuteron [63]). We then restrict the data to where UV corrections are negligible (see Section 2.2).

Figure 2.6: Ground-state energies versus (top), (middle), and (bottom) for the Entem-Machleidt 500 MeV NLO potential [66]. The horizontal dotted lines mark the exact energy .

Figure 2.6 shows that the criterion of a continuous curve with the smallest spread of points clearly favors .

Analytical derivation of

In the asymptotic region (this is the region where IR cutoff is imposed), the potential is negligible and the only relevant part of the Hamiltonian is the kinetic energy or the operator. In what follows, we analytically compute the smallest eigenvalue of in a finite oscillator basis and will see that . In the remainder of this Subsection, we set the oscillator length to one. We focus on -waves and thus consider wave functions that are regular at the origin, i.e. the radial wave functions are identical to the odd wave functions of the one-dimensional harmonic oscillator.

The localized eigenfunction of the operator with smallest eigenvalue is


We employ the -wave oscillator functions


with energy . Here, denotes the Laguerre polynomial, and it is convenient to rewrite this function in terms of the Hermite polynomial . We expand the eigenfunction in Eq. 2.12 as


Before we turn to the computation of the expansion coefficients , we consider the eigenvalue problem for the operator . We have


where and denote the annihilation and creation operator for the one-dimensional harmonic oscillator, respectively. The matrix of is tridiagonal in the oscillator basis. For the matrix representation, we order the basis states as . Thus, the eigenvalue problem becomes a set of rows of coupled linear equations. In an infinite basis, the eigenvector identically satisfies every row of the eigenvalue problem for any value of . In a finite basis , however, the last row of the eigenvalue problem


can only be fulfilled for certain values of , and this is the quantization condition. To solve this eigenvalue problem we need expressions for the expansion coefficients for . Those can be derived analytically as follows.

We rewrite the eigenfunction in Eq. 2.12 as a Fourier transform


and expand the sine function in terms of oscillator functions as


Thus, the expansion coefficients in Eq. 2.14 are given in terms of the Fourier transform as


So far, all manipulations have been exact. We need an expression for for and use the asymptotic expansion


which is valid for , see [67]. Using this approximation, one finds (making use of Fourier transforms)


with due to Eq. 2.12.

Let us return to the solution of the quantization condition in Eq. 2.16. We make the ansatz


and must assume that . This ansatz is well motivated, since the naive semiclassical estimate yields . We insert the expansion coefficients of Eq. 2.21 into Eq. 2.16 and consider its leading-order approximation for and . This yields


as the solution. Recalling that a truncation of the basis at corresponds to the maximum energy , we see that we must identify . Thus, is the lowest momentum in a finite oscillator basis with basis states (and not as stated in Ref. [65]). It is clear from its very definition that is also (a very precise approximation of) the infrared cutoff in a finite oscillator basis, and that (and not as stated in Refs. [68, 69]) is the radial extent of the oscillator basis and the analog to the extent of the lattice in the lattice computations [70].

The derivation of our key result is based on the assumption that the number of shells fulfills . Table 2.1 shows a comparison of numerical results for in different model spaces. We see that is a very good approximation already for , with a deviation of about 1%.

0 1.2247 1.1874 1.8138
2 0.9586 0.9472 1.1874
4 0.8163 0.8112 0.9472
6 0.7236 0.7207 0.8112
8 0.6568 0.6551 0.7207
10 0.6058 0.6046 0.6551
12 0.5651 0.5642 0.6046
14 0.5316 0.5310 0.5642
16 0.5035 0.5031 0.5310
18 0.4795 0.4791 0.5031
20 0.4585 0.4582 0.4791
Table 2.1: Comparison between the lowest momentum , , and for model spaces with up to oscillator quanta.

Note that this approach can be generalized to other localized bases. The (numerical) computation of the lowest eigenvalue of the momentum operator yields the box size corresponding to the employed Hilbert space.

EFT-like approach

We mentioned that the relevant operator for IR truncation is . To get a better understanding of the correspondence between the HO truncation and a hard wall at , let’s compare the spectrum of in the two cases.

Figure 2.7: Eigenfunctions of in the truncated HO basis compared to those in a box of size .

In Fig. 2.7, we compare the low-lying eigenfunctions of in the truncated HO basis to the eigenfunctions in a box of size . As we will see later, the asymptotic or near the wall difference between the two eigenfunctions are high-momentum effects irrelevant for the long-wavelength physics of the bound states.

Another way to look at this is to compute the number of (-wave) states up to a momentum . We find


Here, we apply the semiclassical approximation and write the trace as a phase-space integral. We assume , perform the integrations and use 333For the sake of brevity we replace by . This yields


where is the oscillator length. Figure 2.8 shows a comparison between the quantum mechanical staircase function and the semiclassical estimate of Eq. 2.25 for . For sufficiently small values of , the number of -wave momentum eigenstates grows linearly, and inspection of Eq. 2.25 shows that the slope at the origin is semiclassically. The linear growth of the number of eigenstates of with clearly demonstrate that — at not too large values of — the spectrum of in the oscillator basis is similar to the spectrum of in a spherical box.

Figure 2.8: The staircase function of the states of the operator in a finite oscillator basis with (black) compared to its semiclassical estimate (smooth red curve). denotes the number of states of the operator with eigenvalues .

As a final example of the correspondence between the HO truncation and the hard wall at , we look at the ground state wave functions of a square well in the two bases in Fig. 2.9.

Figure 2.9: Ground-state wave functions for a square well potential of depth (see Eq. 2.6; lengths are in units of and energies in units of with ) from solving the Schrödinger equation with a truncated harmonic oscillator basis with and (dashed) and with a Dirichlet bc at given from Eq. 2.113 (solid). The coordinate-space radial wave functions in a) exhibit a difference at near 1.5, but the Fourier-transformed wave functions in b) are in close agreement at low , showing that the differences are high-momentum modes.

The binding momentum in this case is (in units of ). The Fourier-transformed wave functions differ at much larger momentum and this difference is irrelevant for the long-wavelength physics of bound states. Thus the use of Dirichlet bc to take into account the HO truncation is similar in spirit to the use of contact interactions to describe the effect of unknown short-ranged forces on long-wavelength probes.

2.1.2 Cashing in on the hard wall correspondence

We have so far focused on establishing how HO truncation is analogous to putting the system in a spherical box of a specified radius. Now let’s see how this correspondence helps us in getting the exact energy .

Linear energy method

Our first approximation to the IR correction is based on what is known in quantum chemistry as the linear energy method [71]. Given a hard-wall bc at beyond the range of the potential, we write the energy compared to that for as


We seek an estimate for , which is assumed to be small, based on an expansion of the wave function in . Let be a radial solution with regular bc at the origin and energy . For convenience in using standard quantum scattering formalism below, we choose the normalization corresponding to what is called the “regular solution” in Ref. [72], which means that and the slope at the origin is unity for all . We denote the particular solutions and . Then there is a smooth expansion of about at fixed , so we approximate [71]


for . By evaluating Eq. 2.27 at with the bc , we find


which is the estimate for the IR correction.

Figure 2.10: Testing the linear energy approximation Eq. 2.27 for (a) deep () and (b) shallow () Gaussian potential well Eq. 2.8 (). The solid lines are the exact solutions for energies and , respectively, whose zero crossings determine the corresponding values for .

We can check the accuracy of the linear energy approximation (Eq. 2.27) by numerically solving the Schrödinger equation with a specified energy. This determines as the radius at which the resulting wave function vanishes. Then we compare this wave function for to the right side of Eq. 2.27, with the derivative calculated numerically. Figure 2.10 shows representative examples for a deep and shallow Gaussian potential. In these examples and other cases, the approximation to the wave function is good, particularly in the interior. The estimates for using the right side of Eq. 2.28 are within a few to ten percent: 0.68 versus 0.70 and 0.050 versus 0.055 for the two cases.

The good approximation to the wave function suggests that for the calculation of other observables the linear energy approximation will be useful. For observables most sensitive to the long distance (outer) part of the wave function, such as the radius, this has already been shown to be true [73]. But the good approximation to the wave function at small means that corrections for short-range observables should also be controlled, with the dominant contribution in an extrapolation formula coming from the normalization.

Next we derive an expression for the derivative in Eq. 2.28. To start with we assume we have a single partial-wave channel. For general , the asymptotic form of the radial wave function for greater than the range of the potential is


with for . We take the derivative of Eq. 2.29 with respect to energy, evaluate at using and , to find


We now evaluate at and anticipate that the term dominates:


Substituting Eq. 2.31 into Eq. 2.28, we obtain


Note that this result is independent of the normalization of the wave function.

To calculate the derivative explicitly, we turn to scattering theory, following the notation and discussion in Ref. [72]. In particular, the asymptotic form of the regular scattering wave function for orbital angular momentum and for positive energy is given in terms of the Jost function  [72],


where the functions (related to Hankel functions) behave asymptotically as


The ratio of the Jost functions appearing in Eq. 2.33 gives the partial wave -matrix :


which is in turn related to the partial-wave scattering amplitude by


We will restrict ourselves to for simplicity; the generalization to higher is straightforward and will be considered later.

To apply Eq. 2.33 to negative energies, we analytically continue from real to (positive) imaginary . So,


where is the range of the potential. Upon comparing to Eq. 2.29 we conclude that