. 17
( 24 .)


series, which converge faster and also are much more resistant to round-
o¬ error. It is useful to illustrate how easily these expansions can be
derived to replace the standard divergent asymptotic series.

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

For example, the Stieltjes function S( ) satis¬es the ordinary di¬er-
ential equation
+ (1 + )S = 1 (138)
For paper-and-pencil or symbolic language calculation, the simplest
method is the “Lanczos „ -method”. He observed [177, 284] that if we
perturb the right-hand side of Eq.(138) by a polynomial of degree M ,
multiplied by an unknown constant „ , we can then solve this perturbed
equation exactly by a polynomial of degree M . (Instead of the usual
strategy of approximately solving the exact di¬erential equation, the
„ -method exactly solves an approximate di¬erential equation.) If the
perturbing polynomial is M , then the „ -method yields the ¬rst M
terms of the usual divergent power series in .
However, this is actually a rather stupid choice if the goal is uniform
accuracy on some interval ∈ [0, z] where z is a complex number. The
power function M is extremely non-uniform ” very small near the
origin, but increasing very rapidly away from it. The polynomial of
degree M (and leading coe¬cient of one) which is most uniform on

[0, z] is the shifted Chebyshev polynomial TM ( /z) [177, 178, 55].
Table VI is a short code in the symbolic manipulation language
Maple to solve the ordinary di¬erential equation
2dS —
+ (1 + )S = 1 + „ TM ( /z) (139)
through a Chebyshev „ ’ method. The most important feature of the
table is simply its brevity: all the necessary algebra is performed in
exact, rational arithmetic under the control of only nine lines of code!
The choice of Maple is arbitrary; the same calculation could be per-
formed with equal brevity in any of the other widely used symbolic alge-
bra languages including Mathematica, MACSYMA and Reduce. Boyd
[55, 63] gives many examples of problem-solving via spectral methods
in algebraic manipulation languages. Note that because all operations
involve polynomials, not transcendentals, the code also executes very
The Chebyshev-„ result is rather messy: polynomial in , rational
in z. However, the Chebyshev (or Chebyshev-like) expansion will obvi-
ously converge most rapidly when the expansion interval [0, z] is small
as possible for a given . This implies that for best results, one should
choose z = . Making this substitution not only optimizes accuracy for
a given , but also simpli¬es the result to a rational function of alone.
The M = 8 approximation, which is a polynomial of degree 7 over a
polynomial of degree 8, is
S8 ≡ 16 {4 + 124 + 1336 + 33 7 }/
2 3 4 5 6
+ 6168 + 12173 + 8955 + 1737

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

{256 + 8192 + 93184 2 + 473088 3 + 1108800 4
+1128960 5 + 423360 6 + 40320 7 + 315 8 } (140)

Fig. 18 shows that for small positive , this simple rational approxima-
tion is on the whole a lot more useful than either the superasymptotic
or hyperasymptotic series.
Chebyshev polynomial approximations are usually polynomials rather
than rational functions and are optimized for a particular line segment
in the complex plane. By computing symbolically, we have obtained an
approximation that is more complicated (because it is rational rather
than polynomial) but has the great virtue of being as accurate, for a
given , as the standard Chebyshev approximation of degree M along
the segment [0, ] even when is complex-valued.
In the previous section, we have already given the ordinary asymp-
totics of the Chebyshev coe¬cients of the Stieltjes function. However,
the comparison of convergent Chebyshev series with divergent power
series is not completely straightforward. The asymptotic series uses a
number of terms N which is inversely proportional to . What happens
if we compare the N -term Chebyshev series on the interval ∈ [0, 1/N ]
with the N -term optimally truncated power series for = 1/N ?
Through an elementary steepest descent analysis of the usual inner
product integrals for the coe¬cients of an orthogonal series,, one ¬nds
that the N -th Chebyshev coe¬cient for the series on [0, 1/N ] is (pre-
viously unpublished)

aN ∼ 2.98 N exp(’2.723N )
∼ 2.98 exp(’2.723/ ) (141)

(One can show that the error in truncating the Chebyshev series after
N terms is proportional to aN [55, 73].) Intriguingly, the errors for Pad´
approximants and for hyperasymptotics are of this same form:

|f ’ fN | ¤ exp(’q/ ), N ∼ O(1/ ) (142)

where the constant q > 0 depends on the precise Chebyshev, Pad´, e
or hyperasymptotic scheme used. There are likely deep connections
between these di¬erent families of approximations which are now only
dimly understood [51].
One can make more entertaining approximations by using other
spectral basis sets. For example, the rational Chebyshev functions are a
good basis set for the semi-in¬nite interval, x ∈ [0, ∞]. Boyd [54] gives
three examples in which the usual pair of series “ divergent series in 1/x
for large x and convergent power series for small x “ can be replaced
by a single expansion over the entire semi-in¬nite range. The examples

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

range from the K1 Bessel function, which has a pole at the origin, to
the J0 Bessel function, in which separate series multiply the sine and
cosine in the uniform approximations, to the ground state eigenvalue
of the quantum quartic oscillator as a function of the coupling constant
, which is a Stieltjes function with a factorially divergent power series
about = 0 [47, 19]. These uniform approximations are much compli-
cated and converge more slowly than the pair of Chebyshev series they
replace, but have the advantage of avoiding a conditional statement,
which is needed in the traditional approach to switch between large
and small x approximations.
Lastly, one must not overlook non-series alternatives. Schulten, Ander-
son and Gordon [277] have developed an e¬cient subroutine to evaluate
the Airy functions at arbitrary points in the complex plane. Instead of
using an asymptotic approximation for large | z |, they use a clever
optimized Gaussian quadrature to directly evaluate the integral repre-
sentations for Ai and Bi, even on Stokes™ lines. Their double precision
code, which is accurate to at least 11 decimal places for all | z | (with
use of the power series about z = 0 near the origin) employs a maximum
of just six quadrature points!
Detailed comparisons between high order hyperasymptotics and oth-
er methods of numerical approximation have not yet been carried out.
Still, the examples and illustrations above show that the comparison,
except perhaps for special cases, is likely to be unfavorable to high
order asymptotics.
Although hyperasymptotics look comparable to Chebyshev and Pad´ e
schemes when N ∼ O( ), the Chebyshev and Pad´ have the profound
advantage of converging as N ’ ∞ for ¬xed . Furthermore, these
approximations are built from ordinary polynomials whereas hyper-
asymptotic approximations are series of hyperterminants, which in turn
are approximated by series of hypergeometric functions.
There may be a few exceptions: problems where no alternatives
are available. In quantum chaology, Gutzwiller™s divergent series for
the quantum spectrum has been summed using resurgence[182, 39].
The practical result has been greatly improved energies for quantum
mechanics odd-shaped billiard tables “ idealized but popular for test-
ing theories [182, 39]. Another application is computing the zeros of the
Riemann zeta function, where resurgence has proved to be much better
than the best available competitor, the (un-resurgent) Riemann-Siegel
formula [34].
Berry™s amusing quote,“I am not expecting an early call”, is a frank
admission that most extensions of exponential asymptotics beyond the
lowest nontrivial term are arithmurgically useless. The proper use of
exponential asymptotics is to give insight. A sensible application is to

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

All Methods ev
sy Chebyshev Only
10 a

Higher Order
Chebyshev Only
0 0.1 0.2 0.3 0.4 0.5
Figure 18. A comparison of the rational „ -Chebyshev approximation S8 versus the
superasymptotic and hyperasympotic approximations for the Stieltjes function, S( ).
The three solid curves plot the errors for each of the three methods versus . If
acceptable accuracy for a given is a point in the unshaded region in the upper
left corner, all three methods are satisfactory. In the unshaded lower right region,
none of the three approximations is su¬ciently good (although such tiny errors
< 1 — 10E ’12 can be achieved by simply using a Chebyshev „ approximation of
higher order). The vertical shading “ most of the graph “ shows where the Chebysehv
approximation S8 , a polynomial of degree 7 divided by degree 8, is successful, but
the asymptotic series fail because their minimum error is larger than the required
tolerance. In the vertically-and-horizontally shaded area, both hyperasympotics and
S8 ( ) are successful. Finally, there is a tiny region of horizontal shading where only
hyperasymptotics is successful (though a Chebyshev approximation of higher order
than S8 ( ) would succeed). The hyperasymptotic errors were calculated using the
Berry-Howls scheme (where the errors are O(exp(’2.386/ )), but employing the
more accurate hyperasymptotic methods of later authors such as Olde Daalhuis
would not change the theme of the graph.

compute a small term that is also the leading term to approximate
some crucial feature of a problem, perhaps the lifetime of a quantum
bound state or a nonlocal solitary wave.

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

20. Hybridizing Asymptotics with Numerics

The hyperasymptotic scheme of Boyd [68] and the PKKS method
[171, 172] are both blends of analysis and numerics in the sense that
the ¬nal step, the determination of the proportionality constant which
multiplies the exponential of 1/ , requires a computation. However, the
prior analysis has reduced the problem to a very small calculation that
returns an answer as the product of a number with an analytical fac-
tor. This is far di¬erent than a brute force calculation that requires a
hundred times as much computer time to return only a number.
The ¬‚ow past a sphere or cylinder at small Reynolds number Re[264,
159, 161, 160, 98, 283]. has frustrated ¬‚uid dynamicists for over forty
years, but there has been, very recently, a partial breakthrough by
means of a hybrid numerical-asymptotic method [170]. The source
of pain is that these expansions are double series in powers of Re
and 1/ log(Re) or, de¬ning = 1/ log(Re), in powers of exp(’1/ )
and . Formally, one should include an in¬nite number of logarithmic
corrections to the drag coe¬cient before computing the ¬rst correc-
tion proportional to Re. For the ¬‚ow past the circle, however, Re >
(1/ log(3.70/Re))4 for all Re > 1/12000. (Real ¬‚uid ¬‚ows are typical-
ly at much larger Re.). A systematic scheme for the transcendentally
small terms is still an open problem. For the sphere, which is probably
the easier of the two, Chester and Breach conclude sadly “the expan-
sion is of practical value only in the limited range Re < 1/2 and that
in this range there is little point in continuing the expansion further”.
Kropinski, Ward and Keller [170] made the crucial observation that
if the outer (“Oseen”) problem is solved numerically, the numerical
solution will implicitly incorporate an in¬nite number of logarithmic
corrections. Better yet, the outer solution is independent of whether the
body is a cylinder, an elliptic cylinder, or some other smooth shape: a
single numerical solution provides a good answer to a whole spectrum
of body shapes. The inner solution di¬ers from shape to shape, but
is easy to calculate analytically. They have successfully applied this
same idea, outer-numerical/inner-analytical, to other problems with


. 17
( 24 .)