. 14
( 24 .)


order. However, at each stage, the grid spacing h is halved. In the lim-
it that 1/h and the order of Richardson extrapolation simultaneously
tend to in¬nity, the quadrature scheme converges.
Lyness and Ninham [190] noted this exponential dependence on 1/h
in quadrature errors nearly thirty years ago. Lyness has emphasized
that the Euler-MacLaurin series for the trapezoidal rule is closely relat-
ed to a general formula for the coe¬cients of a Fourier series. His
“FCAE” [Fourier Coe¬cient Asymptotic Expansion] takes the form
[188, 189]
’ g 2k’1 (0)
2k’1 (1)
n+k’1 g
an = (’1) (120)
(’1)M +1 1
g (2M +1) (x) sin(2π n[x ’ 1/2] )dx
(2πn)2M +1 0

plus a similar series for the sine coe¬cients. It is derived by integration-
by-parts. As noted by Lyness, it is usually divergent.
As explained in the books by Boyd [55] and Canuto et al. [91], spec-
tral methods usually employ a basis set so that the error is exponential-
ly small in 1/h or equivalently, in the number of degrees of freedom M .
Fourier series, for example, are restricted to periodic functions. Cheby-
shev polynomials give exponential convergence even for non-periodic
functions, provided only that f(x) is free of singularities on x ∈ [’1, 1].
These polynomials are de¬ned by Tn (cos[θ]) ≡ cos(n θ) so that the

ActaApplFINAL_OP92.tex; 21/08/2000; 16:16; no v.; p.63
64 John P. Boyd

expansion of f (x) as a Chebyshev series is identical, under this change
of variable, with the Fourier expansion of f (cos(θ)). The transformed
function is always periodic in θ, so the error in truncating a Chebyshev
series after M terms decreases exponentially fast with M . For Cheby-
shev and Fourier spectral methods, the power series in h is always the
trivial one with zero coe¬cients.
It follows that all asymptotic approximations to Fourier, Chebyshev,
and other spectral coe¬cients for large M are implicitly hyperasymp-
totic (Table V). One might object that this catalogue of asymptotics for
orthogonal series is out of place here because it is “beyond all orders”
[in h] only because the coe¬cients of all powers of h are zero. The
main tools for the work of Elliott, Miller, Luke, Weideman and Boyd
were steepest descents and the calculus of residues ” no explicit use
of hyperasymptotic thinking at all.
Nevertheless, Table V makes several important points. First, much
of hyperasympotics is steepest descent and the calculus of residues. In
discussing the Stieltjes function earlier, for example, we noted that one
could go beyond the superasympotic approximation by applying steep-
est descent to the error integral of the optimally-truncated power
series. Similarly, the heart of the PKKS technique of matched asymp-
totic expansions in the complex plane is the notion that the singular-
ities or other critical points closest to the real axis control the hyper-
asymptotic behavior. The asymptotic behavior of the coe¬cients of
orthogonal expansions is likewise controlled by complex singularities.
In both cases, the constant q inside the exponential factor, exp(’q/ )
or exp(’q/h), is simply the distance from the dominant singularity to
the real x-axis.3
Second, an important hyperasymptotic strategy is to isolate the
exponentially small contributions. For the large degree behavior of
Chebyshev and Fourier coe¬cients, this isolation is a free gift, the result
of choosing the sensible spectral basis ” Chebyshev or Legendre poly-
nomials for non-periodic problems, spherical harmonics for problems
on the surface of a sphere and so on. For less trivial problems, the key
to isolation is to subtract an optimally truncated asymptotic expansion.
This is the justi¬cation for applying steepest descent to the error inte-
gral for Stieltjes function, for applying Euler™s method or other sum
acceleration scheme, for Dingle™s universal terminants and for Berry™s
smoothing of Stokes phenomenon. The di¬erence between a quantity
and its superasymptotic approximation is as trivially isolated as the

For Chebyshev series, “distance from the real axis” means distance in the trans-
formed coordinate θ, the argument of the equivalent Fourier series, rather than the
argument of the Chebyshev polynomials, x = cos(θ).

ActaApplFINAL_OP92.tex; 21/08/2000; 16:16; no v.; p.64
Exponential Asymptotics

spectral coe¬cients of a Fourier series, for which the superasymptotic
approximation is zero.
Third, the asymptotics of spectral series and other numerical pro-
cesses is a relatively under-cultivated area. Can the ideas reviewed here
lead to optimal order Richardson extrapolation for numerical quadra-
ture, various classes of di¬erential equations, and so on?
Fourth, there have been some limited but important excursions
beyond conventional asymptotics in the analysis of the convergence
of spectral series. For example, Boyd™s 1982 article on the optimization
of Chebyshev methods on an unbounded domain notes that there are
two di¬erent species of contributions to the asymptotic spectral coef-
¬cients: (i) saddle point contributions that depend on how rapidly the
function f (x) being expanded decays as | x |’ ∞ and (ii) contribu-
tions from the poles and other singularities of f (x). In the limit that
the degree n of the rational Chebyshev coe¬cient tends to in¬nity for
¬xed value of the “map parameter” L, one type of contribution will be
exponentially large compared to the other, and it is inconsistent (in the
Poincar´ sense of asymptotics) to retain the other. Boyd points out that
to optimize numerical e¬ciency, one should allow L to vary with the
truncation N of the Chebyshev series. Convergence is optimized when
N and L simultaneously tend to in¬nity in a certain way so that both
the pole and saddle point contributions are of equal order. This sort
of analysis does not explicitly use exponentially-improved asymptotics
of the Dingle-Berry-Olver sort. Nevertheless, hyperasymptotic thinking
” considering the role of terms that at ¬rst glance are exponentially
small compared to the dominant terms ” is absolutely essential to this
kind of numerical optimization.
A few other interesting studies of the role of exponential smallness
in numerical analysis have already been made. For example, the usual
second order di¬erential equation for the nonlinear pendulum, qtt +
sin(q) = 0, can be written as the equivalent system
pt = sin(q), qt = p (121)
where q is the angle of de¬‚ection of the pendulum with q = 0 when the
pendulum is standing (unstably) on its head and p is the momentum.
Hakim and Mallick [136] show that when the ¬rst equation is discretized
by a forward di¬erence and the second equation by a backward di¬er-
ence, the system (121) becomes what in dynamical systems theory, with
a slight change in notation, is called the “standard mapping”:
pn+1 = pn + „ sin(qn ), qn+1 = qn + „ pn+1 (122)
The usual numerical analysis description begins and ends with the
statement that this algorithm is second order accurate, that is, has

ActaApplFINAL_OP92.tex; 21/08/2000; 16:16; no v.; p.65
66 John P. Boyd



π 2π
0 π 2π 0 q
Figure 14. Phase plane for the nonlinear pendulum. Left panel: trajectories for the
exact solution. Right panel: trajectories when the di¬erential equation is integrated
by a ¬nite di¬erence scheme, i. e, the “standard mapping”, with „ = 1. The cross-
hatched area is the region spanned by a single chaotic trajectory. To avoid printer
overload, the individual iterates were erased and replaced by a uniform texture.

an error which is proportional to „ 2 . 4 Hakim and Mallick point out
that there are also changes which are exponentially small in 1/„ ” and
qualitatively di¬erent from the e¬ects of ¬nite di¬erencing which are
proportional to powers of the timestep.
These changes are easiest to explain by examining the trajectories
of the di¬erential equation and the di¬erence system in the phase plane
(Fig. 14). The nonlinear pendulum is an exactly integrable system, and
all trajectories are periodic. The closed curves in the phase plane rep-
resent side-to-side, small amplitude oscillations of the pendulum. The
open curves at top and bottom show trajectories in which the pendu-
lum swings through complete loops like a propeller. These two species of
trajectories are divided by the “separatrix”, which is a trajectory that
passes through q = 0, the unstable equilibrium, with zero momentum
p. The separatrix and trajectories near it are super-sensitive to per-
turbations because only a tiny additional amount of momentum will
su¬ce to push a large amplitude oscillation over the top and thereby
convert it into a propeller-like motion.
The di¬erence system, alias “standard mapping”, is not integrable
and has chaotic solutions. When the time step „ is very small, however,
one would expect that in some sense the di¬erence and di¬erential sys-
The one-sided di¬erences are only ¬rst order, but by eliminating p, one can
show that the system is equivalent to applying centered, 2d order di¬erences to the
second order di¬erential equation.

ActaApplFINAL_OP92.tex; 21/08/2000; 16:16; no v.; p.66
Exponential Asymptotics

tems would be close to one another. Indeed, trajectories away from the
separatrix are not drastically altered by the discretization; ¬nite di¬er-
ences give a good approximation to these trajectories of the nonlinear
pendulum. The neighborhood of the separatrix, however, dissolves into
chaotic motion. The incoming and outgoing separatrices from the equi-
librium split with a splitting angle at the point where the separatrices
cross of approximately [136], [180]

φ∼ exp ’ (123)
„3 „

The width of the region of chaos around the separatrices has a similar
exponential dependence.
The bland statement “the method is second order” or even a formal
expansion of the errors in powers of the timestep „ completely misses
the spawning of this region of chaos. Such exponentially small quali-
tative changes are perhaps less important in the real world, where the
original di¬erential equation probably has regions of chaotic motion
anyway, than in the idealized world of exactly integrable systems such
as the nonlinear pendulum, the Korteweg-deVries equation and so on.
Still, it reiterates the theme that hyperasymptotics is important to
numerical analysis.
Hakim and Mallick observe that the discretized system can be inter-
preted in two ways: (i) a second order accurate approximation to the
pendulum system or (ii) a fourth order accurate approximation to
a nonlinear di¬erential equation which is obtained by modifying the
pendulum equation by the addition of a higher derivative with a „ -
dependent coe¬cient. Although the second interpretation seems rather
arti¬cial, it is also illuminating. Weakly nonlocal solitary waves and
hydrodynamic boundary layers arise in this same way through addi-
tion of a higher derivative with a coe¬cient proportional to the small
parameter. The result of such a singular perturbation is that the power
series in the small parameter is divergent, and there are e¬ects which
depend on the exponential of the reciprocal of the small parameter.

17. Numerical Methods for Exponential Smallness or:
Poltergeist-Hunting by the Numbers, I: Chebyshev &
Fourier Spectral Methods

Because of the messiness of hyperasymptotic methods even for clas-
sical special functions and ordinary di¬erential equations, numerical
algorithms are important both as checks and as alternatives to hyper-
asymptotics. The exponential dependence on 1/ for << 1 cries out

ActaApplFINAL_OP92.tex; 21/08/2000; 16:16; no v.; p.67
68 John P. Boyd

for numerical schemes whose error also falls exponentially with the
number of degrees of freedom M . Fortunately, Chebyshev and Fourier
spectral methods [55, 91] and also Pad´ approximants [9, 19] have this
property. In this section, we shall discuss spectral methods while Pad´ e
algorithms are described in the following section.
However, when the unknown function f ( ) has only a divergent pow-
er series and also has contributions that lie beyond all orders in , both
the rate of convergence and (sometimes) the methodology are altered.
For example, when a function f (x) which is analytic on x ∈ [’1, 1] is
expanded in a Chebyshev series, the error decreases geometrically with
M, the truncation of the Chebyshev series. In other words, the error
EM = O(exp(’µM )) as M ’ ∞.
When f ( ) has only a divergent power series about = 0, the Cheby-
shev or other spectral series on an interval that includes this point
will lack geometric convergence. However, as long as the function is
in¬nitely di¬erentiable (with bounded derivatives), it is easy to prove
by integration-by-parts that the error must decrease faster than any
¬nite power of M . “C ∞ ” singularities do not defeat spectral methods,
but merely slow them down.
By using the method of steepest descents, one can often show that


. 14
( 24 .)