Integrals whose saddle-point expansions terminate

Formulations in mathematical physics have occasionally been reverse-engineered to provide new solutions, or new insights into old ones. Examples are: determining potentials, or refractive-index profiles, that are reflectionless [1]; finding time-dependent Hamiltonians that do not induce transitions between eigenstates of a related system [2–4]; and getting new solutions of the Schrödinger equation by fixing the Madelung–Bohm potential [5].

Here, in a similar spirit, I reverse-engineer the saddle-point approximation [6, 7], leading to a class of canonical integrals for which the technique is exact. In its simplest version, saddle-point integration concerns integrals of the form

$$\begin{equation}I\left(k\right)=\sqrt{\frac{k}{\pi }}\underset{-\infty }{\overset{\infty }{\int }}\mathrm{d}x\enspace \mathrm{exp}\left(-kf\left(x\right)\right),\end{equation} \tag{ 1.1 }$$

in which k is a large parameter and f(x) a smooth function, increasing away quadratically from a minimum conveniently located at x = 0. The saddle-point (=steepest-descent = stationary-phase with 45° rotated contour) approximation is central to theoretical physics: in statistical mechanics, where k is inverse temperature, and in wave mechanics, where k is imaginary and proportional to wavenumber or inverse Planck's constant. The technique generates an expansion of I(k) as a series of inverse powers, whose coefficients involve derivatives of f(x) at the origin. To save writing, we define the minimum value to be zero, and derivatives

$$\begin{equation}{f}_{n}\equiv {\partial }_{x}^{n}f{\left(x\right)}_{x=0},\quad {f}_{0}={f}_{1}=0,\quad {f}_{2}{ >}0.\end{equation} \tag{ 1.2 }$$

In the resulting formal series for I(k),

$$\begin{equation}I\left(k\right)=\sqrt{\frac{2}{{f}_{2}}}\sum _{n=0}^{\infty }\frac{{I}_{n}}{{k}^{n}},\end{equation} \tag{ 1.3 }$$

the coefficients I_n can be found by standard expansion methods; the first few are listed in table 1 in the appendix.

Table 1. Coefficients in the saddle-point expansion, from p 119 of [8].

n	In
0	1
1	$\begin{gathered}{c}\frac{1}{24{f}_{2}^{3}}\left(5{f}_{3}^{2}-3{f}_{2}{f}_{4}\right)\end{gathered}$
2	$\begin{gathered}{c}\frac{1}{1152{f}_{2}^{6}}\left(385{f}_{3}^{4}-630{f}_{2}{f}_{3}^{2}{f}_{4}+105{f}_{2}^{2}{f}_{4}^{2}+168{f}_{2}^{2}{f}_{3}{f}_{5}-24{f}_{2}^{3}{f}_{6}\right)\end{gathered}$
3	$\begin{gathered}{c}\frac{1}{414\enspace 720{f}_{2}^{9}}\left(425\enspace 425{f}_{3}^{6}-1126\enspace 125{f}_{2}{f}_{3}^{4}{f}_{4}+675\enspace 675{f}_{2}^{2}{f}_{3}^{2}{f}_{4}^{2}\right.\\ -51975{f}_{2}^{3}{f}_{4}^{3}+360\enspace 360{f}_{2}^{2}{f}_{3}^{2}{f}_{5}-249\enspace 480{f}_{2}^{3}{f}_{3}{f}_{4}{f}_{5}+13608{f}_{2}^{4}{f}_{5}^{2}\\ \left.-83160{f}_{2}^{3}{f}_{3}^{2}{f}_{6}+22680{f}_{2}^{4}{f}_{4}{f}_{6}+12960{f}_{2}^{4}{f}_{3}{f}_{7}-1080{f}_{2}^{5}{f}_{8}\right)\end{gathered}$

Usually, the series (1.3) is infinite, and factorially divergent in ways now well understood [8, 9]; the divergence depends on saddles and other singularities of f(x) away from the real axis. But for special forms of f(x) the series can terminate after N + 1 terms: I_n = 0 if n > N; then saddle-point integration is exact for all k. These special forms will be determined in section 3. This phenomenon, in which a series terminates for special values of a parameter, is familiar in the different context of the power series in the wavefunction for the quantum harmonic oscillator and the hydrogen atom): the series is a finite (Hermite, Laguerre) polynomial when the energy is one of the eigenvalues.

The saddle-point expansion for exponent functions f(x) that increase quadratically away from their critical point is considered in section 2 because it is the simplest and most familiar situation. But the same procedure for calculating the canonical integrals applies, mutatis mutandis, when the increase of f(x) from the saddle is any power of x. The generalisation is outlined in section 3, which also includes some concluding remarks.

To find the exponent functions f(x) for which the series terminates, we start from the standard procedure in which the integration variable x is changed to X:

$$\begin{equation}f\left(x\right)={X}^{2}{\Rightarrow}I\left(k\right)=\sqrt{\frac{k}{\pi }}\underset{-\infty }{\overset{\infty }{\int }}\mathrm{d}X\left(\frac{\mathrm{d}x}{\mathrm{d}X}\right)\mathrm{exp}\left(-k{X}^{2}\right).\end{equation} \tag{ 2.1 }$$

The series is now determined by expanding the Jacobian dx/dX in powers of X. We require that this series terminates:

$$\begin{equation}\left(\frac{\mathrm{d}x}{\mathrm{d}X}\right)\equiv {G}^{\left(N\right)}\left(X\right)=\sum _{n=0}^{2N}{G}_{n}{X}^{n}.\end{equation} \tag{ 2.2 }$$

This is a polynomial of degree 2N, which must have no real zeros in order for the mapping X(x) to be nonsingular. This restricts the coefficients; in particular G₀ > 0 and G_2N > 0. The result is a finite series for I, depending on k and the chosen coefficients G_n :

$$\begin{equation}I\left(k;\left\{{G}_{2n}\right\}\right)=\sum _{n=0}^{N}\frac{{G}_{2n}}{{k}^{n}}\frac{{\Gamma}\left(n+\frac{1}{2}\right)}{\sqrt{\pi }}={G}_{0}+\frac{{G}_{2}}{2k}+\frac{3{G}_{4}}{4{k}^{2}}+\cdots \frac{{\Gamma}\left(2N\right){G}_{2N}}{{\Gamma}\left(N\right){2}^{2N-1}{k}^{N}}.\end{equation} \tag{ 2.3 }$$

Note that although the odd coefficients will contribute to f(x), only the even coefficients G_2n contribute to the value of the integral.

The mapping X(x) is found by integrating the Jacobian and then inverting the 2N + 1 degree polynomial:

$$\begin{equation}x=\underset{0}{\overset{X}{\int }}\mathrm{d}{X}^{\prime }\enspace {G}^{\left(N\right)}\left({X}^{\prime }\right)=\sum _{n=0}^{2N}\frac{{G}_{n}}{n+1}{X}^{n+1}{\Rightarrow}X={X}^{\left(N\right)}\left(x\right).\end{equation} \tag{ 2.4 }$$

Of the 2N + 1 solutions, it is necessary to choose the one for which X^(N)(x) ∝ x at the origin, in order for the mapping to achieve its aim of standardising the behaviour at the saddle-point. For small and large x, the mapping is

$$\begin{align}\hfill {X}^{\left(N\right)}\left(x\right)& =\frac{x}{{G}_{0}}\left(1+O\left(x\right)\right)\hfill \\ \hfill {X}^{\left(N\right)}\left(x\right)& =\mathrm{sgn}\left(x\right){\left(\frac{\left(2N+1\right)\left\vert x\right\vert }{{G}_{2N}}\right)}^{\frac{1}{2N+1}}\left(1+O\left({x}^{-\frac{1}{2N+1}}\right)\right).\hfill \end{align} \tag{ 2.5 }$$

From (2.1), the exponent function determining the integral whose series truncates after N + 1 terms is

$$\begin{equation}{f}^{\left(N\right)}\left(x\right)={\left({X}^{\left(N\right)}\left(x\right)\right)}^{2}.\end{equation} \tag{ 2.6 }$$

For N = 0, the procedure is trivial:

$$\begin{equation}{X}^{\left(0\right)}\left(x,{G}_{0}\right)=\sqrt{{f}^{\left(0\right)}\left(x,{G}_{0}\right)}=\frac{x}{{G}_{0}},\end{equation} \tag{ 2.7 }$$

reproducing the fact that saddle-point integration is exact for Gaussian integrals. N = 1 is non-trivial; solving the cubic equation (2.4) gives

$$\begin{equation}\begin{aligned}\hfill {X}^{\left(1\right)}\left(x,{G}_{0},{G}_{1},{G}_{2}\right)& =\sqrt{{f}^{\left(1\right)}\left(x,{G}_{0},{G}_{1},{G}_{2}\right)}\hfill \\ \hfill & =-\frac{{G}_{1}}{2{G}_{2}}+{\left(\frac{3}{2{G}_{2}}\right)}^{1/3}\hfill \\ \hfill & \quad {\times}\left[{\left(\sqrt{{\xi }^{2}+{y}^{3}}+\xi \right)}^{1/3}-{\left(\sqrt{{\xi }^{2}+{y}^{3}}-\xi \right)}^{1/3}\right],\hfill \\ \hfill \mathrm{w}\mathrm{h}\mathrm{e}\mathrm{r}\mathrm{e}\;y& \equiv \frac{4{G}_{0}{G}_{2}-{G}_{1}^{2}}{{\left(12{G}_{2}^{2}\right)}^{2/3}},\quad \xi =x+\frac{{G}_{1}\left(6{G}_{0}{G}_{2}-{G}_{1}^{2}\right)}{12{G}_{2}^{2}}.\hfill \end{aligned}\end{equation} \tag{ 2.8 }$$

The condition for the mapping X(x) to be nonsingular is |G₁| < 2√(G₀ G₂), i.e. y > 0. Although the integrand exp(−kf(x)) depends on G₁, the integral does not:

$$\begin{equation}I\left(k;{G}_{0},{G}_{2}\right)={G}_{0}+\frac{{G}_{2}}{2k}.\end{equation} \tag{ 2.9 }$$

For N > 1, determining the mapping, and hence the exponent functions f(x), involves finding roots of quintic and higher polynomials. Although this cannot be carried out explicitly, the procedure is well-defined, and specifies the canonical integrals uniquely. The number of parameters can be reduced by 2, by the scalings

$$\begin{equation}\alpha ={\left(\frac{{G}_{0}}{{G}_{2N}}\right)}^{\frac{1}{2N}},\quad X=\alpha {X}^{\prime },\end{equation} \tag{ 2.10 }$$

and the transformations

$$\begin{equation}{k}^{\prime }=\alpha k,\quad \left\{{{G}^{\prime }}_{0}=1,\quad {{G}^{\prime }}_{0{< }2n{< }2N}={\alpha }^{2n}\frac{{G}_{2n}}{{G}_{0}},\quad {{G}^{\prime }}_{2N}=1\right\}.\end{equation} \tag{ 2.11 }$$

resulting in

$$\begin{equation}I\left(k;\left\{{G}_{2n}\right\}\right)={G}_{0}I\left({k}^{\prime };{{G}^{\prime }}_{2n}\right).\end{equation} \tag{ 2.12 }$$

Therefore the integral depends on 2N–1 essential parameters, in addition to k. For N = 1 this leaves G₁ as the only parameter, and mappings illustrated in figure 1.

Zoom In Zoom Out Reset image size

From the forms of f(x), it is not obvious that the original saddle-point expansion is truncated, i.e. I_n>N = 0. But it is, as will now be demonstrated for N = 1. Expansion of the exponent function f(x) given by (2.8) for G₀ = G₂ = 1 gives the derivatives

$$\begin{equation}\begin{gathered}{c}{f}^{\left(1\right)}\left(x,1,{G}_{1},1\right)=\sum _{2}^{\infty }\frac{{f}_{n}{x}^{n}}{n!},\quad {f}_{2}=2,\enspace {f}_{3}=-6{G}_{1},\enspace {f}_{4}=-16+30{G}_{1}^{2},\\ {f}_{5}=240{G}_{1}-210{G}_{1}^{3},\enspace {f}_{6}=560-3360{G}_{1}^{2}+1890{G}_{1}^{4},\dots .\end{gathered}\end{equation} \tag{ 2.13 }$$

Substitution into table 1 in the appendix gives I₁ = 1/2 (cf (2.9), independent of G₁ as claimed, and I₂ = I₃ = 0.

In the method described here, an asymptotic technique that usually generates a divergent infinite series has been reverse-engineered to find integrals whose expansions terminate. The exponent functions f(x) for which this is achieved are very special. An explicit test that would determine whether a given f(x) is a member of this class would be useful but I have not been able to find one.

The procedure is easily generalised, in two ways. First to consider exponent functions f(x) initially increasing as x^M , rather than x² as in the saddle-point method; second, to consider integrals over the positive real axis rather than the whole real line. This class of integrals, generalising (1.1), is

$$\begin{equation}I\left(k\right)=\frac{{k}^{\frac{1}{M}}}{{\Gamma}\left(\frac{1}{M}\right)}\underset{0}{\overset{\infty }{\int }}\mathrm{d}x\enspace \mathrm{exp}\left(-kf\left(x\right)\right),\quad {f}_{n{< }M}\left(0\right)=0,\quad {f}_{M}{ >}0.\end{equation} \tag{ 3.1 }$$

The natural mapping is now f(x) = X^M . The same Jacobian occurs, and when represented by the same polynomial (2.2), with 2N replaced by N, it leads to the same mapping X(x) (equation (2.4)) as before, also with 2N replaced by N, resulting in the exact truncated series, replacing (2.4),

$$\begin{equation}I\left(k\right)=\sum _{0}^{N}\frac{{\Gamma}\left(\frac{n+1}{M}\right)}{{\Gamma}\left(1M\right)}\frac{{G}_{n}}{{k}^{nM}}.\end{equation} \tag{ 3.2 }$$

The exponent function for which this truncated series exactly represents (3.1) is (2.6) with the exponent 2 replaced by M.

A technically more difficult generalisation would be to find truncated versions of uniform asymptotic expansions [10, 11], where saddle-points within the integration range coalesce as a parameter (additional to k) varies.

The parameter k in (1.1) need not be real. If it is imaginary, the analysis still applies, and gives truncated series for the method of stationary phase that is fundamental in diffraction theory [12] as the way to connect wave physics with ray optics or classical mechanics.

The exponent functions f(x) in the canonical integrals possess branch points away from the real axis, at the values x mapped by the zeros of G(X). In familiar asymptotics, such singularities contribute to the divergence of the expansion; for example, $f \left(x\right)=\sqrt{{x}^{2}+{y}^{2}}-y$ in (1.1) generates the Bessel function $I \left(k\right)=2\sqrt{k/\pi }y\enspace \mathrm{exp}\left(ky\right){K}_{1}\left(ky\right)$ [13], whose large y expansion, generated by $G\left(X\right)=2\left({X}^{2}+y\right)/\sqrt{{X}^{2}+2y}$, diverges. By contrast, the contributions from the branch-points must cancel for the canonical integrals considered here, because their series terminate.

In order to obtain terminated expansions, we have considered cases where the Jacobian G(X) is a polynomial in X. But any function G(X) for which the integral ${\int }_{-\infty }^{\infty }\mathrm{d}X\enspace G\left(X\right)\enspace \mathrm{exp}\left(-{X}^{2}\right)$ can be evaluated analytically leads to integrals that can be evaluated exactly. A class of such integrals corresponds to G(X) = 1 + H(X), where H(X) is any odd function. These integrals are not of the form (1.1) because their exponent functions f(x) depend on k.

I thank Professor C J Howls and Professor Pragya Shukla for helpful comments. This research started in 1980 while I enjoyed the hospitality of the Gleb Wataghin Institute of Physics at UNICAMP, Brazil.

With computer algebra, it is not difficult to generate the coefficients I_n in (1.3), in terms of the derivatives f_n . Table 1 shows the first few. The book [8] conveniently lists more: not only the coefficients I_n for n ⩽ 4, but also their dependence on the derivatives of a k-independent prefactor in the integrand, and the coefficients of half-integer powers of k that occur if the integration in (1.1) is restricted to the positive real axis (as in (3.1) with N = 2).