Coulomb expectation values in $D=3$ and $D=3-2\epsilon$ dimensions

TL;DR

This paper provides detailed methods and formulas for calculating expectation values in Coulomb bound states across three and near-three dimensions, aiding dimensional regularization in quantum electrodynamics corrections.

Contribution

It introduces new techniques and formulas for evaluating expectation values in D=3 and D=3-2ε dimensions, including divergent and finite cases, enhancing dimensional regularization methods.

Findings

Derived formulas for expectation values involving Laguerre polynomials.

Presented perturbative expansions in ε for bound state energies.

Provided regularized expressions for divergent expectation values.

Abstract

We explore the quantum Coulomb problem for two-body bound states, in $D=3$ and $D=3-2\epsilon$ dimensions, in detail, and give an extensive list of expectation values that arise in the evaluation of QED corrections to bound state energies. We describe the techniques used to obtain these expectation values and give general formulas for the evaluation of integrals involving associated Laguerre polynomials. In addition, we give formulas for the evaluation of integrals involving subtracted associated Laguerre polynomials--those with low powers of the variable subtracted off--that arise when evaluating divergent expectation values. We present perturbative results (in the parameter $\epsilon$) that show how bound state energies and wave functions in $D=3-2\epsilon$ dimensions differ from their $D=3$ dimensional counterparts and use these formulas to find regularized expressions for divergent expectation values such as $\big \langle \bar V^3 \big \rangle$ and $\big \langle (\bar V')^2 \big \rangle$ where $\bar V$ is the $D$-dimensional Coulomb potential. We evaluate a number of finite $D$-dimensional expectation values such as $\big \langle r^{-2+4\epsilon} \partial_r^2 \big \rangle$ and $\big \langle r^{4\epsilon} p^4 \big \rangle$ that have $\epsilon \rightarrow 0$ limits that differ from their three-dimensional counterparts $\big \langle r^{-2} \partial_r^2 \big \rangle$ and $\big \langle p^4 \big \rangle$. We explore the use of recursion relations, the Feynman-Hellmann theorem, and momentum space brackets combined with $D$-dimensional Fourier transformation for the evaluation of $D$-dimensional expectation values. The results of this paper are useful when using dimensional regularization in the calculation of properties of Coulomb bound systems.

Equations

$$\bar{V}(\vec{k})=\frac{-4\pi Z\alpha {\bar{\mu }}^{2ϵ}}{\vec{k}{}^2},$$

$${\overline{\mu }}^2={\mu }_{\overline{MS}}^2=\frac{{\mu }^2{e}^{{\gamma }_E}}{4\pi },$$

$$\bar{V}(r)=\int {\text{mathchar}22d}^{D}k{e}^{i\vec{k}\cdot \vec{x}}\bar{V}(\vec{k})=-\frac{\Gamma \left(D/2-1\right){\bar{\mu }}^{2ϵ}Z\alpha }{{\pi }^{D/2-1}{r}^{D-2}}=-\frac{\beta }{r^{1-2ϵ}},$$

$$-\vec{\nabla }{}^2\bar{V}(r)=-4\pi Z\alpha {\bar{\mu }}^{2ϵ}{\delta }^D(\vec{x}).$$

$$\frac{\vec{p}{}^2}{2m_r}\psi (\vec{x})+\bar{V}(r)\psi (\vec{x})=\bar{E}\psi (\vec{x})$$

$$\psi (\vec{x})=\bar{R}(r)Y(\hat{x})$$

$${\vec{\nabla }}^2={\partial }_r^2+\frac{D-1}{r}{\partial }_r-\frac{L^2}{r^2}$$

$$L^2=\sum _{i<j}{L}_{ij}^2\text{with}{L}_{ij}=-i\left(x_i{\partial }_j-x_j{\partial }_i\right).$$

$$\frac{1}{2m_r}\left\{-{\partial }_r^2-\frac{D-1}{r}{\partial }_r+\frac{\ell (\ell +D-2)}{r^2}\right\}{\bar{R}}_{n\ell }(r)+\bar{V}(r){\bar{R}}_{n\ell }(r)={\bar{E}}_{n\ell }{\bar{R}}_{n\ell }(r).$$

$$\left\{{\partial }_r^2+\frac{D-1}{r}{\partial }_r-\frac{\ell (\ell +D-2)}{r^2}\right\}{\bar{R}}_{n\ell }(r)=0.$$

$$m(m-1)+(D-1)m-\ell (\ell +D-2)=(m-\ell )(m+\ell +D-2)=0,$$

$${\bar{R}}_{n\ell }(r)\to r^{\ell }.$$

$$\frac{1}{2m_r}\left\{-{\partial }_r^2\right\}{\bar{R}}_{n\ell }(r)=E_{n\ell }{\bar{R}}_{n\ell }(r);$$

$${\bar{R}}_{n\ell }(r)\to e^{-{\bar{\gamma }}_{n\ell }r}$$

$${\bar{E}}_{n\ell }=-\frac{{\bar{\gamma }}_{n\ell }^2}{2m_r},{\bar{\gamma }}_{n\ell }=\sqrt{-2m_r{\bar{E}}_{n\ell }}.$$

$${\bar{R}}_{n\ell }(r)={\bar{\varphi }}_{n\ell }{\Omega }_{D-1}^{1/2}{\left(\frac{(n+\ell )!}{n(n-\ell -1)!}\right)}^{1/2}\frac{{\rho }^{\ell }}{(2\ell +1)!}{e}^{-\rho /2}{L}_{n\ell }(\rho )$$

$${\bar{\varphi }}_n\equiv {\bar{\varphi }}_{n0}=\underset{r\to 0}{\lim }{\bar{\psi }}_{n0}(r)$$

$${\Omega }_n=\frac{2{\pi }^{\frac{n+1}{2}}}{\Gamma \left(\frac{n+1}{2}\right)}$$

$$\underset{\rho \to 0}{\lim }{L}_{n\ell }(\rho )=1.$$

$$\underset{ϵ\to 0}{\lim }{L}_{n\ell }(\rho )={\left(\genfrac{}{}{0ex}{}{n+\ell }{2\ell +1}\right)}^{-1}{L}_{n-\ell -1}^{2\ell +1}(\rho )\text{and}\underset{ϵ\to 0}{\lim }{\bar{\varphi }}_{n\ell }={\left(\frac{(m_{r}Z\alpha {)}^3}{\pi n^3}\right)}^{1/2}.$$

$$\left\{{\partial }_r^2+\left(\frac{2\ell +D-1}{r}-2\bar{\gamma }\right){\partial }_r-\frac{(2\ell +D-1)\bar{\gamma }}{r}-2m_r\bar{V}(r)\right\}{L}_{n\ell }(\rho )=0.$$

$$\frac{-2m_r}{(2{\bar{\gamma }}_{n\ell }{)}^2}\bar{V}(r)=\frac{2m_r\Gamma (1/2-ϵ){\bar{\mu }}^{2ϵ}Z\alpha }{(2{\bar{\gamma }}_{n\ell }){\pi }^{1/2-ϵ}(2{\bar{\gamma }}_{n\ell }r)}{\left(\frac{\rho }{2{\bar{\gamma }}_{n\ell }}\right)}^{2ϵ}=\frac{{\bar{n}}_{n\ell }{\rho }^{2ϵ}}{\rho },$$

$${\bar{n}}_{n\ell }\equiv \frac{m_{r}Z\alpha }{{\bar{\gamma }}_{n\ell }}\frac{\Gamma (1/2-ϵ)}{{\pi }^{1/2-ϵ}}{\left(\frac{\bar{\mu }}{2{\bar{\gamma }}_{n\ell }}\right)}^{2ϵ},$$

$$\left\{{\partial }_{\rho }^2+\left(\frac{2(\ell +1-ϵ)}{\rho }-1\right){\partial }_{\rho }-\frac{\ell +1-ϵ}{\rho }+\frac{{\bar{n}}_{n\ell }{\rho }^{2ϵ}}{\rho }\right\}{L}_{n\ell }(\rho )=0.$$

$$L_{n\ell }(\rho )$$

$$a_{jk}=\frac{a_{j-1,k}(j+\ell +ϵ(2k-1))-a_{j-1,k-1}}{(j+2ϵk)(j+2\ell +1+2ϵ(k-1))}.$$

$$A_j=\sum _{k=0}^j{a}_{jk}{n}^k=\frac{(-1{)}^j{\Pi }_{k=1}^j(n-\ell -k)}{j!{\Pi }_{k=1}^j(2\ell +1+k)}=\frac{(-1{)}^j(n-\ell -1)!(2\ell +1)!}{j!(n-\ell -j-1)!(2\ell +1+j)!},$$

$$H^{\prime}=\bar{V}(r)-\tilde{V}(r)=-\frac{2Z\alpha}{r}\big{(}\ln(\mu r)+\gamma_{E}\big{)}\epsilon+O(\epsilon^{2}),$$

$${\bar{E}}_{n\ell }=E_n\left\{1+ϵ\left[4\ln \left(\frac{\mu n}{2m_{r}Z\alpha }\right)+4{\text{H}}_{n+\ell }+\frac{2}{n}\right]+O({ϵ}^2)\right\},$$

$${\bar{\varphi }}_{n\ell }$$

Full text

I Introduction

The Coulomb force provides for the binding of all atoms. The nature of the Coulomb force is most clearly displayed in the structure of two-body atoms such as hydrogen, hydrogen-like ions (He+, Li*++*, etc.), positronium ($e^{+}e^{-}$), muonium ($\mu^{+}e^{-}$), true-muonium ($\mu^{+}\mu^{-}$), muonic hydrogen ($p^{+}\mu^{-}$), and many more. Because of its great importance, the quantum mechanics of the two-body Coulomb problem has been extensively studied. In the course of the calculation of relativistic and QFT corrections to the energy levels of two-body Coulombic bound states, expectation values of various operators in quantum bound states are eventually required. Many of these expectation values, for instance that of $1/r^{3}$, are singular due to short-distance divergences, and some regularization is required. An extremely convenient scheme for obtaining finite quantities is provided by dimensional regularization. In this work we discuss the dimensionally regularized Schrödinger equation for the Coulomb potential and show how the $D$-dimensional solution can be used to find regularized values for otherwise divergent quantities. Another goal of this work is to present an extensive list of expectation values–both finite and divergent–that actually arise in practical calculations. To obtain the finite ones, we developed a number of general integral formulas involving associated Laguerre polynomials, logarithms, and powers, and explored their consequences. To obtain the divergent ones, we required a detailed understanding of the short-distance behavior of the $D$-dimensional wave functions.

This paper is organized as follows. In Sec. II we discuss the quantum mechanical Coulomb problem in $D=3-2\epsilon$ dimensions, showing how to find the wave functions and energy levels and exploring their behavior near $\epsilon=0$. In Sec. III we explore the general problem of finding expectation values in $D$ dimensions with a focus on divergent expectation values for which dimensional regularization is used. We give an example of the evaluation of a divergent expectation value and show when the $\epsilon\rightarrow 0$ limit of a finite expectation value calculated in $D=3-2\epsilon$ dimensions might be different from the same finite expectation value calculated directly in three dimensions. We also discuss some implications of the $D$-dimensional recursion relations and the Feynman-Hellmann theorem. In Sec. IV we illustrate the use of divergent Coulomb expectation values in a practical calculation. In Sec. V we develop general formulas for the evaluation of integrals involving one or two associated Laguerre polynomials of the sort that arise in calculations involving Coulomb bound states. Sec. VI is a brief conclusion. The main results of this paper, the expectation values themselves, appear in Appendix A, along with some of the closely related momentum space brackets. The formula for $D$-dimensional Fourier transforms is given in Appendix B. We have included a review of the three-dimensional quantum Coulomb problem in Appendix C. This has a focus on Coulomb bound states, and is included since the notation for associated Laguerre polynomials is not standardized, and because it is helpful to have a listing of useful results in one place and in a consistent notation. Properties of the special functions used in this work are given in Appendix D. Useful formulas for gamma, polygamma, zeta, beta, and hypergeometric functions are included, as is a discussion of the harmonic and “diharmonic” numbers that appear in quantum Coulomb calculations. Finally, we have given in Appendix E a discussion of a specific momentum space bracket, \big{\langle}\!\ln q\big{\rangle}_{\vec{p}_{2},\vec{p}_{1}}, that is needed in calculations but is a special case.

Throughout this paper we use units where $\hbar=c=1$. Standard units can be recovered at every stage by inserting the factors of $\hbar$ and $c$ needed to restore conventional dimensions.

The paper contains a mixture of new, known but obscure, and well-known results. Results of all types are included with the goal of giving a unified and comprehensive presentation and useful tables. Most of the $D=3-2\epsilon$ results are new although many of the momentum space brackets (see Sec.  III) for $n=1$ were obtained in Czarnecki et al. (1999a). The application of the generalized power series solution for the $D=3-2\epsilon$ wave function was introduced in Adkins (2018a) but is applied more systematically here. Section V involves integrals of associated Laguerre polynomials times powers of $x$ and $\ln x$. Various expressions for the basic integral $\int_{0}^{\infty}dx\,e^{-x}x^{s}L_{n}^{k}(x)L_{n^{\prime}}^{k^{\prime}}(x)$ have appeared in the literature. Our formula for this integral differs from others we have seen and is quite convenient for our purposes. General formulas for the integrals involving logarithms are, to our knowledge, new, as are the formulas for integrals involving subtracted associated Laguerre polynomials. We have not seen the summation formulas derived in Section V elsewhere. The definition and derived properties of the diharmonic sums given in Appendix D are, we believe, new.

II The quantum Coulomb problem in $D=3-2\epsilon$ dimensions

In this section we review the quantum mechanics of the Coulomb problem in $D$ dimensions. Dimensional regularization has been extremely useful for atomic physics calculations Brown (1992); Pineda and Soto (1998a, b, c); Czarnecki et al. (1999b, a); Czarnecki et al. (2005); Jentschura et al. (2005); Wundt and Jentschura (2008); Adkins (2018a, b) and for quark-antiquark bound states involving heavy quarks (reviewed, for example in Hoang (2002); Brambilla et al. (2005); Pineda (2012)). The dimensionally regularized Coulomb problem serves as an extension of the usual three-dimensional case that allows us to evaluate quantities that would be divergent in three dimensions. Central to our discussion is the Coulomb potential energy between an electron and a nucleus consisting of $Z$ positive charges. In momentum space, this interaction energy is

[TABLE]

which can be obtained directly from the Coulomb-gauge Feynman rules for the interaction between the two charges with the exchange of a Coulomb photon. The $\overline{MS}$ mass parameter that was introduced to keep $\alpha$ dimensionless can be written as

[TABLE]

where $\mu$ is another parameter with units of mass and $\gamma_{E}$ is the Euler-Mascheroni constant. The corresponding coordinate space potential energy is

[TABLE]

where ${\mathchar 22\relax\mkern-11.0mud}^{n}k$ is the $n$-dimensional momentum space integration measure ${\mathchar 22\relax\mkern-11.0mud}^{n}k=\frac{d^{n}k}{(2\pi)^{n}}$ and $\beta=\Gamma(D/2-1)\bar{\mu}^{2\epsilon}Z\alpha/\pi^{D/2-1}$. (For the $D$-dimensional Fourier transform, see Samko (1978) and Appendix B.) This potential energy satisfies the $D$-dimensional Gauss’ Law

[TABLE]

The Schrödinger-Coulomb equation in $D$ dimensions is Andrew and Supplee (1990); Morales (1996); Dong (2011); Bureš (2015)

[TABLE]

where $m_{r}=m_{1}m_{2}/(m_{1}+m_{2})$ is the reduced mass. Since the Coulomb potential has no angular dependence, it is possible to separate radial from angular variables as in three dimensions. We write

[TABLE]

and use

[TABLE]

where

[TABLE]

We choose $Y_{\ell}(\hat{x})$ to be an eigenfunction of $L^{2}$ with eigenvalue $\ell(\ell+D-2)$ Avery (1989, 2010). For example, the first few of these $D$-dimensional angular eigenfunctions are $Y_{0}(\hat{x})\propto 1$, $Y_{1,i}(\hat{x})\propto\hat{x}_{i}$, $Y_{2,ij}(\hat{x})\propto\hat{x}_{i}\hat{x}_{j}-\frac{1}{D}\delta_{ij}$. Consequently, the radial function satisfies

[TABLE]

The quantum numbers listed for $\bar{R}$ include the angular momentum quantum number $\ell$ and a principal quantum number $n$ that is an integer $n=1$, 2, 3, etc. labeling the states with a given value of $\ell$ ordered by energy. The energy $\bar{E}_{n\ell}$ for $D\neq 3$ depends on $\ell$ as well as $n$.

It is useful to factor out the leading short and long distance behavior of the wave function as in three dimensions. The leading short distance terms are all part of $\vec{\nabla}^{2}$. The $1/r^{2}$ term in $\vec{\nabla}^{2}$ dominates the $1/r^{D-2}$ of the potential and the constant in the energy term. So, at short distances, $\bar{R}_{n\ell}$ satisfies

[TABLE]

This differential equation is of Euler (or Cauchy-Euler) form and is satisfied by a power $\bar{R}_{n\ell}(r)=r^{m}$. Insertion of this trial solution into the equation gives

[TABLE]

so $m=\ell$ or $m=-(\ell+D-2)$. For $\ell\geq 0$ and $D$ near three the second solution is negative and consequently unacceptable. The implied behavior of the radial function near $r=0$ is

[TABLE]

For long distances the leading terms in (II.9) are the first and the last:

[TABLE]

the solution that falls to zero at large distances is

[TABLE]

where

[TABLE]

We define a new function $L_{n\ell}(\rho)$ according to

[TABLE]

where $2\bar{\gamma}_{n\ell}r\equiv\rho$, $\bar{\phi}_{n\ell}$ is a normalization constant such that

[TABLE]

is the S-state wave function at the origin (at “contact”), $\bar{\phi}_{n\ell}$ with $\ell>0$ is related to the $\ell^{\text{th}}$ derivative of the wave function at contact, and

[TABLE]

is the surface area of a unit $n$-sphere. That is, $\Omega_{D-1}$ is the total solid angle of $D$-dimensional space. The new function $L_{n\ell}(\rho)$ is the analog of the associated Laguerre polynomial $L_{n-\ell-1}^{2\ell+1}(\rho)$ but is normalized differently:

[TABLE]

In the three-dimensional limit one has

[TABLE]

When expressed in terms of $L_{n\ell}(\rho)$, the radial equation (II.9) becomes

[TABLE]

It seems natural to change the independent variable from $r$ to the dimensionless variable $\rho$. In doing so we will be dividing equation (II.21) by $(2\bar{\gamma})^{2}$. The potential term becomes

[TABLE]

where

[TABLE]

which would be an integer $n$ with $n\geq 1$ in three dimensions. The radial equation with a dimensionless independent variable becomes

[TABLE]

A bit of trial and error reveals that a power series expansion around the origin is more complicated than usual and takes the form

[TABLE]

with $a_{00}=1$ and $a_{jk}=0$ unless $0\leq k\leq j$. We insert $L_{n\ell}(\rho)$ into (II.24) and assume that all powers $\rho^{j+2\epsilon k}$ are independent. We isolate the coefficient of $\rho^{j-2+2\epsilon k}$, which is required to vanish, and obtain the following recursion relation

[TABLE]

The recursion relation, along with the boundary conditions listed above, allows us to find all of the coefficients $a_{jk}$. For example, we find that $a_{10}=\frac{1}{2}$, $a_{11}=\frac{-1}{2(1+\ell)(1+2\epsilon)}$, etc. It would be easy to write a routine to calculate as many of these coefficients as desired.

As a quick consistency check we consider the three-dimensional limit and calculate the total coefficient $A_{j}$ of a given power $\rho^{j}$. This coefficient is

[TABLE]

where we have used $\epsilon\rightarrow 0$, $\bar{n}_{n\ell}\rightarrow n$. We see that when $n=\ell+k$ where $k$ is a positive integer, the series terminates with a maximum power of $\rho^{n-\ell-1}$. By comparison with (C.20), we see that $L_{n\ell}$ in this limit is just proportional to the usual associated Laguerre polynomials.

The differential equation (II.24) can be solved numerically (as in Adkins (2018a)) to find the wave function and eigenvalues $\bar{n}_{n\ell}$ (which give values for $\bar{\gamma}_{n\ell}$ using (II.23) and then for the energies using (II.15)), or a perturbative approach can be employed. The perturbative method employs a $D$-dimensional lowest-order problem with potential $\tilde{V}(r)=-Z\alpha/r$, for which an exact solution is known Alliluev (1958); Nieto (1979). Building a perturbation scheme on the basis of this lowest order problem, with

[TABLE]

one can obtain expressions for the energies

[TABLE]

where $E_{n}=-\frac{m_{r}(Z\alpha)^{2}}{2n^{2}}$ represents the usual Bohr energies, and for the “wave function at the origin” factor

[TABLE]

where $n_{r}=n-\ell-1$ is the radial quantum number, $\gamma_{n}\equiv\frac{m_{r}Z\alpha}{n}$, and the harmonic, generalized harmonic, and diharmonic numbers $\text{H}_{n}$, $\text{H}^{(2)}_{n}$, and $\text{diH}_{\pm}(n,m)$ are defined in Appendix D.

III Coulomb expectation values in $D=3$ and $D=3-2\epsilon$ dimensions

Expectation values in $D$ dimensions are defined in the usual way:

[TABLE]

where $M$ is an operator that may involve derivatives as well as functions of position. Many of the expectation values we will need are finite in three dimensions, but others, such as \big{\langle}\bar{V}^{3}\big{\rangle}\propto\big{\langle}r^{-3+6\epsilon}\big{\rangle} and \big{\langle}(\bar{V}^{\prime})^{2}\big{\rangle}\propto\big{\langle}r^{-4+4\epsilon}\big{\rangle}, are divergent for S states. These divergences are controlled by dimensional regularization, and working in $D$ dimensions we can find their explicit values. In this section we will give some details of the evaluation of Coulomb expectation values–with a specific focus on the evaluation of divergent ones. Situations can also arise for finite expectation values where taking the $D\rightarrow 3$ limit and taking the expectation value do not commute:

[TABLE]

We will give an example of this subtlety and describe when it will occur. Some special techniques are useful for finding and/or testing expectation values, such as recursion relations and the Feynman-Hellmann theorem. We will explore some consequences of these methods in three and $D$ dimensions. Finally, we will discuss momentum space brackets and their relation to expectation values.

Many innovative procedures for the evaluation of Coulomb expectation values have appeared in the literature, including those of Schrödinger (1926); Waller (1926); Van Vleck (1934); Pasternack (1937); Bethe and Salpeter (1957); Epstein and Epstein (1962); Srivastava et al. (2003); Suslov and Trey (2008); Gonzalez et al. (2017), as well as useful lists of expectation values: Swainson and Drake (1991); Titard and Ynduráin (1994); Czarnecki et al. (1999a); Zatorski (2008). Our lists overlap with these but are more extensive.

The expectation value \big{\langle}(\bar{V}^{\prime})^{2}\big{\rangle}_{n\ell} is divergent for $\ell=0$ with a divergence proportional to \big{\langle}r^{-4}\big{\rangle}_{n0} in three dimensions, so it must be regularized in order to be evaluated. When $\ell>0$ this expectation value is finite due to the extra factors of $r^{\ell}$ in the wave function, and we can use standard three-dimensional methods as described in Sec. V. After performing the angular integrals, the S state expectation value can be written as

[TABLE]

where $\bar{R}_{n0}(r)=\bar{\phi}_{n0}\Omega_{D-1}^{1/2}e^{-\rho/2}L_{n0}(\rho)$ with $2\bar{\gamma}_{n\ell}r=\rho$. The $L_{n0}(\rho)$ function is a generalized series

[TABLE]

Only the first three terms in $L_{n0}(\rho)$–the ones shown–are sensitive at short distances, so we isolate the short distance part $\hat{L}_{n0}\equiv 1+\rho/2-\bar{n}_{n0}\rho^{1+2\epsilon}/(2[1+2\epsilon])$ for separate treatment. The remaining terms in the series, starting at $O(\rho^{2})$, we define as ${{}^{2}\!}L_{n0}(\rho)=L_{n0}(\rho)-\hat{L}_{n0}(\rho)$. We write the expectation value with $D\rightarrow 3-2\epsilon$ as

[TABLE]

where $\bar{V}(r)=-\beta/r^{1-2\epsilon}$. Only the first term involving \big{(}\hat{L}_{n0}\big{)}^{2} is divergent (when $\epsilon=0$), but for $\text{Re}(\epsilon)>1/2$ its $r$ integral is convergent and can be expressed in terms of gamma functions. We analytically continue the result to a neighborhood of $\epsilon=0$ and expand in a Laurent series to obtain a series with a $1/\epsilon$ pole expressing the divergence plus a series with non-negative powers of $\epsilon$. The remaining terms containing ${{}^{2}\!}L_{n0}$ are finite and can be evaluated with $\epsilon\rightarrow 0$ using the methods of Sec. V. In all we find

[TABLE]

The expectation value \big{\langle}\bar{V}^{3}\big{\rangle}_{n0} is simpler because only the first term in the series for $L_{n0}$ must be included in $\hat{L}_{n0}$. For other expectation values as many terms should be included in $\hat{L}_{n0}$ as needed.

Extra care is required for some finite expectation values for which the operations of taking the expectation value and the $D\rightarrow 3$ limit do not commute. An example is \big{\langle}r^{-2+4\epsilon}\partial_{r}^{2}\big{\rangle}_{n0}. When taking the second derivative of the radial part of the wave function, second derivatives of $L_{n0}$ occur, which look like

[TABLE]

Negative powers of $\rho$ generated in this way can lead to divergences in the radial integral that would not have been present in three dimensions where the associated Laguerre is a regular polynomial. The radial integral gives a $1/\epsilon$ divergence, which along with the explicit factor of $\epsilon$ produced by the differentiation, leads to an additive constant that would not have appeared in three dimensions. The same separation of $L_{n0}$ into $\hat{L}_{n0}+{{}^{2}\!}L_{n0}$ allows us to evaluate this expectation value safely, finding

[TABLE]

which differs from \big{\langle}r^{-2}\partial_{r}^{2}\big{\rangle}_{n0} evaluated directly in three dimensions (as tabulated in (A.5d)) only in the $\delta_{\ell=0}$ term.

Expectation values involving the momentum squared positioned next to a wave function (or adjoint wave function) may be simplified by use of the Schrödinger equation. We write (II.5) in the form $p^{2}\psi=2m_{r}(H-\bar{V})\psi$, so that

[TABLE]

for any operator $M$.

Recursion relations (Pasternack (1937); Epstein and Epstein (1962)) allow us to express one expectation value in terms of others, which may already be known. Recursion relations can be found by noting that for any operator $A$, the expectation value of $B\equiv[H,A]$ vanishes: \big{\langle}B\big{\rangle}=0. When we deal with radial operators only, we can use the radial Hamiltonian

[TABLE]

Using $A\rightarrow r^{s}$ we find 0=s(s+1-2\epsilon)\big{\langle}r^{s-2}\big{\rangle}+2s\big{\langle}r^{s-1}\partial_{r}\big{\rangle}, or

[TABLE]

Now taking $A\rightarrow r^{s}\partial_{r}$ and using (III.11) to write $\partial_{r}^{2}$ in terms of $H$ and eliminating \big{\langle}r^{s-2}\partial_{r}\big{\rangle} using (III.12), we find

[TABLE]

Additional relations could be found by making other choices for $A$. We can now derive some consequences. Use of (III.12) with $s=-2+2\epsilon$ gives

[TABLE]

(except when $\ell=0$, since \big{\langle}r^{-D}\big{\rangle} is a special case not regulated by dimensional regularization as noted below (III.36); the $\ell=0$ result for \big{\langle}\bar{V}^{\prime}\partial_{r}\big{\rangle} can be obtained using the Fundamental Theorem of Calculus and is given in (A.15f)). Use of (III.12) with $s=-3+4\epsilon$ gives

[TABLE]

Use of (III.13) with $s=0$ and $\langle 1\rangle=1$ gives

[TABLE]

Use of (III.13) with $s=-2+4\epsilon$ gives

[TABLE]

These results are exact in $D$ dimensions.

The Feynman-Hellmann theorem (Weinberg, 2015, p.195) gives additional relationships involving expectation values. For $\lambda$ a parameter of the Hamiltonian, the theorem states that

[TABLE]

Using the radial Hamiltonian of (III.11), the relevant parameters are $m_{r}$, $\beta$, and $\ell$. The equation for $m_{r}$ is

[TABLE]

and the equation for $\beta$ is

[TABLE]

The first order differential equations implied by (III.16), (III.19) and (III.20) lead to

[TABLE]

which can also be seen using (II.15) and (II.23). The equation for $\ell$ leads to a value for \big{\langle}1/r^{2}\big{\rangle} in three dimensions, but its utility in $D$ dimensions is unclear to us.

Momentum space brackets Czarnecki et al. (1999a) are defined according to

[TABLE]

where ${\mathchar 22\relax\mkern-11.0mud}^{D}p=d^{D}p/(2\pi)^{D}$. Momentum space brackets are useful because Feynman rules are often expressed in momentum space and such brackets emerge naturally in the initial stages of a calculation. Fourier transformation can be employed to give an equivalent coordinate space version

[TABLE]

where $M(\vec{x}_{2},\vec{x}_{1})$ is the Fourier transform

[TABLE]

When $M(\vec{p}_{2},\vec{p}_{1})$ can be expressed entirely in terms of the relative momentum $\vec{q}=\vec{p}_{2}-\vec{p}_{1}$, one can write $M(\vec{p}_{2},\vec{p}_{1})=M(\vec{q}\,)$ and $M(\vec{x}_{2},\vec{x}_{1})$ takes the form

[TABLE]

In that case, the momentum space bracket is equal to a corresponding expectation value

[TABLE]

where $M(\vec{x}\,)=\int{\mathchar 22\relax\mkern-11.0mud}^{D}q\,e^{i\vec{q}\cdot\vec{x}}M(\vec{q}\,)$ (see Appendix B) is the Fourier transform of $M(\vec{p}_{2},\vec{p}_{1})=M(\vec{q}\,)$. Brackets involving components of $\vec{p}_{2}$ or $\vec{p}_{1}$ in the numerator can be dealt with by replacing $p_{1i}$, for example, by $i\partial_{1i}$ acting on the exponential $e^{-i\vec{p}_{1}\cdot\vec{x}_{1}}$ followed by a partial integration to let the derivative act on the wave function $\psi(\vec{x}_{1})$. One finds, for instance, that the momentum space bracket of $p_{2i}N(\vec{q}\,)p_{1i}$ can be written as an expectation value according to

[TABLE]

where $p_{i}=-i\partial_{i}$ and $N(\vec{x}\,)$ is the Fourier transform of $N(\vec{q}\,)$.

The momentum space Schrödinger equation

[TABLE]

where $\bar{V}(\vec{p}-\vec{k}\,)=-4\pi Z\alpha\bar{\mu}^{2\epsilon}/(\vec{p}-\vec{k}\,)^{2}$ can be extremely useful for simplifying momentum space brackets directly in momentum space. We find that

[TABLE]

with an extra integral to do that may lead to an overall simplification. An analogous relation holds for \big{\langle}\vec{p}_{2}^{2}M(\vec{p}_{2},\vec{p}_{1})\big{\rangle}_{\vec{p}_{2},\vec{p}_{1}}.

A number of expectation values and momentum space brackets are presented in Appendix A. Many of the three-dimensional ones appear elsewhere, but we thought it would be useful to have a complete list containing both finite and divergent expectation values in one place. The finite ones, where possible, have been checked using numerical integration in three dimensions for many values of the quantum numbers, and we have verified that all of them are consistent with the various recursion relations and other identities.

A few of the expectation values in the lists are singular when evaluated directly in three dimensions but are finite in the $D\rightarrow 3$ limit of the $D$-dimensional result. One of these is \big{\langle}p^{6}\big{\rangle}. After use of the Schrödinger equation on the left and right, we find

[TABLE]

The final expectation value here is divergent in three dimensions as can best be seen in momentum space:

[TABLE]

where the $\vec{k}$ integral is divergent for large $|\vec{k}\,|$ in three dimensions. In $D$ dimensions, for appropriate values of $D$, the $\vec{k}$ integral is well-defined and can be performed by use of a Feynman parameter. We find that

[TABLE]

which has a finite $D\rightarrow 3$ limit. We write this in coordinate space using the Fourier transform formula of Appendix B as

[TABLE]

and note that the three dimensional expectation value \big{\langle}p_{i}V^{2}p_{i}\big{\rangle} has the finite value (A.13i) that can be obtained from (A.3e). The final regularized value of \big{\langle}p^{6}\big{\rangle} is given in (A.3). Another exceptional case is the momentum space bracket \big{\langle}1/\vec{q}\,^{4}\big{\rangle}_{\vec{p}_{2},\vec{p}_{1}}, which is divergent in three dimensions, but when evaluated in $D$ dimensions has a finite $D\rightarrow 3$ limit. We use the $D$-dimensional Fourier transform of Appendix B to see that

[TABLE]

One particular bracket bears special mention–that of $\ln q$ where $q=|\vec{q}\,|$. This bracket appears in calculations, but its evaluation has particular difficulties. The bracket can be expressed as an expectation value using the Fourier transform (FT) of $\ln q$

[TABLE]

We write $\ln q$ as $\ln q=\displaystyle{\lim_{s\rightarrow 0}}\frac{d}{ds}q^{s}$. If we assume that the operation of taking the limit of a derivative commutes with the operation of taking a Fourier transform (done using the formulas of Appendix B), we find that the transform of $\ln q$ is

[TABLE]

This function $1/r^{D}$ is problematic because it gives divergent S state expectation values that are not regulated by dimensional regularization: $\int_{0}^{\infty}dr\,r^{D-1}\bar{R}_{n0}(r)r^{-D}\bar{R}_{n0}(r)$ has a short distance divergence. It seems that our usual approach of evaluating momentum space brackets by Fourier transformation followed by calculating a coordinate space expectation value must fail. We could fall back on momentum space calculations, which work for individual states, but it seems difficult to find a general formula (for all states) by using momentum space methods. Instead, we evaluate the bracket of $\ln q$ by keeping $s$ as the regulating parameter and hold off on taking the $s$ derivative and $s\rightarrow 0$ limit until after the expectation value is done. This procedure gives a finite result for \big{\langle}\!\ln q\big{\rangle}_{\vec{p}_{2},\vec{p}_{1}} that agrees with direct evaluations of momentum space brackets for particular states. Details of the calculation are given in Appendix E.

IV Evaluation of an example of an energy level contribution

We give a sample calculation of a dimensionally regularized contribution to Coulomb energy levels as an illustration of the usefulness of Coulomb expectation values of the sort studied and tabulated in this work. We base our calculation on the dimensionally regularized low-energy effective quantum field theory Non-Relativistic QED (NRQED) Caswell and Lepage (1986). We consider the contribution illustrated in Fig. 1, which shows a process involving exchange of a Coulomb photon having a Coulomb vertex on one end and a higher-order vertex on the other. In particular, we will calculate the effect of the “$CX1$” higher-order vertex, which comes from the term

[TABLE]

in the NRQED Lagrangian (as given in Hill et al. (2013)). The factor $c_{X1}$ is a Wilson coefficient (i.e. a matching coefficient) and lends its name to this contribution. The fermion fields $\psi^{\dagger}$ and $\psi$ are two-component Pauli spinors representing the electron or the “nucleus” (assumed also to be a spin-1/2 fermion) with charge $q$ and mass $m$. The photon field is contained in the electric field $\vec{E}=-\vec{\nabla}A^{0}-\partial_{0}\vec{A}$, and $\vec{D}=\vec{\nabla}-iq\vec{A}$ is the gauge covariant derivative. The Feynman rule from $\delta{\cal L}$ for the interaction of a fermion with a Coulomb photon is $c_{X1}\frac{-iq}{m^{4}}\left(\vec{p}_{2}^{2}-\vec{p}_{1}^{2}\right)^{2}$. A formalism for using NRQED to calculate bound state energy levels is given in Adkins (2018b). It makes use of the NRQED Bethe-Salpeter equation and a perturbative scheme based on a lowest-order problem that can be reduced to the usual Schrödinger equation with a ($D=3-2\epsilon$ dimensional) Coulomb potential. The Bethe-Salpeter wave functions in this formalism are $\bar{\Psi}(p)$ and $\Psi(p)$, and the energy shift coming from the two diagrams of Fig. 1 is

[TABLE]

where $d=4-2\epsilon$ is the dimension of spacetime. The $c_{X1}$ matching coefficient might be different for different particles. For elementary spin 1/2 particles such as the electron and muon (and their antiparticles) it is $c_{X1}=\frac{5}{128}+O(\alpha)$ Hill et al. (2013).

The reduction of (IV.2) to an expectation value and then to a number happens through several steps. We perform the energy integrals by use of a convenient property of the lowest-order NRQED Bethe-Salpeter wave functions: $\int{\mathchar 22\relax\mkern-11.0mud}p_{0}\,\bar{\Psi}(p)=\psi^{\dagger}(\vec{p}\,)$ and $\int{\mathchar 22\relax\mkern-11.0mud}p_{0}\Psi(p)=\psi(\vec{p}\,)$, where $\psi(\vec{p}\,)$ is the spatial wave function discussed in Sec. II. The $CX1$ energy correction simplifies to

[TABLE]

where the integrals now are over the $D=3-2\epsilon$ spatial dimensions, the charge product is $q_{1}q_{2}=-4\pi Z\alpha\bar{\mu}^{2\epsilon}$ as in (II.1), and we define $\vec{q}\equiv\vec{p}_{2}-\vec{p}_{1}$. We expand the numerator factor $(\vec{p}_{2}^{2}-\vec{p}_{1}^{2})^{2}=\vec{p}_{2}^{4}-2\vec{p}_{2}^{2}\vec{p}_{1}^{2}+\vec{p}_{1}^{4}$. We use the Fourier transform as discussed in Sec. III and Appendix B to turn our expression for the energy correction into a standard (D-dimensional) expectation value:

[TABLE]

The required expectation values have been tabulated as (A.13k) and (A.13). Or, we could use the identity of (A.15e) to write \big{\langle}\vec{p}\,^{4}\bar{V}-2\vec{p}\,^{2}\bar{V}\vec{p}\,^{2}+\bar{V}\vec{p}\,^{4}\big{\rangle} as -4m_{r}\big{\langle}(\bar{V}^{\prime})^{2}\big{\rangle}, and so reduce our result for the $CX1$ energy correction to

[TABLE]

where we have used the value for \big{\langle}(\bar{V}^{\prime})^{2}\big{\rangle} tabulated in (A.13e). It is extremely convenient to have results for needed expectation values readily available.

V Integrals and matrix elements involving associated Laguerre polynomials

Many expectation values in $D=3-2\epsilon$ dimensions are finite and reduce to their three-dimensional form in the limit $\epsilon\rightarrow 0$. In this section we develop the techniques for evaluating the types of integrals involving standard three-dimensional associated Laguerre polynomials that occur in the quantum mechanical Coulomb problem. Where possible, general formulas are obtained, and a number of special cases most useful for atomic physics calculations are discussed in detail. Some useful sums that follow from the integration formulas are given as well.

We begin by finding general expressions for the integrals

[TABLE]

(where $s$ must satisfy $s>-1$ in order to avoid a divergence for small $x$) and some related sums. It is also occasionally necessary to find values for integrals like those above but with negative values of $s$ for which the integrals are divergent. In that case subtractions are required–we choose to modify one of the associated Laguerre polynomials in the integrand by removing its first $p$ powers of $x$ where $p+s>-1$. If we write

[TABLE]

for the associated Laguerre function where $c_{nkr}=\frac{(-1)^{r}}{r!}\binom{n+k}{n-r}$ as in (C.20), then we can define

[TABLE]

The subtracted integrals are convergent for $p+s>-1$

[TABLE]

We can immediately obtain a result for $I_{s}(n,k)$ using the series expansion for $L_{n}^{k}(x)$ as

[TABLE]

However, a more convenient form without a sum is possible through use of integration by parts. We use the Rodrigues formula (C.24) for $L_{n}^{k}(x)$ and define $D\equiv d/dx$ to write

[TABLE]

Now we integrate by parts $n$ times using the assumption that the real part of $s$ is sufficiently large that all integrated terms vanish at $x=0$. They vanish for $x\rightarrow\infty$ due to the exponential factor. One finds that

[TABLE]

Results for $I_{s}(n,k)$ with smaller values of $s$ are obtained by analytical continuation and possibly the appropriate limiting case as $s$ approaches the desired value. For example, the integral $I_{0}(n-1,1)$ with $n\geq 1$ has the value

[TABLE]

where (D.7d) was used to see that

[TABLE]

Other integrals, similarly obtained, are

[TABLE]

and

[TABLE]

A general result is

[TABLE]

where we assume that $s$, $k$, and $n$ are all non-negative integers.

We can obtain results for the subtracted integrals ${}^{p}\!I_{s}(n,k)$ starting from (V.5). We find

[TABLE]

The sum can be done, for example by use of Mathematica Wolfram Research, Inc. (2017), giving

[TABLE]

(The hypergeometric functions are defined and some properties are given in Appendix D.) For specific values of $p$ (with $p\geq 0$) this can be simplified. As examples, one has

[TABLE]

where values having $s$ an integer in the range $-p\leq s\leq-1$ must be interpreted as limits. We use (D.39) to rewrite these as

[TABLE]

interpreting the formulas as limits as appropriate. For the once subtracted integral with integral $s$, $s\geq-p=-1$ and $n\geq 1$, one has

[TABLE]

and with $n\geq 2$

[TABLE]

For the twice subtracted integral with integral $s$, $s\geq-p=-2$, one has

[TABLE]

all of which vanish automatically when $n=1$ or $n=2$ since ${}^{2}\!L_{n-1}^{1}(x)$ vanishes in those cases.

The $J_{s}(n,k)$ integrals can be obtained from the $I_{s}(n,k)$ ones using

[TABLE]

and

[TABLE]

One finds immediately that

[TABLE]

where as in $I_{s}(n,k)$ a limiting procedure for $s$ is often required. The result for $J_{s}(n,k)$ can be simplified somewhat using the digamma identity (D.16b):

[TABLE]

As an example we calculate $J_{0}(n-1,1)=\lim_{s\rightarrow 0}J_{s}(n-1,1)$:

[TABLE]

and more generally one has

[TABLE]

Some of the $J_{s}(n-3,5)$ integrals are also useful, and one has

[TABLE]

A more general result is

[TABLE]

where we assume that $s$, $k$, and $n$ are all non-negative integers.

A comparison of (V.5) with (V.7) and a similar comparison for the logarithmic integrals results in a number of useful summation formulas. One finds that

[TABLE]

for non-negative integral values of $n$ and $k$ and any $s$, perhaps with a limit required. For example, with $n\rightarrow n-1$, $k\rightarrow 1$, and $s\rightarrow 0$, one has

[TABLE]

For the logarithmic case, we return to the evaluation of $J_{s}(n,k)$ by the direct method:

[TABLE]

So one has the summation formula

[TABLE]

valid for non-negative integral values of $n$ and $k$, and $s>-1$, perhaps with a limit required. For example, with $n\rightarrow n-1$, $k\rightarrow 1$, $s\rightarrow 0$, one has

[TABLE]

Also, with $n\geq 1$, $k\rightarrow 0$, and $s\rightarrow 0$, one finds

[TABLE]

The denominator gamma factor is singular for all such $n$, so the whole term would vanish except if there is a pole in the curly bracket. All terms there are finite except when $r=1$, in which case they go as $1/s$, and (V.9) leads to

[TABLE]

For the evaluation of $K_{s}(n,k;n^{\prime},k^{\prime})$ we use the expansion (V.2) for one of the associated Laguerre functions and the Rodrigues formula (C.24) for the second and find

[TABLE]

The integration by parts above was done under the assumption that $Re(s)$ is large enough so that the integrated terms vanish for $x\rightarrow 0$ as in the derivation of (V.7). Integrals for smaller values of $s$ are obtained through analytic continuation and, if necessary, a limiting procedure for $s$.

The subtracted $K$ integrals differ from the above only by having a different lower limit for the $r$ summation:

[TABLE]

and converge whenever $s+p>-1$.

We can use (V.35) to prove the Laguerre orthogonality relation (C.26):

[TABLE]

The gamma function in the denominator is infinite for every r in the sum when $m>n$, so $K_{k}(n,k;m,k)$ vanishes in this case. By symmetry under $n\leftrightarrow m$, $K_{k}(n,k;m,k)$ also vanishes for $m<n$. When $m=n$, the only nonvanishing contribution comes from $r=m=n$, with the result

[TABLE]

The orthogonality relation (C.26)

[TABLE]

follows immediately.

Another integral useful for normalizing the Coulomb wave functions is

[TABLE]

Because of the function $\Gamma(r+2-n)$ in the denominator, the only nonvanishing terms are those with $r=n-1$ and $r=n$. One finds that

[TABLE]

Another set of useful $K$ integrals is $K_{1}(n-1,1;n-a,a)$ with $a$ an integer $1\leq a\leq n$, for which we use (V.35):

[TABLE]

For $0\leq r\leq a-1$ there is a divergent gamma function in both the numerator and the denominator (with a finite limit), while for $a\leq r\leq n-1$ there is a divergent gamma function only in the denominator. It follows that

[TABLE]

A similar calculation leads to

[TABLE]

We also find

[TABLE]

and

[TABLE]

As a summary of results for the $K$ integrals one has

[TABLE]

Some values for the subtracted $K$ integrals are

[TABLE]

A general expression for the evaluation of ${}^{p}\!K_{s}(n,k;n^{\prime},k^{\prime})$ can be obtained by considering three possible regions in the sum over $r$. Our general assumptions are that $n,k,n^{\prime},k^{\prime},p$ are all non-negative integers and that the convergence condition $s+p>-1$ holds. While the general formula (V.39) has no restriction that $s$ be an integer, we now assume that $s$ is an integer as well. We note that the factors

[TABLE]

in the sum over $r$ are all well defined and finite for $r$ in the summation range $p\leq r\leq n$. For $p\leq r\leq k^{\prime}-s-1$, both $\Gamma(s+r+1-k^{\prime})$ in the numerator and $\Gamma(s+r+1-k^{\prime}-n^{\prime})$ in the denominator are divergent. The ratio is evaluated as a limit in which $s\rightarrow\epsilon+s$ approaches its integral value $s$ when $\epsilon\rightarrow 0$. One has

[TABLE]

For values of $r$ in the range $k^{\prime}-s\leq r\leq n^{\prime}+k^{\prime}-s-1$, the divergent gamma function in the denominator is not canceled by a numerator divergence and the contribution vanishes. For $r$ in the range $n^{\prime}+k^{\prime}-s\leq r\leq n$, all gamma functions are finite and can be written as factorials. The result for the subtracted $K$ integral is

[TABLE]

This formula works as well when $p=0$, giving the standard $K$ integrals of (V.35). We note that one or both of the summation regions might be void, giving a contribution of zero for that sum. Also, the factorial functions in the denominators enforce the overall limits of summation, keeping $r$ in the range $0\leq r\leq n$.

An associated integral with a log function can be obtained by way of the same method of differentiation by $s$ that gave $J_{s}(n,k)$ from $I_{s}(n,k)$:

[TABLE]

Some examples include:

[TABLE]

The subtracted $L$ integral is defined as

[TABLE]

where $\psi(n,z)$ is defined in (D.20). As an example, one has

[TABLE]

A second differentiation by $s$ gives integrals involving $\ln^{2}x$:

[TABLE]

We find that

[TABLE]

As an example, we evaluate $M_{1}(n-1,1;n-1,1)$. Firstly, by explicit integration, one has the values $\gamma_{E}^{2}-2\gamma_{E}+\zeta(2)$ and $2\gamma_{E}^{2}-6\gamma_{E}+2\zeta(2)+4$ for $n=1$ and $n=2$ respectively. In general one has

[TABLE]

For general $n\geq 2$ one breaks the sum into two parts: $0\leq r\leq n-2$ and $r=n-1$ and works them out separately. The gamma function in the denominator is singular in the sum, so only the singular parts of the numerator contribute and a limiting process is required for their evaluation. The result is

[TABLE]

where the sum is void for $n=1$ and $\text{H}_{0}\equiv 0$. The sums can be done using the formulas of Appendix D, resulting in

[TABLE]

More generally, one has

[TABLE]

where $\text{diH}_{-}(n,m)$ is the falling diharmonic number defined in (D.29).

VI Discussion

We have given an extensive list of expectation values that are useful for the calculation of corrections to the properties of Coulombic two-body bound systems. We used dimensional regularization to regulate divergences and have found the values for a number of expectation values using this regulatory scheme. In addition to exploring the use of dimensional regularization in the Coulomb problem, we gave convenient formulas for a number of integrals involving the standard associated Laguerre polynomials. These integration formulas involved one or two associated Laguerre polynomials, arbitrary powers of the radial variable $r$, and possibly various powers of $\ln r$. We also obtained formulas for analogous integrals involving subtracted associated Laguerre polynomials that were required for the evaluation of divergent expectation values. Most of the tabulated expectation values, Laguerre polynomial integrals, and associated sums described in Sec. V were used in a calculation of corrections to the energy levels of positronium at order $m\alpha^{6}$ using dimensionally regularized nonrelativistic quantum electrodynamics (NRQED) Adkins (2018b); Adkins et al. (2019), and would also appear in calculations of energy levels for any Coulombic two-body bound state at this order when using dimensional regularization.

Appendix A Coulomb expectation values

We define coordinate space expectation values in the usual way:

[TABLE]

Some finite expectation values in $D=3$ dimensions are given here:

[TABLE]

A set of finite expectation values containing the momentum operator $\vec{p}=-i\vec{\nabla}$ follows:

[TABLE]

Expectation values involving the radial derivative $\partial_{r}=\hat{x}_{i}\partial_{i}=i\hat{x}\cdot\vec{p}$ include

[TABLE]

Some expectation values involving $\partial_{r}^{2}$ are

[TABLE]

and also

[TABLE]

Some expectation values involving $\partial_{r}^{3}$ are

[TABLE]

In the following we also use the adjoint radial derivative $\partial_{r}^{\dagger}=-\partial_{i}\hat{x}_{i}=-i\vec{p}\cdot\hat{x}$. The general relation holds:

[TABLE]

from which follows

[TABLE]

Some expectation values involving both $\partial_{r}^{\dagger}$ and $\partial_{r}$ are

[TABLE]

We also find

[TABLE]

Some expectation values involving logs are

[TABLE]

Expectation values involving the Coulomb potential include the following, where we note that S state expectations involving $\bar{V}$ must be evaluated in $D$ dimensions:

[TABLE]

Only the leading terms are shown for S state expectation values involving $\bar{V}$–terms of $O(\epsilon)$ are ignored. All expectation values for $\ell>0$ were evaluated in the $D\rightarrow 3$ limit and the results given are exact in that limit, as are expectation values involving the three-dimensional Coulomb potential $V$.

Second and higher order derivatives in coordinate space expectation values can lead to unexpected results. Some examples of expectation values that do or might seem to lead to unexpected results are shown here:

[TABLE]

All of these formulas are correct only through $O(\epsilon^{0})$–there are uncalculated $O(\epsilon)$ corrections.

Some exact $D$-dimensional relations among expectation values involving the potential $\bar{V}=-\beta/r^{1-2\epsilon}$ where $\beta=\Gamma(1/2-\epsilon)Z\alpha\bar{\mu}^{2\epsilon}\pi^{-1/2+\epsilon}$ include (III.17) and

[TABLE]

Some momentum-space brackets involving $\vec{q}=\vec{p}_{2}-\vec{p}_{1}$ are

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

We can keep higher order terms if necessary, for example

[TABLE]

Appendix B The Fourier transform in $D$ dimensions

The $D$-dimensional Fourier transform is given by Samko (1978)

[TABLE]

where the $Y^{D}_{\ell m}$ are $D$-dimensional spherical harmonics of angular momentum $\ell$, while $m$ represents the indices necessary to distinguish one spherical harmonic with angular momentum $\ell$ from another. Formula (B.1) is singular when the arguments of the gamma functions are non-positive integers. The Fourier transform contains delta functions and their derivatives when $(\ell-n)/2$ is a non-positive integer and logs when $(\ell+D+n)/2$ is a non-positive integer. The three-dimensional special case is standard and is discussed, for example, in Adkins (2016). Some useful transforms include

[TABLE]

where $p=|\vec{p}\,|$. The transforms (B.2a) and (B.2b) are immediate since the transformed functions $1/p^{\alpha}$ and $p_{i}/p^{\alpha}$ have $\ell=0$ and $\ell=1$, respectively, while $p_{i}p_{j}/p^{\alpha}$ is a combination of $\ell=2$ and $\ell=0$, $p_{i}p_{j}p_{k}/p^{\alpha}$ is a combination of $\ell=3$ and $\ell=1$, and $p_{i}p_{j}p_{k}p_{\ell}/p^{\alpha}$ is a combination of $\ell=4$, 2, and 0.

Appendix C Review of the bound-state Coulomb problem in three dimensions

In the appendix we assemble the main formulas relevant to the three-dimensional bound-state Coulomb problem. We use the conventions for associated Laguerre polynomials found in Abramowitz and Stegun (1972); Gradshteyn and Ryzhik (1980); Galindo and Pascual (1990); Arfken and Weber (2001); Olver et al. (2010) and that are built into the computer algebra systems Mathematica Wolfram Research, Inc. (2017) and MatLab MAT (2019). Other references, including many textbooks on quantum mechanics, use a variety of other conventions–a partial tabulation is given in Liboff (Liboff, 1998, Table 10.3).

The coordinate-space Schrödinger-Coulomb equation is

[TABLE]

where $\vec{p}=-i\vec{\nabla}$, $m_{r}=m_{1}m_{2}/(m_{1}+m_{2})$ is the reduced mass and $V(r)=-Z\alpha/r$. For the bound state Coulomb wave functions we have

[TABLE]

where

[TABLE]

with $\rho=2r/(an)$, and $a$, the Bohr radius for a system with positive charge $Z$ and reduced mass $m_{r}$, given by

[TABLE]

in our units. For $Z=1$ (as for hydrogen) with the assumption of no recoil (nucleus of infinite mass, $m_{r}=m_{e}$), the usual Bohr radius is

[TABLE]

We have separated out the short distance behavior (given by $\rho^{\ell}$) and the long distance behavior (controlled by the exponential $e^{-\rho/2}$). As examples, the S state radial functions are

[TABLE]

and the P state radial functions are

[TABLE]

The radial Schrödinger equation is

[TABLE]

and in terms of the dimensionless radial variable $\rho$ it is

[TABLE]

The solutions to the radial equation can be expressed in terms of associated Laguerre polynomials. The Laguerre polynomials are the conventionally normalized solutions to the equation

[TABLE]

and have the generating function

[TABLE]

Rodrigues’ formula

[TABLE]

recursion relations

[TABLE]

and orthogonality relation

[TABLE]

An explicit formula is

[TABLE]

An integral formula for them can be derived easily by use of a Laplace transform to solve the defining differential equation:

[TABLE]

where the integration contour circles the origin in a counterclockwise direction Merzbacher (1961), and they can be expressed in terms of hypergeometric functions as

[TABLE]

The associated Laguerre polynomials

[TABLE]

are the regular solutions to the equation

[TABLE]

An explicit formula for them is

[TABLE]

We use the general definition of the binomial symbol

[TABLE]

which works even when $a$ and $b$ aren’t integers, to extend the definition of $L_{n}^{k}(x)$ to non-integral $k$. We define

[TABLE]

The associated Laguerre polynomials have the generating function

[TABLE]

Rodrigues’ formula

[TABLE]

recursion relations

[TABLE]

and orthogonality relation

[TABLE]

An integral formula for the associated Laguerre polynomials is

[TABLE]

and their explicit expression in terms of hypergeometric functions is

[TABLE]

The Schrödinger-Coulomb equation in momentum space has the form

[TABLE]

where ${\mathchar 22\relax\mkern-11.0mud}^{n}p=\frac{d^{n}p}{(2\pi)^{n}}$. The momentum space Coulomb wave functions are Fourier transforms of the coordinate space functions:

[TABLE]

where

[TABLE]

The constants in (C.31) are the coordinate space wave function at contact

[TABLE]

and

[TABLE]

and we have used $D_{n}=p^{2}+\gamma_{n}^{2}$ and $\overline{D}_{n}=p^{2}-\gamma_{n}^{2}$. The momentum space version of the Coulomb wave function can be obtained by Fourier transformation of the coordinate space wave function Podolsky and Pauling (1929), by direct solution of the momentum-space Schrödinger-Coulomb equation Fock (1936), or from the the momentum-space Green function Schwinger (1964); Lieber (1989); Adkins and Cammerata (1992). Some values for the Gegenbauer polynomials are $C_{0}^{\lambda}(\beta)=1$, $C_{1}^{\lambda}(\beta)=2\lambda\beta$, $C_{2}^{\lambda}(\beta)=2\lambda(\lambda+1)\beta^{2}-\lambda$, etc. They satisfy the Gegenbauer differential equation

[TABLE]

the Rodrigues formula

[TABLE]

and have the series expansion

[TABLE]

where the floor $\lfloor n/2\rfloor$ is the greatest integer in $n/2$. The Gegenbauer generating function is

[TABLE]

The natural interval is $-1\leq\beta\leq 1$, corresponding to $0\leq p<\infty$. The Gegenbauer polynomials satisfy the orthonormality condition

[TABLE]

Various derivative identities exist, including

[TABLE]

Appendix D Special functions: gamma, polygamma, beta, zeta, hypergeometric, and harmonic and diharmonic numbers

In this Appendix we give definitions and useful formulas for a number of special functions. In particular, we consider the gamma, polygamma, beta, zeta, and hypergeometric functions. We also discuss properties of the harmonic numbers that arise, for example, when expanding the gamma function near a non-positive integer. We define “diharmonic numbers” as the natural generalization of harmonic numbers for use in related circumstances. Our definitions follow the standards set, for example, in Abramowitz and Stegun (1972); Gradshteyn and Ryzhik (1980); Olver et al. (2010).

The gamma function is defined by

[TABLE]

so that

[TABLE]

The gamma function satisfies $\Gamma(1)=1$ and the recursive formula

[TABLE]

so that $\Gamma(n)=(n-1)!$ if $n\geq 1$ is an integer. It also satisfies the reflection identity

[TABLE]

and the duplication formula

[TABLE]

Some useful values for integral and half-integral arguments are

[TABLE]

Some identities that prove useful in the evaluation of $\langle r^{p}\rangle$ for various values of $p$ and elsewhere are

[TABLE]

where $n$ is a positive integer and $\text{H}_{n}$ is the harmonic number

[TABLE]

We also define the generalized harmonic number $\text{H}_{n}^{(\alpha)}$ as

[TABLE]

which has the form of a truncated version of the sum used to define the Riemann zeta function

[TABLE]

Expansion of the gamma function near integral arguments can be facilitated by the formula

[TABLE]

The beta function $B(a,b)$ is defined as

[TABLE]

Many useful integrals can be done in terms of the beta function, including

[TABLE]

Derivatives of the gamma function are given in terms of the polygamma functions. Specifically, the digamma function $\psi(z)$ is the logarithmic derivative of $\Gamma(z)$:

[TABLE]

As a direct consequence of $\Gamma(z+1)=z\Gamma(z)$ one has

[TABLE]

and the derived formulas

[TABLE]

Some useful values for integral and half-integral arguments are

[TABLE]

where the Euler, or Euler-Mascheroni, constant is $\gamma_{E}\approx 0.57721\,56649$. Near the origin $\psi(z)$ has a pole and the expansion

[TABLE]

and also, for $N$ a non-negative integer,

[TABLE]

The higher polygamma functions are

[TABLE]

Some useful values and relations involving the polygamma function are

[TABLE]

The functions $\psi(n,z)$ have a pole of order $n+1$ near the origin and the expansion

[TABLE]

and satisfy identities of the sort

[TABLE]

The Hurwitz generalized zeta function $\zeta(n,x)$ is defined as

[TABLE]

and satisfies

[TABLE]

with $\text{H}_{0}^{(s)}\equiv 0$. The polygamma function satisfies a corresponding identity

[TABLE]

In particular, one has

[TABLE]

Some useful harmonic number sums are

[TABLE]

where the latter sum comes from equation (24) of Spieß Spieß (1990).

The “diharmonic” sums shown above are of the form $\sum_{\text{region}}1/(ij)$ where “region” is some region in the plane of integer pairs $(i,j)$ with $1\leq i,j<\infty$. Diharmonic sums for rectangular regions are easy to work out since the sums factorize. For the rectangle $a\leq i\leq b$, $c\leq j\leq d$, the diharmonic sum is $\sum_{a\leq i\leq b,\;c\leq j\leq d}1/(ij)=(\text{H}_{b}-\text{H}_{a-1})(\text{H}_{d}-\text{H}_{c-1})$. More generally, in analogy to the dilogarithm function, we define the rising and falling “diharmonic numbers”

[TABLE]

The regions in the $(i,j)$ plane over which the sums of $1/(ij)$ are taken are illustrated in Fig. 2. For the purposes of these definitions, we assume that $\text{H}_{n}=0$ for $n\leq 0$ and that the summations vanish when the upper limit is less than the lower. The diharmonic sums allow us to find the values of the sum of $1/(ij)$ over any polygonal region of $(i,j)$–space with $1\leq i,j<\infty$ having sides that are either vertical, horizontal, or diagonal at $45$ degrees. Some special cases come from (D.28) above

[TABLE]

and

[TABLE]

The diharmonic numbers are much more general, as the special cases both come from regions with a vertex at $i=j=1$ and with one side along the horizontal axis. Some more general sums appear in Spieß, but apparently only for triangular regions with the diagonal on or close to the line $i=j$. The diharmonic sums for arbitrary polygonal regions in the $(i,j)$ plane having all sides either horizontal, vertical, or diagonal at $45$ degrees can be obtained as linear combinations of $\text{diH}_{\pm}(n,m)$ and the sums for rectangular regions. Closed forms for diharmonic numbers of nearby regions can be obtained using the recursion relations

[TABLE]

Values for diharmonic numbers of nearly-triangular regions include:

[TABLE]

The diharmonic numbers also satisfy reflection identities

[TABLE]

found by consideration of the regions in $(i,j)$ space covered by the corresponding sums.

The hypergeometric function ${}_{2}F_{1}(a,b;c;x)$ is defined by

[TABLE]

and more generally by

[TABLE]

where the rising and falling factorials are

[TABLE]

(The rising factorial is also known as the “Pochhammer symbol”: $a^{\bar{n}}=(a)_{n}$. Other notations for the rising and falling factorials are $a^{\bar{n}}=a^{(n)}$ and, confusingly, $a^{\underline{n}}=(a)_{n}$. One must always verify definitions when dealing with rising and falling factorials. The advantage of the notation chosen here is that it is unambiguous and suggestive of its meaning.) Rising and falling factorials satisfy identities reminiscent of the binomial formula:

[TABLE]

The hypergeometric series terminates when $a$ or $b$ is a negative integer. A hypergeometric function identity that is useful when doing sums involving factorials and binomial coefficients is given as (5.93) of Graham, Knuth, and Patashnik Graham et al. (1994):

[TABLE]

One often requires series expansions of ${}_{2}F_{1}(a,b;c;x)$ when $a$, $b$, and $c$ differ from integers by terms of order $\epsilon$ with one of those integers in $a$ or $b$ being negative. For example, when $a\rightarrow-n+a\epsilon$, $b\rightarrow k+b\epsilon$, $c\rightarrow-n+1+a\epsilon$, and $x\rightarrow-1$, with $n$ an integer $n\geq 1$ and $k$ an integer $k\geq 0$, we define

[TABLE]

For the case $n\geq 1$, $k\geq 1$, one has

[TABLE]

where sums are taken to be zero when the upper limit is smaller than the lower limit. The first sum here is easy to evaluate for given $n$ and $k$. The second term can be expanded in $\epsilon$ using

[TABLE]

The final infinite sum evaluates to

[TABLE]

which, at least for integers $n$, $k$ with $0\leq n,k\leq 16$, evaluates to a form $\alpha\ln 2+\beta$ for some $\alpha$, $\beta$. So in all one has

[TABLE]

Other cases require separate analyses, which proceed along similar lines. For instance, with $n\geq 1$, $k=0$ one finds

[TABLE]

Appendix E Evaluation of the bracket \big{\langle}\!\ln q\big{\rangle}_{\vec{p}_{2},\vec{p}_{1}}

The momentum space bracket \big{\langle}\!\ln q\big{\rangle}_{\vec{p}_{2},\vec{p}_{1}} is easy to evaluate for any specific state via momentum space integration. This is a finite bracket with no need for regularization, so we calculate in three dimensions. In order to obtain a general formula for this bracket, it is best to work with coordinate space expectation values. However, the three dimensional Fourier transform of $\ln q$ is $\text{FT}[\ln q](\vec{x}\,)=-1/(4\pi r^{3})$, which has a divergent expectation value for S states. Dimensional regularization fails to regulate the Fourier transform of $\ln q$, which is proportional to $1/r^{D}$ in $D$ dimensions, so we stick to three dimensions and regulate differently. We note that $\frac{d}{ds}q^{s}=q^{s}\ln q$ and

[TABLE]

so we use $s$ as a regulating parameter and hold off on the derivative and $s\rightarrow 0$ limit until after taking the expectation value. That is, we evaluate \big{\langle}\!\ln q\big{\rangle}_{\vec{p}_{2},\vec{p}_{1}} as

[TABLE]

The three-dimensional Fourier transform of $q^{s}$, from Appendix B, is

[TABLE]

with the expectation value

[TABLE]

after use of the S-state wave function of (C.6) and the integration formula (V.35). Application of an $s$ derivative and the limit $s\rightarrow 0$ yields the S-state bracket of $\ln q$:

[TABLE]

We can simplify the expression in curly brackets slightly by use of $\psi(-s-1+r)-\psi(-s+r-n)=\sum_{j=2}^{n}\frac{1}{r-s-j}$, which has a linear divergence for $s\rightarrow 0$ when $j=r$. In the $s\rightarrow 0$ limit $\Gamma(-s/2)$ has a linear divergence, as does $\Gamma(-s-1+r)$ for $r=0$ and $r=1$, $\Gamma(-s+r)$ for $r=0$, and $\Gamma(-s+r-n)$ for $0\leq r\leq n-1$. Near $s=0$, the expression in curly brackets is $\ln\left(\frac{2}{an}\right)+H_{n}+O(s)$ for $r=0$, $1/s+O(s^{0})$ for r=1, and $O(1/s)$ for $2\leq r\leq n-1$. We see that the only terms in the sum over $r$ that contribute to the $s\rightarrow 0$ limit have $r=0$ and $r=1$, and we find

[TABLE]

for S states. We have checked this formula for \big{\langle}\!\ln q\big{\rangle}_{\vec{p}_{2},\vec{p}_{1}} by explicit integration in momentum space for $1\leq n\leq 36$. For $\ell>0$ the bracket \big{\langle}\!\ln q\big{\rangle}_{\vec{p}_{2},\vec{p}_{1}}=-\frac{1}{4\pi}\big{\langle}r^{-3}\big{\rangle} is finite and can be evaluated in the usual way. Our result for \big{\langle}\!\ln q\big{\rangle}_{\vec{p}_{2},\vec{p}_{1}} agrees with that of Titard and Ynduráin Titard and Ynduráin (1994).

Acknowledgements.

We acknowledge the support of the National Science Foundation through Grant No. PHY-1707489 and of the Franklin & Marshall College Hackman Scholars Program.