ARGONNE NATIONAL LABORATORY

9700 South Cass Avenue Argonne, IL 60439

———————————-

ANL/MCS-TM-180 ———————————-

ASYMPTOTIC ANALYSIS

Working Note #2

APPROXIMATION OF INTEGRALS

by

Marc Garbey^{∗} and Hans G. Kaper

Mathematics and Computer Science Division

Technical Memorandum No. 180

July 1993

This work was supported by the Office of Scientific Computing, U.S. Department of Energy, under Contract W-31-109-Eng-38.

∗Address: Université Claude Bernard Lyon 1, Laboratoire d’Analyse Numérique, 69622 Villeurbanne cedex, France

Preface

For some time, we have been interested in the development and application of asymptotic methods for the numerical solution of boundary value problems with critical parameters— that is, parameters that determine the nature of the solution in some critical way. We are thinking, for example, of fluid flow (viscosity), combustion (Lewis number), and superconductivity (Ginzburg-Landau parameter) problems. Their solution may remain smooth over a wide range of parameter values, but as the parameters approach critical values, complicated patterns may emerge. Boundary layers may develop, or the region over which the solution extends may take on the appearance of a patchwork of subregions; on each subregion, the solution is smooth, but between subregions the solution undergoes dramatic changes over very short distances. Shock layers in fluid flow are a visible manifestation of this type of behavior.

Boundary value problems with critical parameters pose some of the most challenging problems in computational science, and much effort is being spent on developing new techniques for their numerical solution. Some of the most useful techniques, in particular on parallel computing architectures, are based on domain decomposition. In a domain decomposition method, one partitions the domain into subdomains, approximates the solution on each subdomain, and assembles these solutions to obtain an approximate solution on the entire domain. Many criteria, involving considerations from linear algebra to computer architecture, go into the design of a useful domain decomposition method. Our aim is to explore the use of asymptotic methods.

Asymptotic analysis, in particular singular perturbation theory, is the study of boundary value problems involving critical parameters. It provides a methodology to identify and characterize boundary layers, transition layers, and initial layers; hence, our idea to use asymptotic methods in the design of domain decomposition algorithms.

We have organized two workshops on the subject of asymptotic analysis and domain decomposition: a workshop at Argonne, jointly sponsored by the Department of Energy and the National Science Foundation (February 1990), and a NATO Advanced Research Workshop in Beaune, France (May 1992). Proceedings of these workshops have been published (Asymptotic analysis and the numerical solution of partial differential equations, edited by H. G. Kaper and M. Garbey, Lecture Notes in Pure and Applied Mathematics – Vol. 130, Marcel Dekker, Inc., New York, 1991; Asymptotic and numerical methods for partial differential equations, edited by H. G. Kaper and M. Garbey, NATO ASI Series C: Mathematical and Physical Sciences – Vol. 384, Kluwer Academic Publishers, Dordrecht, Neth., 1993).

We currently have plans to develop a full-length book on the subject. To formulate our thoughts before final publication, we intend to produce a series of Working Notes on various relevant topics. Some of the notes will contain new material; others may offer new presentations of existing material. We certainly expect the notes to evolve in time; the notes may or may not appear eventually as chapters of the book. The notes are intended for our own use, but we will be happy to supply copies to interested colleagues.

Marc Garbey, Lyon, France Hans G. Kaper, Argonne, Illinois, USA

Working Note #1: Asymptotic Analysis—Basic Concepts and Definitions, ANL/MCS-TM-179 (July 1993)

Contents

Abstract 1 1 Integration by Parts 1 2 Laplace Integrals 3 3 Laplace’s Method 5 4 Fourier Integrals 8 5 Stokes’ Method of Stationary Phase 11 References 16 Index 17

ASYMPTOTIC ANALYSIS

Working Note #2

APPROXIMATION OF INTEGRALS

by

Marc Garbey and Hans G. Kaper

Abstract

In this note we discuss the approximation of integrals that depend on a parameter. The basic tool is simple, namely, integration by parts (Section 1). Of course, the power of the tool is evidenced in applications. The applications are many; they include Laplace integrals (Section 2), generalized Laplace integrals (Section 3), Fourier integrals (Section 4), and Stokes’ method of stationary phase for generalized Fourier integrals (Section 5). These results illustrate beautifully Hardy’s concept of applications of mathematics, that is, “certain regions of mathematical theory in which the notation and the ideas of the [method of integration by parts] may be used systematically with a great gain in clearness and simplicity” [G. H. Hardy, Orders of infinity, App. II, Cambridge University Press, Cambridge, England (1910)].

The notation differs slightly from Working Note #1, for reasons that are mainly historical. The asymptotic analysis of integrals originated in complex analysis, where the (real or complex) parameter, usually denoted by x, is usually introduced in such a way that the interesting behavior of the integrals occurs when x → ∞ in some sector of the complex plane. As there is nothing sacred about notation, and historical precedent is as good a guide as any, we follow convention and denote the parameter by x, focusing on the behavior of integrals as x → ∞ along the real axis or, if x is complex, in some sector of the complex plane. The connection with the notation of Working Note #1 is readily established by identifying the small parameter ϵ with |x|^{−1}.

1 Integration by Parts

One of the simplest ways of finding the asymptotic expansion of a function defined by a definite integral is the method of integration by parts. The successive terms of the expansion are produced by repeated integration by parts, and the asymptotic character of the series is then proved by examining the remainder, which is in the form of a definite integral. Although the field of applications of the method of integration by parts is somewhat limited, the method is very powerful when it works. Precise theorems of sufficient generality are difficult to formulate; examples can be found in [1, Section 2.1] and [2, Chapter 3].

Because of the limited usefulness of those theorems, we prefer to illustrate the method with a very simple example to get the idea across. The example concerns the incomplete gamma function,

γ(a, x) = ∫_{0}^{x} e^{−t}t^{a−1} dt.

We will assume that a and x are both real and positive, in which case the integral is certainly well defined. We are interested in its behavior for large positive values of x.

The asymptotic behavior of γ(a, x) as x → ∞ is most conveniently studied by means of the related function

Γ(a, x) = Γ(a) − γ(a, x) = ∫_{x}^{∞} e^{−t}t^{a−1} dt.

Integrating by parts once, we get Γ(a, x) = e^{−x}x^{a−1} + (a − 1)Γ(a − 1, x). Repeating this process m − 1 more times, we find

Γ(a, x) = e^{−x}x^{a−1} ∑_{n=0}^{m−1} \frac{Γ(a)}{Γ(a − n)}x^{−n} + \frac{Γ(a)}{Γ(a − m)} ∫_{x}^{∞} e^{−t}t^{a−m−1} dt.

Now,

∫_{x}^{∞} e^{−t}t^{a−m−1} dt ≤ x^{a−m−1} ∫_{x}^{∞} e^{−t} dt = e^{−x}x^{a−m−1},

provided m > a − 1. Thus we obtain the asymptotic expansion of the Γ-function,

Γ(a, x) = e^{−x}x^{a−1} ∑_{n=0}^{∞} \frac{Γ(a)}{Γ(a − n)}x^{−n} as x → ∞.

The infinite sum actually converges in the usual sense, so in this case there is no need to distinguish between the asymptotic approximation and the function itself.

Exercises

1. Verify that the asymptotic expansion of γ(a, x) remains valid as x → ∞ in the complex plane in any sector S_{Δ} = {x ∈ C : | arg(x)| ≤ Δ} with Δ < 3π/2.

2. The function γ^{∗}(a, x) = x^{−a}γ(a, x)/Γ(a) is a single-valued analytic function of a and x, possessing no finite singularities. Use Mathematica to draw the graph of γ^{∗} on the rectangle {(a, x) : |a| ≤ 5, |x| ≤ 4}. Compare the result with [3, Figure 6.3].

3. Use the asymptotic expansion of the Γ-function obtained above to derive the asymptotic expansion of the complementary error function,

erfc(x) = \frac{2}{√π} ∫_{x}^{∞} e^{−t^{2}} dt = \frac{e^{−x^{2}}}{π} ∑_{n=0}^{∞} Γ(n + \frac{1}{2} ) \frac{(−1)^{n}}{x^{2n+1}} as x → ∞.

4. Prove that

∫_{x}^{∞} e^{it}t^{−a} dt = \frac{ie^{ix}}{x^{a}} ∑_{n=0}^{∞} \frac{Γ(a + n)}{Γ(a)(ix)^{n}} as x → ∞,

for any complex a with Re a > 0. Deduce from this result the asymptotic expansions of the Fresnel integrals C(x) = ∫_{0}^{x} cos( \frac{1}{2} πt^{2}) dt and S(x) = ∫_{0}^{x} sin( \frac{1}{2} πt^{2}) dt as x → ∞.

5. Use the method of integration by parts to derive Euler’s asymptotic expansion ∑_{n}(−1)^{n}n!x^{n} from the integral ∫_{0}^{∞} e^{−t}(1 + xt)^{−1}dt.

2 Laplace Integrals

In the example of the preceding section, it was possible to obtain a complete asymptotic expansion by integrating by parts an indefinite number of times. Often, it is possible to integrate only a finite number of times, and the process then leads to a finite expansion. The simplest such case occurs for Laplace-type integrals, ∫_{a}^{b} e^{−xt}f(t) dt, where the interval (a, b) is a finite segment of the real axis. We investigate the asymptotic behavior of these integrals as x → ∞.

Theorem 1 If f ∈ C^{m}([a, b]), then the integral L_{a,b}(f)(x) = ∫_{a}^{b} e^{−xt}f(t) dt has an m-term asymptotic expansion in powers of x^{−1},

EL_{a,b}(f)(x) = e^{−ax} ∑_{n=0}^{m−1} f^{(n)}(a)x^{−n−1} as x → ∞.

The result remains true if b → ∞, provided f(t) = O(e^{αt}) as t → ∞ for some constant α.

Proof. Let h_{0}(t) = e^{−xt}, and define h_{−1}, …, h_{−m} recursively,

h_{−n−1}(t) = − ∫_{t}^{∞} h_{−n}(s) ds = (−1)^{n+1}e^{−xt}x^{−n−1}, n = 0, …, m − 1.

After m partial integrations, we have

\begin{aligned}L_{a,b}(f)(x) & = ∑_{n=0}^{m−1}(−1)^{n}[h_{−n−1}(b)f^{(n)}(b) − h_{−n−1}(a)f^{(n)}(a)] + R_{m}(x) \\ & = e^{−ax} ∑_{n=0}^{m−1} f^{(n)}(a)x^{−n−1} + B_{m}(x) + R_{m}(x), \end{aligned}

where B_{m}(x) = −e^{−bx} ∑_{n=1}^{m−1} f^{(n)}(b)x^{−n−1} and R_{m}(x) = (−1)^{m} ∫_{a}^{b} h_{−m}(t)f^{(m)}(t) dt.

Now, e^{ax}B_{m}(x) = e^{−(b−a)x}O(x^{−m}) = o(x^{−m}) as x → ∞, for any positive integer m. Furthermore, f^{(m)} is bounded, and |h_{−m}(t)| ≤ e^{−xt}x^{−m} on [a, b], so e^{ax}R_{m}(x) = O(x^{−m−1}) as x → ∞. The theorem follows.

If a = 0 and the limit of ∫_{0}^{b} exists as b → ∞, then L_{a,b}(f)(x) is actually the Laplace transform of f,

L_{0,∞}(f)(x) = ∫_{0}^{∞} e^{−xt}f(t) dt.

For L_{0,∞}(f)(x) to exist for x sufficiently large, we must require that f(t) = O(e^{αt}) as t → ∞ for some constant α. The asymptotic expansion of L_{0,∞}(f)(x) is obtained by setting a = 0 in Theorem 1,

EL_{0,∞}(f)(x) = ∑_{n=1}^{m−1} f^{(n)}(0)x^{−n−1} as x → ∞.

A significant extension of Theorem 1 may be based on the following lemma.

Lemma 1 Let f and g be two given functions, and let L_{a,b}(f)(x) = ∫_{a}^{b} e^{−xt}f(t) dt and L_{a,b}(g)(x) = ∫_{a}^{b} e^{−xt}g(t) dt. Suppose that g(t) > 0 for all t ∈ (a, b). If e^{δx}L_{a,b}(g)(x) → ∞ as x → ∞ for all δ > 0 and if f(t) = o(g(t)) as t ↓ a, then L_{a,b}(f)(x) = o(L_{a,b}(g)(x)) as x → ∞.

Proof. For any ϵ > 0, there exists a δ > 0 such that |f(t)| ≤ ϵg(t) for all t ∈ [a, a + δ]. With this δ, we write L_{a,b}(f)(x) = A(x)+B(x), where A(x) = ∫_{a}^{a+δ} e^{−xt}f(t) dt and B(x) = ∫_{a+δ}^{b} e^{−xt}f(t) dt. Then |A(x)| ≤ ϵ ∫_{a}^{a+δ} e^{−xt}g(t) dt ≤ ϵL_{a,b}(g)(x) and |B(x)| ≤ Ce^{−δx} as x → ∞, so

\frac{|L_{a,b}(f)(x)|}{L_{a,b}(g)(x)} ≤ ϵ + \frac{C}{e^{δx}L_{a,b}(g)(x)}.

The second term becomes vanishingly small as x → ∞, so the upper bound is less than 2ϵ for all sufficiently large x.

The lemma enables us to analyze the behavior of asymptotic sequences and asymptotic expansions under Laplace-type transforms.

Theorem 2 Let f_{n} (n = 0, 1, …) be given positive-valued functions on the interval (a, b), and let L_{a,b}(f_{n})(x) = ∫_{a}^{b} e^{−xt}f_{n}(t) dt. (i) If {f_{n}(t) : n = 0, 1, …} is an asymptotic sequence for t ↓ a and e^{δx}L_{a,b}(f_{n})(x) → ∞ as x → ∞ for all δ > 0 and for each n, then {L_{a,b}(f_{n})(x) : n = 0, 1, …} is an asymptotic sequence for x → ∞. (ii) Furthermore, if a function f has an asymptotic expansion of the form Ef(t) = ∑_{n} a_{n}f_{n}(t) as t ↓ a, then the integral L_{a,b}(f) = ∫_{a}^{b} e^{−xt}f(t) dt has an asymptotic expansion of the form EL_{a,b}(f)(x) = ∑_{n} a_{n}L_{a,b}(f_{n})(x) as x → ∞.

Proof. The assertion (i) is an immediate consequence of Lemma 1; (ii) follows likewise from Lemma 1 if one replaces f by f^{(m)} = f − ∑_{n=0}^{m−1} a_{n}f_{n} and g by f_{m−1}. Here, m can be any integer.

The most important application of this theorem is for the Laplace transform (a = 0 and b = ∞) and f_{n}(t) = t^{λ_{n}−1}, where 0 < λ_{0} < λ_{1} < · · ·. If f has an asymptotic approximation

Ef(t) = ∑_{n=0}^{m−1} a_{n}t^{λ_{n}−1} as t ↓ 0,

then the asymptotic approximation of its Laplace transform L_{0,∞}(f) = ∫_{0}^{∞} e^{−xt}f(t) dt is

EL_{0,∞}(f)(x) = ∑_{n=0}^{m−1} a_{n}Γ(λ_{n})x^{−λ_{n}} as x → ∞.

The results of this section can be generalized in many ways. For example, results similar to those of Theorems 1 and 2 hold if x is a complex parameter and x → ∞ in the complex plane within any sector S_{Δ} = \{x ∈ C : | arg(x)| ≤ Δ\} with Δ < \frac{1}{2}π. We refer the reader to [1, Section 2.2] for details.

3 Laplace’s Method

The essence of Theorem 2 is that, under certain circumstances, the asymptotic behavior of a Laplace-type integral as x → ∞ is determined by the asymptotic behavior of its integrand as t ↓ a. This result extends to more general situations. For example, consider the integral ∫_{a}^{b} g(t)e^{xh(t)} dt, where (a, b) is a finite segment of the real axis, x is a large positive parameter, and h is real-valued. (The function g may be complex-valued.) If h(t) has a maximum at t = τ and h(t) < h(τ) for all t ≠ τ, then the modulus of the integrand will have a maximum at a point near τ, and most of the contribution to the integral will arise from the immediate vicinity of this maximum. The integral can then be evaluated approximately by expanding both g(t) and h(t) in the neighborhood of t = τ. This idea is central in Laplace’s method.

Without loss of generality, we may assume that h(t) reaches its maximum at one of the endpoints and at no other point of the interval; if necessary, we break up the interval of integration in a finite number of subintervals. Accordingly, we assume that h(t) reaches its maximum at a and that h(t) < h(a) for all t ∈ (a, b].

Theorem 3 Let h be real-valued and continuous at a, and let g and h be such that the integral ∫_{a}^{b} e^{xh(t)}g(t) dt is well defined for all sufficiently large positive x. Suppose that there exist ϵ > 0 and η > 0, such that (i) g(t) = (t − a)^{λ−1}g_{1}(t) for some λ > 0, where g_{1} is m times continuously differentiable on [a, a + η]; (ii) h′(t) = −(t − a)^{ρ−1}h_{1}(t) for some ρ > 0, where h_{1} is m times continuously differentiable on [a, a+η] and h_{1}(t) > 0 for all t ∈ [a, a+η]; and (iii) h(t) < h(a)−ϵ for all t ∈ (a+η, b). Then ∫_{a}^{b} e^{xh(t)}g(t) dt has an m-term asymptotic approximation with respect to the asymptotic sequence {x^{−(n+λ)/ρ} : n = 0, 1, …},

E ∫_{a}^{b} e^{xh(t)}g(t) dt

= e^{xh(a)} ∑_{n=0}^{m−1} \frac{1}{n!} [\frac{d^{n}}{dt^{n}} g(t)\frac{(h(a) − h(t))^{1−λ/ρ}}{−h′(t)} ]_{t=a} Γ (\frac{n + λ}{ρ}) x^{−(n+λ)/ρ} as x → ∞.

Proof. We write the integral as a sum of two terms, ∫_{a}^{b} e^{xh(t)}g(t) dt = A(x) + B(x), where A(x) = ∫_{a}^{a+η} e^{xh(t)}g(t) dt and B(x) = ∫_{a+η}^{b} e^{xh(t)}g(t) dt. On the interval [a, a + η], we introduce a new variable of integration,

u = (h(a) − h(t))^{1/ρ}, t ∈ [a, a + η].

Then

u^{ρ} = − ∫_{a}^{t} h′(s) ds = ∫_{a}^{t}(s − a)^{ρ−1}h_{1}(s) ds = (t − a)^{ρ} ∫_{0}^{1} y^{ρ−1}h_{1}(a + y(t − a)) dy.

The last integral is m + 1 times continuously differentiable, positive, and increasing as a function of t. With U = (h(a) − h(a + η))^{1/ρ}, we have thus an m + 1 times continuously differentiable mapping from the interval [a, a + η] to the interval [0, U], whose inverse exists and is also m + 1 times continuously differentiable. Then,

A(x) = e^{xh(a)} ∫_{0}^{U} e^{−xu^{ρ}}u^{λ−1}k(u) du,

where k(u) = g(t)u^{1−λ}dt/du; k is m times continuously differentiable on [0, U].

Let l_{0}(u) = e^{−xu^{ρ}}u^{λ−1}, and let l_{−1}, …, l_{−m} be defined recursively,

l_{−n−1}(u) = − ∫_{u}^{∞} l_{−n}(v) dv, n = 0, …, m − 1.

The general expression for l_{−n−1}(u) is

l_{−n−1}(u) = \frac{(−1)^{n+1}}{n!} ∫_{u}^{∞} e^{−xv^{ρ}}(v − u)^{n}v^{λ−1} dv, n = 0, …, m − 1.

After m partial integrations, we have

A(x) = e^{xh(a)} ∑_{n=0}^{m−1}(−1)^{n}[l_{−n−1}(U)k^{(n)}(U) − l_{−n−1}(0)k^{(n)}(0)] + R_{m}(x),

where the remainder is R_{m}(x) = e^{xh(a)}(−1)^{m} ∫_{0}^{U} l_{−m}(u)k^{(m)}(u) du. Now,

l_{−n−1}(0) = \frac{(−1)^{n+1}}{n!} ∫_{0}^{∞} e^{−xv^{ρ}}v^{n+λ−1} dv = \frac{(−1)^{n+1}}{ρn!} Γ (\frac{n + λ}{ρ}) x^{−(n+λ)/ρ},

while l_{−n−1}(U) = O(e^{−xU^{ρ}}) for n = 0, …, m − 1. Furthermore,

|R_{m}(x)| ≤ \frac{Ce^{xh(a)}}{(m − 1)!} ∫_{0}^{U} ∫_{u}^{∞} e^{−xv^{ρ}}(v − u)^{m−1}v^{λ−1} dv du
≤ \frac{Ce^{xh(a)}}{m!} ∫_{0}^{∞} e^{−xv^{ρ}}v^{m+λ−1} dv = \frac{Ce^{xh(a)}}{ρm!} Γ (\frac{m + λ}{ρ}) x^{−(m+λ)/ρ},

so R_{m}(x) = e^{xh(a)}O(x^{−(m+λ)/ρ}). Since O(e^{−xU^{ρ}}) = o(x^{−(m−1+λ)/ρ}) as x → ∞, for any positive integer m, it follows that A has the m-term asymptotic expansion

EA(x) = e^{xh(a)} ∑_{n=0}^{m−1} \frac{k^{(n)}(0)}{ρn!}Γ (\frac{n + λ}{ρ}) x^{−(n+λ)/ρ}.

We complete the proof by observing that B(x) = e^{xh(a)}O(e^{−ϵx}) and e^{−ϵx} = o(x^{−(m−1+λ)/ρ}) as x → ∞, for any positive integer m.

It follows from Theorem 3 that, if

g(t) ∼ α(t − a)^{λ−1}, h′(t) ∼ −β(t − a)^{ρ−1} as t ↓ a,

then

∫_{a}^{b} e^{xh(t)}g(t) dt ∼ \frac{α}{ρ}Γ (\frac{λ}{ρ}) (\frac{ρ}{βx})^{λ/ρ} e^{xh(a)} as x → ∞.

The integral studied by Laplace corresponds to the special case λ = 1, ρ = 2.

The results of this section remain true if x is a complex variable and x → ∞ in any sector S_{Δ} with Δ < \frac{1}{2}π.

Exercises

1. The logarithmic derivative of the Γ function is given by

ψ(x) = \frac{Γ′(x)}{Γ(x)} = log x + ∫_{0}^{∞} e^{−xt} [\frac{1}{t} − \frac{1}{1 − e^{−t}} ] dt.

Obtain an asymptotic expansion of ψ(x) as x → ∞. Justify that this expansion may be integrated term by term and obtain an asymptotic expansion of log Γ(x) (modulo an arbitrary constant) as x → ∞. Use the duplication formula for the Γ-function, Γ(2x) = (2π)^{−1/2}2^{2x−1/2}Γ(x)Γ(x + \frac{1}{2} ), to fix the arbitrary constant and establish the asymptotic expansion

log Γ(x) = (x − \frac{1}{2} ) log x − x + \frac{1}{2} log(2π) + ∑_{n=1}^{∞} \frac{B_{2n}}{2n(2n − 1)} x^{1−2n} as x → ∞.

Here, B_{n} is the nth Bernoulli number; cf. [3, Section 23.1].

2. When x is positive, the logarithmic integral li(x) is defined by

li(x) = ∫_{0}^{x} \frac{dt}{log t} ,

the integral being a Cauchy principal value (P) when x > 1. If we put x = e^{a} and t = e^{a−v}, the definition becomes

li(x) = e^{a}P ∫_{0}^{∞} e^{−v} \frac{dv}{a − v} ,

when a > 0. Apply the method of integration by parts to establish the asymptotic expansion

e^{−a}li(e^{a}) = ∑_{n=0}^{∞} n! a^{−n+1} as a → ∞.

Verify that the smallest term in the asymptotic expansion is the mth term, where m is the largest integer less than or equal to a. Verify that the difference between e^{−a}li(e^{a}) and the first m terms of the asymptotic expansion (i.e., the remainder after m partial integrations, R_{m}(a)) is given by e^{−a}(α − 1/3)(2π/m)^{1/2} as a → ∞, where α = a − m (0 ≤ α < 1). Show that R_{m−1}(a) and R_{m+1}(a) have opposite signs as a → ∞.

3. The modified Bessel function of order ν has the integral representation

I_{ν}(x) = \frac{( \frac{1}{2} x)^{ν}}{Γ(ν + \frac{1}{2} )Γ( \frac{1}{2} )} ∫_{−1}^{1} e^{−xt}(1 − t^{2})^{ν−1/2} dt, ν > − \frac{1}{2} .

Apply the divide-and-conquer technique to separate the effects of the endpoints and establish the asymptotic expansion

I_{ν}(x) = \frac{e^{x}}{(2πx)^{1/2}} ∑_{n=0}^{∞} \frac{(−1)^{n}Γ(ν + \frac{1}{2} + n)}{2^{n}n! Γ(ν + \frac{1}{2} − n)} x^{−n} as x → ∞.

Cf. [4, Section 7.23].

4 Fourier Integrals

Next, we consider Fourier-type integrals, ∫_{a}^{b} e^{ixt}f(t) dt, where (a, b) is a finite segment of the real axis. If f is integrable, these integrals exist for all real x. We are interested in their asymptotic behavior as x → ∞.

Theorem 4 If f ∈ C^{m}([a, b]), then the integral F_{a,b}(f)(x) = ∫_{a}^{b} e^{ixt}f(t) dt has an m-term asymptotic expansion in powers of x^{−1},

EF_{a,b}(f)(x) = ∑_{n=0}^{m−1} i^{n−1}x^{−(n+1)} [e^{ibx}f^{(n)}(b) − e^{iax}f^{(n)}(a)] as x → ∞.

The result remains true if a → −∞ (b → ∞), provided that f^{(n)}(t) → 0 as t → −∞ (t → ∞) for n = 0, …, m − 1, and provided that f^{(m)} is integrable over (a, b).

Proof. The remainder after m partial integrations is R_{m}(x) = −(ix)^{−m} ∫_{a}^{b} e^{ixt}f^{(m)}(t) dt. Because f^{(m)} is of bounded variation, it follows from the Riemann–Lebesgue lemma [5, Section 9.41] that the integral vanishes as x → ∞. Therefore, R_{m}(x) = o(x^{−m}).

Typically for the method of integration by parts, the result of Theorem 4 relies heavily on the differentiability properties of f. Suppose that f(t) = (t − a)^{λ−1}(b − t)^{μ−1}f_{1}(t), where f_{1} ∈ C^{∞}([a, b]) and λ and μ are not integers. Then Theorem 4 gives us at best that F_{a,b}(f)(x) = O(x^{−n}), where n is the largest integer less than or equal to min(λ, μ); in particular, if λ < 1 and μ < 1, we find that F_{a,b}(f)(x) = O(1) as x → ∞. This result is not very good, as F_{a,b}(f)(x) = o(1) in this case. Fortunately, we can do better, thanks to a method developed by Erdélyi for functions f that are weakly singular near the endpoints a and b; see [1, Section 2.8].

Theorem 5 If f(t) = (t − a)^{λ−1}(b − t)^{μ−1}f_{1}(t), where 0 < λ, μ ≤ 1 and f_{1} ∈ C^{m}([a, b]), then the asymptotic expansion of the integral F_{a,b}(f)(x) = ∫_{a}^{b} e^{ixt}f(t) dt is given by

EF_{a,b}(f)(x) = ∑_{n=0}^{m−1} \frac{Γ(n + μ)}{n!}e^{πi(n−μ)/2}x^{−(n+μ)}e^{ibx} [\frac{d^{n}}{dt^{n}} (t − a)^{λ−1}f_{1}(t)]_{t=b}

− ∑_{n=0}^{m−1} \frac{Γ(n + λ)}{n!}e^{πi(n+λ−2)/2}x^{−(n+λ)}e^{iax} [\frac{d^{n}}{dt^{n}} (b − t)^{μ−1}f_{1}(t)]_{t=a} as x → ∞.

Proof. Divide and conquer. Let χ be an infinitely smooth function satisfying the conditions χ(t) = 1 for all t ∈ [a, a + η] and χ(t) = 0 for all t ∈ [b − η, b], where η < \frac{1}{2}(b − a). Write the integral as a sum of two integrals, F_{a,b}(f)(x) = A(x) + B(x), where

A(x) = ∫_{a}^{b−η} e^{ixt}(t − a)^{λ−1} [χ(t)(b − t)^{μ−1}f_{1}(t)] dt,
B(x) = ∫_{a+η}^{b} e^{ixt}(b − t)^{μ−1} [(1 − χ(t))(t − a)^{λ−1}f_{1}(t)] dt.

Thus, the singularities at the endpoints have been separated; the integrand of A has inherited the singularity at a, but goes to zero smoothly at b−η; the integrand of B has inherited the singularity at b, but goes to zero smoothly at a + η.

We discuss the asymptotic behavior of A(x) as x → ∞. The integral is of the form

A(x) = ∫_{a}^{b−η} e^{ixt}(t − a)^{λ−1}f_{2}(t) dt,

where f_{2} ∈ C^{m}([a, b − η]), f_{2}^{(n)}(a) = [(d^{n}/dt^{n})(b − t)^{μ−1}f_{1}(t)]_{t=a} and f_{2}^{(n)}(b − η) = 0 for n = 0, …, m − 1.

Let h_{0}(t) = e^{ixt}(t − a)^{λ−1}, and let h_{−1}, …, h_{−m} be defined recursively,

h_{−n−1}(t) = − ∫_{t} h_{−n}(z) dz, n = 0, …, m − 1,

where the integral connects the point t ∈ [a, b] with the point at infinity in the complex plane along the directed path {z ∈ C : z = t + iτ, 0 ≤ τ < ∞}. The general expression for h_{−n−1}(t) is

h_{−n−1}(t) = \frac{(−1)^{n+1}}{n!} ∫_{t} e^{ixz}(z − t)^{n}(z − a)^{λ−1} dz, n = 0, …, m − 1.

The integral converges absolutely.

After m partial integrations, we have

A(x) = − ∑_{n=0}^{m−1}(−1)^{n}h_{−n−1}(a)f_{2}^{(n)}(a) + R_{m}(x),

where the remainder is R_{m}(x) = (−1)^{m} ∫_{a}^{b} h_{−m}(t)f_{2}^{(m)}(t) dt. Now,

h_{−n−1}(a) = \frac{(−1)^{n+1}}{n!} ∫_{a}^{a+i∞} e^{ixz}(z − a)^{n+λ−1} dz = (−1)^{n} \frac{Γ(n + λ)}{n!} \frac{e^{πi(n+λ−2)/2}}{x^{n+λ}} e^{iax},

for n = 0, …, m−1. To estimate R_{m}, we need an estimate of h_{−m}(t) for t ∈ [a, b]. It is most readily obtained from the general expression given above. Along the path of integration we have |z − a|^{λ−1} ≤ (t − a)^{λ−1} for all λ ∈ (0, 1], so

|h_{−m}(t)| ≤ \frac{(t − a)^{λ−1}}{(m − 1)!} ∫_{t} |e^{ixz}| |z − t|^{m−1} |dz| ≤ (t − a)^{λ−1}x^{−m}.

Hence, R_{m} = O(x^{−m}) as x → ∞. If λ = 1, we can replace this estimate by R_{m} = o(x^{−m}), by virtue of the Riemann–Lebesgue lemma; cf. the proof of Theorem 4. Thus, for all λ ∈ (0, 1] we find that the asymptotic expansion of A is given by

EA(x) = − ∑_{n=0}^{m−1} \frac{Γ(n + λ)}{n!} \frac{e^{πi(n+λ−2)/2}}{x^{n+λ}} e^{iax}f_{2}^{(n)}(a) as x → ∞.

Similarly,

EB(x) = ∑_{n=0}^{m−1} \frac{Γ(n + μ)}{n!} \frac{e^{πi(n−μ)/2}}{x^{n+μ}} e^{ibx}f_{3}^{(n)}(b) as x → ∞,

where f_{3}^{(n)}(b) = [(d^{n}/dt^{n})(t − a)^{λ−1}f_{1}(t)]_{t=b}.

The essence of Theorems 4 and 5 is that, under certain circumstances, the asymptotic behavior of a Fourier integral as x → ∞ is determined by the behavior of its integrand near the endpoints of the interval of integration. As x increases, the integrand oscillates more and more rapidly, and the contribution from each oscillation tends to zero, except near the endpoints.

5 Stokes’ Method of Stationary Phase

We now consider the integral ∫_{a}^{b} e^{ixh(t)}g(t) dt, where (a, b) is a finite segment of the real axis, x is a large positive parameter, and h is real-valued. (The function g may be complexvalued.) In general, the rapid oscillations of e^{ixh(t)} tend to cancel large contributions to the integral, but this cancellation does not occur near the endpoints and near the stationary points of h—that is, those points where h′ vanishes. Moreover, if there are any stationary points, their contributions tend to be more important than the contributions from the endpoints. Stokes’ method of stationary phase appraises the contribution of the stationary points to the integral.

Assuming that h has only a finite number of stationary points, we may break up the integral in a finite number of integrals, in each of which h is monotonic, and we may assume h(t) to be increasing. Thus, we shall consider integrals in which h′(t) > 0 for all t ∈ (a, b) and where a and b are either ordinary points (where h′(t) is positive) or stationary points (where h′(t) vanishes to some positive, possibly fractional, order).

Theorem 6 Let g and h be such that the integral ∫_{a}^{b} e^{ixh(t)}g(t) dt is well defined for all sufficiently large positive x. Let h be real-valued and differentiable, such that h′(t) = (t − a)^{ρ−1}(b − t)^{σ−1}h_{1}(t), where ρ, σ ≥ 1, h_{1} ∈ C^{m}([a, b]), and h_{1}(t) > 0 for all t ∈ [a, b]. If g(t) = (t−a)^{λ−1}(b−t)^{μ−1}g_{1}(t), where 0 < λ, μ ≤ 1 and g_{1} ∈ C^{m}([a, b]), then the asymptotic expansion of ∫_{a}^{b} e^{ixh(t)}g(t) dt is given by the expression

E ∫_{a}^{b} e^{ixh(t)}g(t) dt = ∑_{n=0}^{m−1} \frac{β_{n}}{n!} Γ (\frac{n + μ}{σ}) e^{−πi(n+μ)/(2σ)}x^{−(n+μ)/σ}e^{ixh(b)}

+ ∑_{n=0}^{m−1} \frac{α_{n}}{n!} Γ (\frac{n + λ}{ρ}) e^{πi(n+λ)/(2ρ)}x^{−(n+λ)/ρ}e^{ixh(a)} as x → ∞,

where

α_{n} = ρ^{n} [(\frac{(h(t) − h(a))^{1−1/ρ}}{h′(t)} \frac{d}{dt})^{n} g(t)\frac{(h(t) − h(a))^{1−λ/ρ}}{h′(t)} ]_{t=a} ,

β_{n} = (−σ)^{n} [(\frac{(h(b) − h(t))^{1−1/σ}}{h′(t)} \frac{d}{dt})^{n} g(t)\frac{(h(b) − h(t))^{1−μ/σ}}{h′(t)} ]_{t=b} .

Proof. Divide and conquer. Separate the effects of the endpoints of the integral as in the proof of Theorem 5: ∫_{a}^{b} e^{ixh(t)}g(t) dt = A(x) + B(x), where

A(x) = ∫_{a}^{b−η} e^{ixh(t)}(t − a)^{λ−1} [χ(t)(b − t)^{μ−1}g_{1}(t)] dt,

B(x) = ∫_{a+η}^{b} e^{ixh(t)}(b − t)^{μ−1} [(1 − χ(t))(t − a)^{λ−1}g(t)] dt.

We discuss the asymptotic behavior of A(x) as x → ∞. The integral is of the form

A(x) = ∫_{a}^{b−η} e^{ixh(t)}(t − a)^{λ−1}g_{2}(t) dt,

where g_{2} ∈ C^{m}([a, b − η]), g_{2}^{(n)}(a) = [(d/dt)^{n}(b − t)^{μ−1}g_{1}(t)]_{t=a} and g_{2}^{(n)}(b − η) = 0 for n = 0, …, m − 1. We introduce a new variable of integration,

u = (h(t) − h(a))^{1/ρ}, t ∈ [a, b − η].

Then,

u^{ρ} = ∫_{a}^{t} h′(s) ds = ∫_{a}^{t}(s − a)^{ρ−1}(b − s)^{σ−1}h_{1}(s) ds
= (t − a)^{ρ} ∫_{0}^{1} y^{ρ−1}(b − a − y(t − a))^{σ−1}h_{1}(a + y(t − a)) dy.

The last integral is m + 1 times continuously differentiable, positive, and increasing as a function of t. With U = (h(b − η) − h(a))^{1/ρ}, we have thus an m + 1 times continuously differentiable mapping from the interval [a, b − η] to the interval [0, U], whose inverse exists and is also m + 1 times continuously differentiable. Then,

A(x) = e^{ixh(a)} ∫_{0}^{U} e^{ixu^{ρ}}u^{λ−1}k(u) du,

where k(u) = (t − a)^{λ−1}u^{1−λ}g_{2}(t)dt/du; k is m times continuously differentiable on [0, U].

Let l_{0}(t) = e^{ixu^{ρ}}u^{λ−1}, and let l_{−1}, …, l_{−m} be defined recursively,

l_{−n−1}(u) = − ∫_{u} l_{−n}(z) dz, n = 0, …, m − 1,

where the integral connects the point u ∈ [0, U] with the point at infinity in the complex plane along the directed path {z ∈ C : z = u + τe^{πi/(2ρ)}, 0 ≤ τ < ∞}. Thus,

l_{−n−1}(u) = \frac{(−1)^{n+1}}{n!} ∫_{u} e^{ixz^{ρ}}(z − u)^{n}z^{λ−1} dz, n = 0, …, m − 1.

After m partial integrations, we have

A(x) = −e^{ixh(a)} ∑_{n=0}^{m−1}(−1)^{n}l_{−n−1}(0)k^{(n)}(0) + R_{m}(x),

where the remainder is R_{m}(x) = e^{ixh(a)}(−1)^{m} ∫_{0}^{U} l_{−m}(u)k^{(m)}(u) du. Now,

l_{−n−1}(0) = \frac{(−1)^{n+1}}{n!} ∫_{u} e^{ixz^{ρ}}z^{n+λ−1} dz = \frac{(−1)^{n+1}}{ρn!} e^{πi(n+λ)/(2ρ)}Γ (\frac{n + λ}{ρ}) x^{−(n+λ)/ρ},

for n = 0, …, m − 1.

Estimating the remainder R_{m} is a bit more involved. First, consider the case 0 < λ < 1. Then |z|^{λ−1} ≤ u^{λ−1} along the path of integration, so

|l_{−m}(u)| ≤ \frac{u^{λ−1}}{(m − 1)!} ∫_{u} |e^{ixz^{ρ}}| |z − u|^{m−1} |dz|.

Also, from the identity

ixz^{ρ} + x|z − u|^{ρ} = iρx ∫_{0}^{u}(s + |z − u|e^{πi/(2ρ)})^{ρ−1} ds,

together with the fact that the imaginary part of the integral is positive, we deduce that Re(ixz^{ρ}) ≤ −x|z − u|^{ρ}, so

|l_{−m}(u)| ≤ \frac{u^{λ−1}}{ρ(m − 1)!}Γ (\frac{m}{ρ}) x^{−m/ρ}.

Hence, R_{m}(x) = O(x^{−m/ρ}) as x → ∞.

If λ = 1 and ρ > 1, we proceed in two steps. First, given ϵ > 0 arbitrarily small, we choose δ > 0 such that

\frac{1}{ρm!}Γ (\frac{m}{ρ}) ∫_{0}^{δ} |k^{(n)}(u)| du < \frac{1}{2}ϵ.

Then,

l_{−m}(u) = \frac{(−1)^{m}}{(m − 1)!} ∫_{u} e^{ixz^{ρ}}(z − u)^{m−1} dz = \frac{(−ue^{πi/(2ρ)})^{m}}{(m − 1)!} ∫_{0}^{∞} e^{ixu^{ρ}(1+τe^{πi/(2ρ)})^{ρ}}τ^{m−1} dτ,

so l_{−m}(u) = u^{m}O((xu^{ρ})^{−m}) as xu^{ρ} → ∞. Consequently, l_{−m}(u) = O(x^{−m}) as x → ∞, uniformly in u when u ≥ δ. Hence,

x^{m/ρ} ∫_{δ}^{U} |l_{−m}(u)| |k^{(m)}(u)| du < \frac{1}{2}ϵ,

for all sufficiently large x, and therefore R_{m}(x) = o(x^{−m/ρ}) as x → ∞.

The only remaining case is λ = 1, ρ = 1. But this case is already covered by Theorem 5. So we find that, for all ρ ≥ 1 and 0 < λ ≤ 1, the asymptotic expansion of A is given by

EA(x) = e^{ixh(a)} ∑_{n=0}^{m−1} \frac{α_{n}}{n!} Γ (\frac{n + λ}{ρ}) e^{πi(n+λ)/(2ρ)}x^{−(n+λ)/ρ} as x → ∞,

where α_{n} = k^{(n)}(0)/ρ. The coefficient α_{n} can be expressed in terms of the values of g and h and their derivatives at a; the result is given in the statement of the theorem.

A similar result holds for B(x). The integral is of the form

B(x) = ∫_{a+η}^{b} e^{ixh(t)}(b − t)^{μ−1}g_{3}(t) dt,

where g_{3} ∈ C^{m}([a + η, b]), g_{3}^{(n)}(a + η) = 0 and g_{3}^{(n)}(b) = [(d/dt)^{n}(t − a)^{λ−1}g_{1}(t)]_{t=b} for n = 0, …, m − 1. We introduce the new variable of integration,

v = (h(b) − h(t))^{1/σ}, t ∈ [a + η, b],

and put l(v) = (b−t)^{μ−1}v^{1−μ}g_{3}(t)dt/dv. In the repeated integrals of e^{ixv^{σ}}v^{μ−1} we integrate along the directed path {z ∈ C : z = v + τe^{−πi/(2σ)}} and obtain by a process similar to that used in the case of A(x) that, for all σ ≥ 1 and 0 < μ ≤ 1, the asymptotioc expansion of B is given by

EB(x) = e^{ixh(b)} ∑_{n=0}^{m−1} \frac{β_{n}}{n!} Γ (\frac{n + μ}{σ}) e^{−πi(n+μ)/(2σ)}x^{−(n+μ)/σ} as x → ∞,

where β_{n} = −l^{(n)}(0)/σ. The coefficient β_{n} can be expressed in terms of the values of g and h and their derivatives at b; the result is given in the statement of the theorem.

Let us consider in more detail the case where g is regular at both endpoints a and b (i.e., λ = μ = 1). If a is a stationary point of order one and b an ordinary point (i.e., h′(t) = (t − a)h_{1}(t) with h_{1}(t) > 0 for all t ∈ [a, b], so h has a minimum at a), then the leading term in the asymptotic expansion given in Theorem 6 comes from the point a and is O(x^{1/2}), while the next term comes from the point b and is O(x^{−1}). From the theorem we obtain

α_{0} = [g(t)\frac{(h(t) − h(a))^{1/2}}{h′(t)} ]_{t=a} = \frac{g(a)}{\sqrt{2h′′(a)}}.

Hence,

∫_{a}^{b} e^{ixh(t)}g(t) dt ∼ \sqrt{\frac{π}{2xh′′(a)}} g(a)e^{ixh(a)+πi/4} as x → ∞.

Similarly, if a is an ordinary point and b a stationary point of order one (i.e., h′(t) = (b − t)h_{1}(t) with h_{1}(t) > 0 for all t ∈ [a, b], so h has a maximum at b), then

β_{0} = [g(t)\frac{(h(b) − h(t))^{1/2}}{h′(t)} ]_{t=b} = \frac{g(b)}{\sqrt{−2h′′(b)}}

and

∫_{a}^{b} e^{ixh(t)}g(t) dt ∼ \sqrt{\frac{π}{2x|h′′(b)|}} g(b)e^{ixh(b)+πi/4} as x → ∞.

In general, if both endpoints a and b are ordinary points, where h′ does not vanish, but h has an extremum (either a minimum or a maximum) at some point τ in the interior of the interval (a, b), then we find, by combining the leading terms in the asymptotic expansions of the integrals over [a, τ] and [τ, b], that the leading term in the asymptotic expansion of ∫_{a}^{b} e^{ixh(t)}g(t) dt is given by

∫_{a}^{b} e^{ixh(t)}g(t) dt ∼ \sqrt{\frac{2π}{x|h′′(τ)|}} g(τ)e^{ixh(τ)+πi/4} as x → ∞.

In this case, the contributions from the endpoints are O(x^{−1}) or smaller as x → ∞. This observation confirms our earlier remark that the contributions from stationary points are generally more important than those from regular endpoints.

Exercises

1. The Bessel function of order ν has the integral representation

J_{ν}(x) = \frac{( \frac{1}{2} x)^{ν}}{Γ(ν + \frac{1}{2} )Γ( \frac{1}{2} )} ∫_{−1}^{1} e^{ixt}(1 − t^{2})^{ν−1/2} dt, ν > − \frac{1}{2} .

Apply the divide-and-conquer technique to separate the effects of the endpoints and establish the asymptotic expansion

J_{ν}(x) = (\frac{2}{πx})^{1/2} [cos (x − \frac{π}{4} − \frac{νπ}{2} ) ∑_{n=0}^{∞} \frac{(−1)^{n}Γ(ν + \frac{1}{2} + 2n)}{2^{2n}(2n)!Γ(ν + \frac{1}{2} − 2n)} x^{−2n}

− sin (x − \frac{π}{4} − \frac{νπ}{2} ) ∑_{n=0}^{∞} \frac{(−1)^{n}Γ(ν + \frac{3}{2} + 2n)}{2^{2n+1}(2n + 1)!Γ(ν − \frac{1}{2} − 2n)} x^{−2n−1}] as x → ∞.

Cf. [4, Section 7.21].

2. Prove that

∫_{0}^{1} e^{ixt^{3}} dt = Γ (\frac{4}{3}) e^{πi/6}x^{−1/3} − ∑_{n=0}^{∞} \frac{Γ(n + 2/3)}{Γ(−1/3)} (ix)^{−n−1}e^{ix} as x → ∞.

3. Prove that

∫_{0}^{∞} e^{ix(t^{3}/3+t)} dt = i ∑_{n=0}^{∞} \frac{(3n)!}{3^{n}n!} x^{−2n−1} as x → ∞.

4. Airy’s integral is defined by

Ai(x) = \frac{1}{π} ∫_{0}^{∞} cos(s^{3}/3 + xs) ds.

With s = x^{1/2}t and y = x^{3/2}, we obtain

Ai(y^{2/3}) = \frac{y^{1/3}}{2π} ∫_{−∞}^{∞} e^{iy(t^{3}/3+t)} dt.

Prove that

Ai(x) = \frac{e^{−(2/3)x^{3/2}}}{2πx^{1/4}} ∑_{n=0}^{∞} \frac{(−1)^{n}Γ(3n + \frac{1}{2} )}{9^{n}(2n)!} x^{−3n/2} as x → ∞.

References

[1] A. Erdélyi, Asymptotic expansions, Dover Publications, Inc., New York (1956).

[2] N. Bleistein and R. A. Handelsman, Asymptotic expansions of integrals, Holt, Rinehart and Winston, New York (1975).

[3] M. A. Abramowitz and I. A. Stegun (eds.), Handbook of mathematical functions, National Bureau of Standards (1965).

[4] G. N. Watson, Theory of Bessel functions, 2nd ed., Cambridge University Press, Cambridge, England (1944).

[5] E. T. Whittaker and G. N. Watson, A course of modern analysis, 4th ed., Cambridge University Press, Cambridge, England (1927).

Index

Airy’s integral, 15 applications of mathematics

according to Hardy, 1

Bernoulli number, 7 Bessel function, 15

Complementary error function, 2

Erdélyi, A., 8

Fourier integrals, 8–10 Fourier-type integrals, 8 Fresnel integral, 2

Gamma function logarithmic derivative of, 7

Hardy, G. H., 1

Incomplete Gamma function, 2 integration by parts, 1

Laplace integrals, 3–5 Laplace transform, 3, 4

asymptotic approximation, 4 Laplace’s method, 5–8 Laplace-type integrals, 3 logarithmic integral, 7

Method of stationary phase, 11 modified Bessel function, 8

Riemann–Lebesgue lemma, 8, 10

Stationary points, 11 Stokes’ method of stationary phase, 10–15